. 72
( 136 .)


imcreases simultaneously as the neighbor-to-neighbor distance decreases.

Fig. 17.1 is a schematic of the optimum relationship between h and N . If we keep h
¬xed, and merely use a larger and larger number of terms in (17.12), our approximation
will converge not to f (y) but to fW (y), and we are stuck with the ¬xed error |f (y) ’ fW (y)|
even with million grid points. Conversely, if we only decrease h and truncate the sum in
(17.12) by omitting all grid points such that |yj | > L for ¬xed L, we would be left with
an O(f [L]) error even if the grid spacing is ridiculously small. Because of the need to

simultaneously increase L and decrease h, the error is O(exp[’p N ]) for a typical function
like f (y) = sech(y): Subgeometric convergence with exponential index of convergence =
There are several alternatives for representing functions on in¬nite or semi-in¬nite in-
tervals, but they do no better. Sub-geometric convergence is a way of life when the interval
is unbounded. The grid of Gaussian quadrature points for the Hermite functions, for ex-
ample, automatically becomes both wider (in total span) and more closely spaced as N
increases as shown in Fig. 17.2. There is always this need to serve two masters ” increas-
ing the total width of the grid, and decreasing the spacing between neighboring points ”
for any method on an unbounded interval. The cardinal function series differs from or-
thogonal polynomial methods only in making this need explicit whereas it is hidden in the
automatic behavior of the interpolation points for the Hermite pseudospectral algorithm.
For more information about sinc numerical algorithms, consult the review article by
Stenger (1981) and the books by Lund and Bowers (1992) and Stenger (1993) or the articles
in Table 17.1.

Table 17.1: A Selected Bibliography of Sinc (Whittaker Cardinal) Applications

References Comments
Whittaker(1915) A few words from our founder: theory only
Stenger (1979) Boundary value problems
Lundin (1980) Nonlinear PDE boundary value problem
Stenger (1981) Review article; theory plus numerical applications
Saffman&Szeto(1981) Vortex stability
Lund&Riley(1984) Sinc with mapping to solve radial Schroedinger equation
Elliott&Stenger(1984) Numerical solution of singular integral equations
Robinson&Saffman (1984) Nonlinear eigenvalue problem in 2 space dimensions
Lund (1986) Symmetrization of sinc-Galerkin matrix
Boyd (1986b) Shows that exponential mapping of Stenger to in¬nite
interval is key to solving problems with
endpoint singularities on [’1, 1];
sinc basis is effective, but other
in¬nite interval bases work just as well
Schaffer&Stenger (1986) Multigrid iteration with sinc
Bowers & Lund(1987) Singular Poisson equations
Eggert&Jarratt&Lund(1987) Finite and semi-in¬nite intervals, too, through mapping;
singularities at end of ¬nite interval
Lewis&Lund&Bowers(1987) Sinc in both x and t for parabolic PDE
McArthur&Bowers&Lund(1987a) Second order hyperbolic PDEs
McArthur&Bowers&Lund(1987a) Numerical methods for PDEs
Bialecki(1989, 1991) Boundary value problems
Lund&Bowers&McArthur (1989) Elliptic PDEs
Smith&Bowers&Lund(1989) Fourth order problems in control theory
Jarratt&Lund&Bowers(1990) Sinc basis; endpoint singularities , ¬nite interval
Lund&Vogel (1990) Inverse problem for parabolic PDE
Boyd (1991f) Sinc solution of boundary value problems with a banded
matrix obtained by Euler sequence acceleration
Smith&Bogar&Bowers&Lund(1991) Fourth order differential equations
Lund&Bowers&Carlson(1991) Boundary feedback stabilization (control theory)
Smith&Bowers(1991) Inverse problem for Euler-Bernoulli beams
Smith&Bowers&Lund(1992) Inverse problem for Euler-Bernoulli beams
Boyd (1992c) Fast Multipole Method provides a fast sinc transform
Boyd (1992d) Fast off-grid sinc interpolation
Lund & Bowers (1992) Monograph
Stenger (1993) Monograph
Yin (1994) Poisson equation through tensor-product sinc basis
Koures (1996) Coulomb Schroedinger equation

Appendix F gives formulas for the derivatives of the Whittaker cardinal functions at the
grid points. The ¬rst derivative matrix is skew-symmetric while the second derivative is
symmetric (Stenger, 1979). These are desirable properties. Skew-symmetry implies that all
eigenvalues of the ¬rst derivative matrix are pure imaginary; thus the sinc discretization
of the one-dimensional advection equation, ut + cux = 0, is nondissipative. Similarly, the
symmetry of the second derivative matrix implies that all eigenvalues are real; the sinc
discretization of the diffusion equation gives a system of ODEs in time whose eigenmodes
all decay monotonically in time.
The sinc expansion is the easiest of all pseudospectral methods to program and to un-
derstand. One weakness, shared with the Hermite functions, is that the Whittaker cardinal
series is useful only when f (y) decays exponentially as |y| ’ ∞. In contrast, the T Bn (y)
will give rapid, exponential convergence for functions which decay algebraically with |y| or
asymptote to a constant, provided certain conditions are met. Thus, the mapped Cheby-
shev polynomials have a wider range of applicability than the sinc expansion.
A more serious weakness of sinc series is that ” unlike the mapped Chebyshev func-
tions ” they cannot be summed via the Fast Fourier Transform. This is especially unfortu-
nate for time-dependent problems and for iterative solution of multi-dimensional bound-
ary value problems.
A third limitation is that sinc functions do not give the exact solution to any classical
eigenvalue problem (to my knowledge). In contrast, the Hermite functions are the ide-
alized or asymptotic solutions to many physical problems. Thus, the Whittaker cardinal
functions, despite their simplicity, have not eliminated the competing basis sets discussed
in the next two sections.
The “Big Sky School” of Stenger, Lund, and their students have applied the sinc expan-
sion to problems on a ¬nite interval with endpoint singularities ” always in combination
with a map that exponentially stretches the bounded domain to an in¬nite interval. As ex-
plained in Chap. 16, Sec. 5, this is a good strategy for coping with endpoint branch points.
Eggert, Jarratt, and Lund (1987) also compute some solutions without endpoint singulari-
ties, but this is not recommended. Chebyshev series are much more ef¬cient when u(x) is
analytic at the boundaries.
In summary, it is best to regard the Whittaker cardinal functions as a good basis for
[’∞, ∞] only. Although it has some important defects in comparison to Hermite and
rational Chebyshev functions, the sinc series is the simplest and it converges just as fast as
the alternatives.

17.4 Hermite functions
The members of this basis set are de¬ned by

ψn (y) ≡ e’0.5 y Hn (y)

where Hn (y) is the Hermite polynomial of degree n. Like the familiar Chebyshev poly-
nomials, the n-th Hermite polynomial is of degree n. The members of the set satisfy a
three-term recurrence relation. The derivative of Hn (y) is simply (2 n Hn’1 ), so it is easy to
use Hermite functions to solve differential equations.
Unlike most other sets of orthogonal functions, the normalization constants for the Her-
mite functions increase exponentially with n, so it is dif¬cult to work with unnormalized
functions without encountering over¬‚ow (unless N is very small). In contrast, normalized
Hermite functions satisfy the bound
|ψn (y)| ¤ 0.816 all n, all real y (17.15)

Table 17.2: Hermite Function Theory & Applications

References Comments
Hille(1939, 1940a,b, 1961) Convergence theorems
Englemann&Feix 1D Vlasov eq. (plasma physics); Fourier in x
& Minardi&Oxenius(1963) Hermite in velocity coordinate u
Grant&Feix(1967) Fourier-Hermite for 1D Vlasov equation
Fokker-Planck collision term to improve Hermite convergence
Armstrong(1967) 1D Vlasov equation; Fourier-Hermite Galerkin
Hermite truncation reduced with increasing time to mitigate
slow convergence of Hermite part of tensor basis
Birkhoff&Fix(1970) Eigenproblems
Joyce&Knorr 1D Vlasov eq.; Fourier-Hermite basis compared with
&Meier(1971) double Fourier & method-of-moments (!)
Anderson (1973) Time-dependent model of the tropical ocean
Moore & Philander (1977) Equatorial waves and jets in the ocean (review)
Shoucri&Gagn´ (1977)
e 2D Vlasov; tensor Hermite basis in velocity
¬nite differences in space; plasma physics
Bain (1978) Convergence theorems
Banerjee (1978) Eigenvalue problem: quantum anharmonic operator
Banerjee et al. (1978) ψn (±y) with variable ±
Boyd (1980b) Rate of convergence theorems
Boyd (1984a) Asymptotic theory for Hermite coef¬cients through
steepest descent and the calculus of residues
Tribbia (1984) Hough-Hermite Galerkin method
Boyd & Moore (1986) Euler sequence acceleration of algebraically-converging
Hermite series with applications to oceanography
Marshall & Boyd (1987) Equatorial ocean waves
Smith(1988) Hermite-Galerkin model of tropical ocean
Funaro & Kavian (1991) Time-dependent diffusion problems
Weideman (1992) Hermite differentiation matrices are well-conditioned
Holvorcem(1992) Summation methods for Hermite series including
Holvorcem&Vianna(1992) Green™s functions, which depend on two variables
Boyd (1992c) Fast Multipole Methods provide fast Hermite transform
Boyd (1993) Hermite methods in symbolic manipulation languages
Tang (1993) Gaussian functions; role of the scale factor
Holloway(1996a,b) Vlasov-Maxwell equations (plasma physics)
Introduces “asymmetric” Hermite basis, which is very
ef¬cient for this application
Tse&Chasnov(1997) Hermite for vertical coordinate
Convection in unstable layer surrounded by
unbounded stable layers; Fourier in x, y
Schumer&Holloway(1998) 1D Vlasov equation; Hermite-Fourier basis
asymmetric Hermite is unstable; standard Hermite is good
with use of anti-¬lamentation ¬ltering

where the overbar denotes that the ψn (y) have been normalized. As for Chebyshev and
Fourier series, this bound (17.15) implies that the error in truncating a Hermite series is
bounded by the sum of the absolute values of all the neglected terms.
The optimum pseudospectral points are the roots of the (N +1)-st Hermite function. The
points are not exactly evenly spaced ” the outermost points are a little farther apart than
the separation between points near y = 0 ” but the asymptotic approximations to the Her-
mite functions show that the Hermite grid does not differ very much from the sinc grid
as shown in Fig. 17.2. In the previous section, we found that it was necessary to simul-
taneously (i) decrease the spacing between adjacent grid points and (ii) increase the span
between the outermost points on the grid if the sinc series was to give increasing accuracy
as N increased. Exactly the same √ thing happens with the Hermite grid: the separation be-

tween grid points decreases as 1/ N while the span of the grid increases as N . It is not
necessary to choose h and the span of the grid explicitly; using the roots of the Hermite
function as the grid points automatically insures the correct behavior.
Fig. 17.3 illustrates the ¬rst six (normalized) Hermite functions.
However, it is still necessary to pick a scaling factor ± since we can can always use
ψn (±y) as the basis set for any ¬nite ±. This freedom does not exist for a ¬nite interval, but
after a change of scale in y, an in¬nite domain is still an in¬nite domain. Banerjee (1978)
and Banerjee et al.(1978), Boyd (1984a), and Tang(1993) have discussed the choice of scale
factor. The theory (Boyd, 1984a) involves heavy use of steepest descent to asymptotically
approximate the Hermite coef¬cient integrals as n ’ ∞. One general conclusion is impor-
tant since it is also true of mapped Chebyshev methods: The optimum ± is a function of
N . A good choice of scaling parameter for N = 5 is a mediocre choice for N = 50; at least
for the cases studied, ± should increase with N .
The rates-of-convergence theory for Hermite functions is now well-developed for a ¬xed
(N -independent) scale factor ±. It is worth summarizing since it indicates the subtleties in-
troduced by the unboundedness of the expansion interval. We must begin with a de¬nition
because the rate of convergence depends not only on the location of singularities, as true of
all expansions, but also upon the rate at which the function f (y) decays as |y| tends to in¬n-

Figure 17.2: The pseudospectral grid for Hermite functions. For all N , the grid points are
the roots of the (N +1)-st Hermite polynomial. The grid spacing varies slightly, but is almost
uniform over the interval. The grid closely resembles the sinc grid: when N is quadrupled,
the span of the grid points is roughly doubled and the grid spacing is halved.

Figure 17.3: Graphs of the ¬rst 6 Hermite functions, ψn (y).

The ORDER of REAL AXIS DECAY k is the least upper bound of j for which

f (y) = O e’p|y|

for some constant p as |y| ’ ∞ along the real axis. The function is said to be “SUB-GAUSSIAN”
or “SUPER-GAUSSIAN” respectively if

[“SUB-GAUSSIAN”] (17.17)

[“SUPER-GAUSSIAN”] (17.18)

(i) The domain of convergence of a Hermite series is the in¬nite STRIP about the real y-axis bounded

|Im(y)| = w (17.19)

(ii) For an ENTIRE function, that is, f (y) which has no singularities except at in¬nity, the width is

∞ if k > 1
if k < 1

(iii) For SINGULAR functions, let „ denote the absolute value of the imaginary part of the location
of that singularity which is closest to the real y-axis, and let p denote the constant in the de¬nition
of real axis decay. Then
if k > 1

smaller of [p, „ ] if k = 1 (17.21)


. 72
( 136 .)