BIB-VERSION:: CS-TR-v2.0
ID:: CORNELLCS//TR93-1337
ENTRY:: 1993-10-14
ORGANIZATION:: Cornell University, Computer Science Department
LANGUAGE:: English
TITLE:: Pseudospectra of the Convection-Diffusion Operator
AUTHOR:: Reddy, Satish C.
AUTHOR:: Trefethen, Lloyd N. 
AUTHOR:: Pathria, Dimpy
DATE:: April 1993
PAGES:: 32
ABSTRACT::
The spectrum of the simplest 1D convection-diffusion operator is a discrete 
subset of the negative real axis, but the pseudospectra are regions in the 
complex plane that approximate parabolas. Put another way, the norm of the 
resolvent is exponentially large as a function of the Peclet number 
throughout a certain parabolic region. These observations have a simple 
physical basis, and suggest that conventional spectral analysis for 
convection-diffusion operators may be of limited value in some applications.
END:: CORNELLCS//TR93-1337
BODY::
Pseudospectra of the
Convection-Diffusion Operator
Satish C. Reddy
Lloyd N Trefethen*
Dimpy Pathria
TR 93-1337
April 1993
Department of Computer Science
Cornell University
Ithaca, NY 14853-7501
*Supp6rted by NSF Grant DMS-9116110
Pseudospectra of the
Convection-Diffusion Operator
Satish C. Reddy
Courant Institute of Mathematical Sciences
251 Mercer St.
New York, NY 10012
reddyOQ?dms.nyu.edu
Ltoyd N. Trefethen*
Department of Computer Science
Cornell University
Ithaca, NY 14853
LNT?cs.corndl.edu
Dimpy Pathria
Program in Computational and Applied Mathematics
Princeton University
Princeton, NJ 08544
dpG?acm.princeton.edu
Abstract. The spectrum of the simplest tD convection-diffusion operator is a
discrete subset of the negative real axis, but the pseudospectra are regions in the
complex plane that approximate parabolas. Put another way, the norm of the
resolvent is exponentially large as a function of the P?clet number throughout a
certain parabolic region. These observations have a simple physical basis, and
suggest that conventional spectral analysis for convection-diffusion operators
may be of limited value in some applications.
*Supported by NSF Grant DMS-911611O.
CONTEXTS
2.
3.
4.
5.
6.
7.
Introduction... 1
Eigenfunctions and pseudo-eigenfunctions... 5
The limiting case			oc... 8
Estimates of the pseudospectra... 11
Symmetrizability and a further estimate... 16
Finite difference and spectral discretizations... 19
Applications... 22
Appendix: formulas for arbitrary convection and diffusion coefficients... 28
References... 30
1. Introduction
Figure 1 illustrates the phenomenon that is the subject of this paper. Let
 denote the convection-diffusion operator
flu			u11 + u',
u(0)=u(d)=0,			(1.1)
acting in the Hilbert space L2[o,d].t For d = 40, the solid dots in the figure
show the first 27 eigenvalues of fi, a set of numbers on the negative real axis.
These eigenvalues are the points A E C where the resolvent of fi, (Al--Hfl)--H', has
a pole, which we may write as
II (Al --H fi)--H' I = 00.
This much is standard, but now consider the behavior of (Al --H fi)--H' at other
points A. Counting from the outside in, the figure plots the following contour
lines in the complex A-plane:
II(AI--Hfl)--H'II = 10', 102 i03,..., lO7.
(Section 6 explains how these curves were computed.) Outside the shaded region
fi bounded by the dashed "critical parabola." II(AI--Hfl)--H' is small. Inside that
region, however, it is hiige, attaining vaiues as great as l0? at points A that
are far from the eigenvaiues. In fact, for any A in the interior of the shaded
region, II(AI--Hfl)--H111 grows exponentially as dH 00. If d were 400 instead of
40, the figure would look much the same but witb the contours representing
1O1o,102o,...,107().
It is common prad-ice in both pure and applied mathematics to analyze
an operator by investigating its spectrum. Ilowever, Figure 1 shows that in the
case of convection-diffusion operators, there may be a great deal of structure
in the surface II(AI--Hfl)--H'II that is not revealed by the points where its height
is infinite. This is a reflection of the fact that these operators are non-normal.
* Convection-diffusion operators are often written with explicit convection and/or diffusion parameters
e.g. flu = vu1, +cu'. The results of this paper carry over to this alternative formulation if d is replaced
by the Peclet number dc/?. Details are given in the Appendix.
tT0 be precise the domain of fi is the subset of L2[o,dl of continuous functions on [o?dj which satisfy
the boundary conditions and possess a second derivative in L2[(), dj. Similar remarks apply to the
operators tO,oo) and fl(--Hoo,co) defined in Section 3.
critical parabola ReA = --H(ImA)2
A
complex A-plane
-5
-4			-3			-2			-1			0			1
Figure 1. Contour plot of the resolvent norm surface (Al--Hfl)-1 in the
complex A-plane (d = 40). The contours represent levels 101,102,.. .,to?,
and the dots are the eigenvalues. At each point A in the interior of the
region II (shaded), II(AI --H )-` grows exponentially as d 00. Equiva-
lently, this can be iitterpreted as a depiction of ?-pseudospectra of fl for
E=10?1,10?2,...,10--H7.
This means that they cannot be unitarily diagonalized, or to put it another way,
their eigenfunctions are not orthogonal. For the convection-diffusion problem,
the degree of non-normality grows exponentially with d. It follows that any
attempt to make quantitative estimates of the behavior of  by means of its
eigeiifunctions or eigenvalues is likely to lead to exponentially large constants.
Such estimates are of little use when d is large, and have no content at all that
is uniformly valid as d cc. Indeed, as far as spectra are concerned, the limit
dH cc is a singular one. For each d< cc, the spedrum of  is a discrete subset
of the negative real axis, but for the analogous operator [o,oc) on a semi-infinite
interval, the spectrum is the entire shaded region II.
This paper is devoted to an alternative analysis of a?nvection-diffusion
operators: the investigation of their pseudospectra, which are robust with respect
2
to the limit d H 00. Here is the definition:
DEFINITION. Let  be a closed linear operator in a Hilbert space TL
and let ? > 0 be arbitrary. The &pseudospectrum A(() is the set
of all AEiC such that
(i) II(AI--H)--H111 > ?--H1
or equivalently,
(ii) For any ?` > ?, A belongs to the spectrum of  + S for some
bounded operator g on ? with ??&? ?
By convention we wnte II(AI--H)--H11 = 00 if (AI--H)--H' is unbounded or nonex-
istent, i.e., if A is in the spectrum A(). An equivalent formulation of condi-
tion (ii) is to consider ??I ? and then take the closure of the resulting set.
In words: the pseudospectrum is the set of all points where the resolvent
norm is large, or equivalently, which are spectral values of a slightly perturbed
Thus Figure 1 can be described as a depiction of the boundaries of
operator.
the pseudospectra
A10?1 (), A10?2(). . . ., A10?7(),
and illustrates that the pseudospectra of ? are approximately parabolas. If 
were normal, the ?-pseudospectrum for each ? would be just the union of the
closed ?-balls about each eigenvalue.
The analysis of linear operators by means of the resoivent and its norm is an
old technique, especially among theoreticians; see, e.g., [24]. The specific idea of
considering the sets that we call pseudospectra was perhaps first introduced by
H. J Landau in 1975 [2(3,27], and has been subsequently and for the most part
independently employed by others including Varah [43] ?Viikinson [45], Demmel
[9], Chatelin [7], and especially Godunov, Kostin, Malyshev, and their colleagues
in Novosibirsk [16,28]. Apart from the work by Landan, these applications have
been primarily concerned with the effect of rounding errors on matrix compu-
tations. Most recently, ti? work in this vein by Trefethen and Reddy and their
colleagues has followed Landan's original view that pseudospectra also provide
information about the behavior of non-normal matrices and operators of a more
fundamental kind, unrelated to rounding errors or indeed to perturbations of
3
any sort. Apptications of this nature have been investigated in numericat lin-
ear algebra [14,32,39], the numerical solution of differential equations [11,23,36],
and fluid dynamics [35,42]. For introductions to these ideas, see [34,37,40]. A
survey of the theory and applications of pseudospectra is in preparation [41].
The results of the present paper can be summarized as follows. Though
the spectrum of  is a subset of the negative real axis, its pseudospectra are
large regions in the left half-plane bounded approximately by parabolas, as
depicted in Figure 1 Any point A at a distance dist (A, II) outside the critical
parabola lies outside all of the pseudospectra A?() with ? < dist (A, II), indepen-
dently of d; eqmvalently, II(AI--H)--H'II < J/dist(A,fl) (Theorem 4). Any point
A inside the parabola is an ?pseudo-eigenvalue for a value of ? that decreases
exponentially as a function of d; equivalently, II(A! --H grows exponen-
tially as dH oc (Theowm 5). Overall, although A() 74 A([0?)) as dH 00,
H Ac([0??)) as d oc for every ? > 0 (Theorem 6). Finally, the condi-
tion number of the basis of eigenfunctions of  grows exponentially as d H oc
(eq. (5.3)). Pseudospectra of ?mte-difference and spectral discretizations of the
convection-diffusion operator are considered in Section 6, and Section 7 discusses
some applications.
Convection-diffusion operators have a long history in applied mathematics
and numerical analysis. They are the canonical examples of non-self-adjoint
operators, or after discretization, of nonsymmetric matrices. For d? 00, (1.1)
becomes the simplest example of a singular perturbation problem, with a so-
lution involving a boundary layer at r = 0 (relative to the length of the full
interval). Such problems pose a natural challenge to analytical methods and to
numerical algorithms. Finite difference and finite element strategies for solving
them, proposed as early as 1955, have been continuously defined and refined
since then, resulting in the almost overwbelming selection of methods available
today [1,2,3,6,10,30]. In recent years convection-diffusion equations have also
assumed a new role in mimerical analysis, emerging as favorite test problems for
nonsymmetric iterative linear solvers and preconditioners [8,12,13,14,15,18,22].
Our results imply that in any of these applications, when convection is domi-
nant, predictions based on the exact spectrum are likely to be misleading, and
numerical methods based on the exact spectrum are likely to be suboptimal.
4
2. Eigenfunctions and pseudo-eigenfunctions
If we ignore for the moment the matter of boundary conditions, the oper
ator  multiplies the function cax by the factor
A = c?+?.			(2.1)
Conversely, for each A E C, there are two functions ??` that  multiplies by A,
namely those given by the roots of (2.1),?
cr'?=--H-21?-21VmA.			(2.2)
In stating these observations we are making implicit use of the factorization
--HAJ = (D--Ha+I)(D--Hcr?i),			(2.3)
where D is the differentiation operator (cf. [33]).
For any A and corresponding c'+ and a__ the function
= ?a+r?6o?x			(2.4)
a+--Ha
satisfies the condition for an eigenfunction with eigenvalue A in the interior of
[0, d] and the boundary condition at x =0. It satisfies the boundary condition
at x = d too provided ea+d = ?a?d, that is, (a+ --H a?)d = 2rin for some nonzero
n E Z. By (2.2), this amounts to the condition d#mA = 2r?n, and upon
squaring we obtain the following eigenvalues:
THEoREM 1. The spectrum of is A()=Un>cfAu?, with
--H $ ?2?2
4			d2
(2.5)
Thus A() is a discrete set of negative real numbers in the interval (--Hoc, 41).t
*In the confluent case a? = 1 the second function becomes ze?X/2 All of our results carry over
continuously to this confluent point so in the remainder of this paper for simphcity we shall ignore
this case and speak as if and 0 are always distinct.
t The calculation above does not prove Theorem 1 completdy but by more careful reasoning it can
be shown that each number A? is a simple eigenvalue of  and that there are no other points in the
spectrum.
5
x
0			2.5			5
Figure 2. The first two eigenfunctions of  (d = 5). The non-normality of
fl can be seen in the fact that u1 and u2 are not orthogonal.
By (2.2), the eigenfundions (2.4) can be written (after rescaling by a con-
stant)
= e?X/2sin(7rnx/d),			n=1,2,3			(2.6)
This is a sine wave times an exponentially decaying envelope, and the cases
n = 1,2 are plotted for d = 5 in Figure 2. The boundary layer at the left of
the domain can be interpreted as the result of leftward convection introduced
by the term u?. Note that this figure confirms that  is non-normal, since the
eigenfunctions are obviously not orthogonal.
If  were purely a diffusion operator, the term --H ? would be absent from
2
(2.2), o+ and a would be negatives of each other, and for every A, one of
ea+x and e??z would be decreasing as a function of x and the other would be
increasing. For our problem, however, there are choices of A for which both a+
and a_ lie in the left half-plane and thus both e?+? and 6??? are decreasing
functions. For the eigenfunctions (2.6), this occurs with Rea+ = Rea = 21.
More generally, it occurs if and only if a belongs to the strip
S=fa?C: --H1<Rea<0f,			(2.7)
since if a is one solution of (2.2), the other is --H1 --Ha. The corresponding region
6
a-plane
A = a2 + a
A-plane
Figure 3. The algebra of the critical parabola P. The value A belongs to II
if and only if both of the corresponding values a+ lie in the left half-plane,
or equivalently, if and only if either one lies in the strip 5.
in the A-plane is the image of 5 under the function A = a+a2, which we denote
by fi:
II = fAEC: A=a2+a, --H1?Wa?OJ.
(2.8)
The "critical parabola" that bounds fi is the image of the boundary of 5 under
the same function, whicb we can represent simply by
P = fAe C: A=a2+a, Rea=OJ
(2.9)
since Rea = --Ht maps onto the same parabola as Rea = 0. See Figure 3.
Suppose now that A is any complex number in the interior of II, so that
Rea+ <0 and R?ea <0. Then ?(x) decreases exponentially with x, so if d
is reasonably large, it follows that the boundary condition u(d) = 0 is nearly
satisfied, with an error of order epd, where
ii			= maxfRea+,Rea?J =			2'+ ??ea?+-21 .			(2.10)
See Figure 4. Thus ?(x) is "nearly an eigenfunction" of , though A may be
far from any of the exact eigenvalues.* This phenomenon is the essence of the
*The term ?pseudo-eigenfunction in the heading of this section is generaJly applied not to a function
like ?(r) that nearly satisfies the prescribed boundary conditions but to a function that satisfies
the boundary conditions exact]y and thus belongs to the domain of the operator under investigation.
It is a straightforward matter to add a small correction A? to ? so as to construct a perturbed
function ?			? of this kind, with I (Al --H )?1?II/Ij?II of order approximately e			is then an
&pseudo-eigenfunction of fl with ? of order approximately c?d.
7
A = --H0.15
_____________ ______________ H
10			20
Figure 4. For any A in the interior of II, the fundion ?x) = (e?+r --H
e ?X)/(a+ --H??) satisfies the eigenvalue condition and the boundary con-
dition at x = 0, and it satisfies the boundary condition at x = d up to an
error that decreases exponentially with d. Thus ?(x) behaves nearly like
an eigenfunction.
fl(ptd)
reason why the pseudospedra of  approximate parabolic regions in the left
half-plane.
Equation (2.10) suggests the back-of-the-envelope estimate
ll?,			--H1???0,			(2.11)
where il? denotes the region of the A-plane bounded by the parabola P? that
is the image of the strip in the &plane bounded by the lines Rea = ? and
Rea = --H1 --H ?. Sections 4 and 5 will make this estimate more prec?se.
3. The limiting case
Before estimating the pseudospectra of  more carefully, it is instructive to
consider the analogous operators defined on the intervals [0, cc) and (--Hcc, cc).
Intuitively speaking, the boundary conditions become ti(0) = u(cc) = 0 and
u(--Hcc) = ti(cc) = 0, respectivAy. Since we are working in L2 spaces, however,
the boundary conditions at icc are taken care of automaticjly and need not
be imposed explidtly. Thus formally, the two operators in question are [o,oo),
acting in L2[0,cc) with boundary condition u(0) = 0, and ?(--H??), acting in
L2(--Hcc, cc) with no boundary conditions.
The operator (--H???) is normal, and its spectrum and pseudospectra are
easily determined. In the following theorem A? denotes the closed disk fA E C:
8
IA < tJ and the sum of two sets A, BE C has the usual meaning: A+B = [A E
C: A=Aat%, AaeA, %cBJ.
THEOREM 2. The spectrum of (-oooc) j?
and for each ?> 0 the ?-pseudospectmm is
Equivalently, for A ? P the resolvent norm satisfies
II(AJ (?ocoo))?1II --H dist(A, P)
(3.1)
(3.2)
(3.3)
In words, the ?-pseudospectrum of (--Hoo,oc) is the parabola P together with a
cushion of thickness E.
Proof of Theorem 2. This result is standard (see e.g. [17]), but since the
proof is easy we shall give it. Suppose (AI--H(????))u = f for some A E C with
u, f ? L2(--Hoo, oo). Then by taking Fourier transforms we obtain
(Atw2--Hiw)u(w) =
The Fourier variable iw corresponds to a in (2.1), and by (2.6), the factor
A+w2 --Hiw = A--Ha2 --Ha is zero for some w C R if and only if A lies on the critical
parabola P. If this is so, then Al --H (--H???) annihilates the corresponding
function e:wx = e?z, and thus A c A((????)). On the other hand for A ? P, we
calculate
--H w5??P? ?(w)
1 _ 1
wS??Pi?IA+w2?iwI --H dist(A,P)
This proves (3.1) and (3.3), and the equivalence of (3.2) and (3.3) follows from
the definition of pseudospectra. I
9
On the semi-infinite interval, P is replaced by II:
THEOREM 3. The spectrum of [()?oo) is
--H II,
and for each ?> 0 the ?-pseudospectmm is
Equivalently, for A ? II the resolvent norm satisfies
(3.4)
(3.5)
I (Al [o?oc))?1 --H dist (A, II)			(3.6)
Moreover, if w([0?)) denotes the numerical range of [Ooo), then
= fi.
(3.7)
The numerical range ??([o?oo)) is defined as the set of all Rayleigh quotients
(u,[0?)u)/(u,u) for u in the domain of [Ooo), and dt.? denotes set closure.
Proof of Theorem 3. For any A in the interior of II, the function 9(x) of
(2.4) satisfies both boundary conditions and belongs to L2[O, oc). Thus ?(x) is
an eigenfunction of [o?oo) with eigenvalue A. Since the spectrum is closed and
contained in the closure of the numerical range [21, problem 214], this implies
II C A([0??)) c dfw([0??))J.
We now establish the converse inclusion ll7([o?oo)) C II, which implies
--H A([0??)) = fi and thus completes the proof of both (3.4) and
(3.7). By standard arguments to be found for example in [24], (3.7) is also
equivalent to (3.6), and as in the proof of Theorem 2, the equivalence of (3.6)
and (3.5) follows from the definition of pseudospectra.
Let u be a function in the domain of [o?oc) normalized by 1Jul = 1. By
integration by parts we calculate
(u, [O,oo)u) --H !o? u(x)(u11 t u') dx = !o? Ju'(x)J2dx + !o? u(x)u'(x)dx.
to
The first term is real, and another integration by parts shows that the second
is imaginary, so by the Cauchy-Schwarz inequality we obtain
Re(u,[0,?)u) = IIu?II2.,
Therefore
Im(u,[ooc)u)I ? lull IIu'II = llu'll.
Re(u,[0??)u) ? --HlIm(u,[0??)u)l2,
or by (2.9), (u,[00o)u)Cil. Thus H7([0?))cfl, as claimed.
Figure 5 summarizes the results of Theorems 2 and 3.
P
Figure 5. Spectra of (-oo?) and [ooo)
4. Estimates of the pseudospectra
Our first theorem (?ncerning the pseudospectra of  asserts that outside
the parabolic region fi, the resolvent nornt of  decreases inverse-linearly with
constant 1 in the numerator?i.e., the "Kreiss constant" is 1 (cf. [11]).
11
THEOREM 4. For each ? > 0, the ?-pseudospectrum of  satisfies
C fitAc
Equivalently, for A ? 11 the resolvent nonn satisfies
II(AI--H)?'II ? dist(1A,fl)'
and the numerical range satisfies
(4.1)
(4.2)
c II.			(4.3)
Proof. By changing the upper limit of integration from oc to d in the proof
of Theorem 3, we obtain a proof of (4.3), and this implies (4.2) and equivalently
(4.1). ?
We now turn to the problem of estin?ating the resolvent norm and hence
the pseudospectra insid?' II. To motivate our estimates, note that the resolvent
is the solution operator for the inhomogeneous problem
(AJ--H)u = f,			u(O)=u(d)=O			(4.4)
for f ? L2[o,d]. For each fixed A and y, the Green's function G(x,y) for this
problem is the solution u(x) couesponding to f(x) = 6(x y), and it can be
written
G(x, y) = a ?(x) t ?([x --H y]+),			(4.5)
where ?(x) is the function (2.4) illustrated in Figure 4, [x --H y]+ denotes 0 for
x--Hy<Oandx--Hyforx-?y>0 and aisaconstant. The term ?[x--Hy]+) would
be a solution to (4.4) if there were no right-hand boundary condition, and the
term a?(x) is a correction added to enforce that boundary condition. Figure 6
sketches G(x, y) for d = 20, y = 10, and thr?? values of A. The figure reveals that
in the case A = --H0.2, this correction term is exponentially large in 0 < x
For general functions f(x) in (4.4), we can make an analogous decomposi-
tion of u(x) into two terms. Assuming that A is in the interior of II, write
(4.6)
= R?+Rd,
12
G(x,y)
A=--HO.2
A=O
A=O.2
0			10
--H x
20
Figure 6. Green's functions G(x,?) for three values of A on the interval
[o,d], d=20, y= 10, normalized by G(y,y)= 1. The value A=--H().2 lies in
the interior of II, and this explains the exponentially large lobe of G(x, y)
for 0< x < ? in the fiTst plot.
where R? is the operator on L2 [0, d] defined as follows: extend the data by
zero to L2(--Hoo,oo), apply (AJ--H(--H??))--H', then restrict the result to [O,dj.
By construction IiR?II < Ii(AI--H(--H??))--H'II, so by (3.3) we have
IIR?II < dist(A,P)			(4.7)
This quantity is small compared with the exponentials of interest, and thus (4.6)
and (4.7) effectively reduce the estimation of II(AJ--H)--H111 to the estimation of
liRdli, which we can carry out with the aid of explicit representations of R? and
13
Rd. First, by integrating the term ?([x--Hy]+) in (4.5) against f(y), we obtain
the following formula for v = R?f:
In particular,
v(x) = $0Z?(x?y)j(y)dy.			(4.8)
v(d) = f??(?--Hy)?(y??y.			(4.9)
Since Rd is designed to enforce the boundary condition (R?f+Rdf)(d) = 0, we
have the following formula for w = R?f:
Combining (4.9) and (4.10) gives
w(x) =			9(x)			(4.10)
liRdil = sup liwil --H II?II Iv(d)I
f Ilfil I?(d)I S??P ???
(4.11)
By the Cauchy-Schwarz inequaiity applied to (4.9), the second supremum in
(4.11) is ???, attained with f(y) = ?(d--Hy), and thns we have
liRdil --H II?II2
I?d)I
(4.12)
Estimates for II(AI--H)--H'II can now be obtained by combining (4.7) and (4.12)
with the inequality
liRd I --H IIR?II < II(Al--H)--H'II < liRdli t IIR?II, (4.13)
a conseqnence of (4.6).
The following theorem presents several estimates of this kind. Note that
(4.16) shows that the back-of-the-envelope estimate (2.11) got the leading ex-
ponential behavior right but missed an algebraic factor of size o(Ic?I?1) =
o(IAI--H1/2) as A?+oo.
14
THEOREM 5. Let A be an arbitrary point in the interior of II, with ?+ and
?(x) defined by (2.2) and (2.4), and write 0+ --H 0 = ?= a + ir. Then
and
II(AJ--H)--H'II N IJ?II2
(4.14)
2
II?II2			(1?a2)(1 +r2)			(4.15)
as d H 00. If in addition A ? (?00,?4?], then ?(d) N e?d/1o+ --H ?--HI and
therefore
AI--H'			2e??d(a2+r2)l/2
--H			N (1 --H a2)(1 + r2) `			(4.16)
where? = maxfReo'+, Reo?J <0 as in (2.10).
Proof. By (4.7), jR?? is bounded independently of d as dH 00, whereas
it follows from (4.12) that li1?d11 grows without bound. Therefore (4.14) follows
from (4.12) and (4.13). Note that if A c (--Hoo, ??`), the terms on either side of
(4.14) oscillate through 00 between large positive and negative values each time
d passes through a zero of v(d) (making A an eigenvalue of ). Nevertheless,
the ratio of the two converges to 1 as d H 00, as the ` N'? requires.
To derive (4.15), we compute
11			dx			?0?Iea+??e??z12dx
= $0?(e(Y+Z --H ea?z)(e?+r --H eQ'?z)dx
= ?0?(6(?++Q-+)r + e(??+a?)z --H e??-+?+)? e(a?+?+)x) dx
1			1			1
+			+
o++&+ a?+o? a?+a% 0?+0+
1			1			1			1
--H 1--Ha + 1+ a --H 1+ ir --H 1--Hjr
15
2			2 2(a2+?2)
1--Ha2 1+r2 = (1--Ha2)(1+r2)
Since 4)(x) = (e?+x --H e??x)/(o+ --H c?) and I?+ --H a 12 = a2 + T2, this establishes
(4.15). From here, the derivation of ?(d) N ePd/ja+ --H a_ and then (4.16) is
trivial. o+
When d is reasonably large, the estimate (4.16) is a very good one. For
example, if Figure 1 is redrawn to plot contours of (4.16) instead of II(AI--H)--H'II,
all the curves except the outermost one (? = 10--H1) remain unchanged to plotting
accuracy.
Theorems 4 and 5 imply as a corollary that for each ? > 0, Ac() behaves
continuously in the limit d H 00. The following theorem to this effect is analo-
gous to Theorem 3.3 of [37].
THEoREM 6. For each A c C,
II(AI--H)--H1II H II(AJ--H[0?))--H1ii as dH00.
(4.17)
Proof. If A is in the interior of 11, the right-hand side of (4.17) is oc by
Theorem 3, and by Theorem 5, the left-hand side converges to 00 as dH 00.
Suppose on the other hand that A is in the exteriorof fi or on the boundary P.
The right-hand side of (4.17) is now 1/dist(A,fl), and by Theorem 4, the left-
hand sideis ?1/dist(A,fi)for alld. Toderiveareverseinequality,forany? >0
and ?>0, letA' beapointintheinterioroffi with jA'--HAj<dist(A,fl)+?, and
taked largeenough sothat II(A'I--H)--H'II >?`. Thenitiseasilyseenthat there
existsaperturbation,with `--H? <IA'--HAI+? <dist(A,fl)+?+? suchthat
A isin the spectrumof`, which implies Ii(AJ--H)--H1II >1/(dist(A,fl)+?+?).
Taking SiH0 and ?H0 completes the proof. ?
5. Symmetrizability and a further estimate
For each lixed d, the operator  is symmetrizable: similar by a diagonal
similarity transformation to a self-adjoint operator with, of course, the same
16
real elgenvalues. This is pointed out in Example III.6.1? of [24]? and is also
suggested by the form of the cigenfunctions (2.6).
To carry out the syinmetrization we define ?i(x) = e--Hz/2v(x), which implies
tit --H e?x2?2lv+vt,			till = e?xI2(4iv?vt+vll)
and therefore
Thus if we define
then we have
ti = till + tit --H c?x/2(vll --H 41 v).
--Hi			_
kv=v't			4v,			w&4v --H e?XI2v			(5.1)
=			(5.2)
Here K is a self-adjoint operator arid ?\4 is a diagonal operator with JIM I = 1,
JIM--H1I = edi2, and consequently
= JIMlI JIM--H1 II --H ?d/2 (5.3)
The notation t:(M) comes from numerical anJysis, where JIMlI JIM--H1 JJ is in-
terpreted as a condition number. On the face of it ?(Ai) is just the condition
number of the symmetrizing transformation Al, but a little thought shows that
it is also equal to the condition number of tlie basis of eigenfunctions (2.6) if
they are normalized to satisfy JJti?JJ = 1, for these eigenfunctions are related
by the same transformation Al to an ortlionormal set of eigenfunctions of k,
namely the suitably normalized sine functions.
Since K is self-adjoint, its resolvent norm is
1 1
JJ(AJ--HK)--H1JJ = dist(A,A(K)) = dist(A,A())
On the other hand the resolvents of  and k are related by
= (AI--HMKM--H1)-1 = M(A1--Hk)?1M--H1
Combining these formulas yields the following new bound:
*A good deal of information about convection-diffusion operators can bc found in [24]. See Section
111.2.3 and Examples 111.5.32, [11.6.11,111.6.20, VIII.l.19.
17
THEOREM 7. Foranyd>O andAEC,
--H1			edI2			ed/2
AJ--H dist(A,A()) --H IlmAt
(5.4)
It is interesting to compare this result with Theorem 5. That theorem,
as well as Figure 1 and our discussion up to this point, suggested that the
pseudospectra of  approximate parabolas. Theorem 7 reveals that although
this may be true in a large part of the A-plane, it cannot be true as lAl H 00.
Instead, for any fixed ? and d, Ac() must be contained in a strip of finite
(though typically very large) width:
Ac(fl) C fAcC: IJmAI<?ed/2t.
(5.5)
With hindsight we can see that this deviation from a parabola was foreshadowed
in the presence of the aigebraic factors O(IAI--H"2) in (4.16), and is visible in
Figure 1.
Because of the equivalence of conditions (i) and (b) in the definition of
pseudospectra, (5.4) can be interpreted as a statement about the sensitivity
of the spectrum A(): if  is perturbed to  + g, each eigenvalue changes by
at most ed/2II?II. In numerical analysis this bound is known as the Baner-
Fike Iheorem. To see how close it is to sharp, one can calculate the (absolute)
condition number of an individual eigenvalue A?, equal to the inverse of the
absolute value of the inner product of the corresponding normalized eigenvectors
of  and its adjoint. Such a calculation shows that the condition numbers are
?ed/2/d for fixed d as 7? H 00, 50 (5.5) is sharp up to a factor asymptotic to d.
One might be tempted to conclude from the above observations that since
fl is symmetrizable, with a spectrum whos? sensitivity to perturbations is finite,
there is no need to pay attention to non-normality of convection-diffusion opera-
tors or to distinguish spectra and pseudospectra. Relatedly, the symmetrization
(5.2) amounts to the observation that  is a self-adjoint operator in a space with
an exponentially weighted inner product, and thus ()flC might argue that since
non-normality is a norm-dependent property, it is of no consequence. However,
arguments like these miss a great deal. Theoretically, they miss the point that
18
the quantities in question are finite only for each fixed d;  is by no means
uniformly symmetrizable as d H 00, nor are the eigenvalue condition numbers
uniformly bounded. Practically, they miss the point that even for fixed d,
these quantities are huge, so that ignoring them is asking for trouble. Norm-
independence may seem an attractive property mathematically, but quantities
that are measurable in the laboratory are rarely norm-independent. Indeed,
one might take that as a principle: if it's norm-independent, it's probably non-
physical.
6. Finite difference and spectral discretizations
Pseudospectra of differential operators can be computed numerically by
finite difference, finite element, or spectral methods. In future work we hope to
prove convergence of some of these numerical approMmations. liere, we shall
just present two examples.
The simplest finite-difference discretizations of (1.1) involve 3-point ap-
proximations to the first and second derivatives. in practice it may be better
to use an irregular grid in the convection-dominated case, but for simplicity,
consider the regular grid
Xj=jh, h=?d j=O,1,...,AT (6.1)
AT'
for some N> 1 and the centered difference approximations
2VjtVj?J			____________
u"(xj) ? Vj+1--Hh2			U'(rj)			Vj+1--HVj?1
3			With this discretization and our			2h
with v			u(x.).			homogeneous boundary con-
ditions,  is approximated by the tridiagonal Toeplitz matrix L(N) of dimension
N --H 1 with entries
L(N) --H
= --H2h--H2			j4?i --Hh--H2i(2h)--H'.			(6.2)
Figure 7 plots numerically computed pseudospectra of L(N) for d =40 and
N =30. The parameters and contour lines are the same as in Figure 1, and the
19
2t
iF
I--H
-2
-5			-4			-3			-2			-1			0
Figure 7. A Jinite-difference approximation to Figure t with N = 30.
2
0
-1
-2
-5			-4			-3			-2			-1			0
Figure 8. A Chebyshev spectral approximation to Figure 1 with N = 30.
Six real eigenvahes lie off-scale to the left of the figure.
20
reader should compare these figures. Evidently A?(L(N)) approximates A(()
well near the origin of the A-plane, corresponding to smooth modes that are
well resolved by the discretization. For larger A the approximation deterio-
rates quickly and it is clear that one would have to take N  200 to achieve
a good match throughout the region displayed. As it happens, quite a bit is
known about the pseudospectra of Toeplitz matrices [37] In particular, the fact
that the pseudospectra A?(L(N)) approximate ellipses is related to the fact that
the symbol of a tridiagonal Toeplitz matrix maps circles about the origin onto
ellipses.
Better approximations to A?() can be achieved with spectral methods. In
Figure 8, the same pseudospectra are plotted for the case of a spectral colloca-
tion approximation LsYN) based on the Chebyshev grid
xj=2Id(ltcos(3))), j=09l,...,i'\ffi
(6.3)
For the details of such approximations see [5]. It is obvious at a glance that the
approximation A((Ls(pN)) A?() is accurate over a larger region of the A-plane
than in the finite-difference case. Indeed it appears to be `spectrally accurate"
over a sizable portion of the plane. Note that this effect happens automatically
when the differential operator is discretized; in no sense is the discretization an
explicit attempt to capture pseudospectra.
The spectral structure of LsYN) is complicated, by contrast, and in fact,
there are three real pairs of additional eigenvalues off-scale to the left of Fig-
ure 8, with absolute values approximately 7, 17, and 95. The properties of
the eigenvalues of spectral discretization matrices are of some importance, for
sometimes the location of the outlying eigenvalues determines the allowable step
size in discretization of time-dependent problems. A comparison of Figure 8
with Figure 1 suggests, however, that the location of the higher eigenvalues in
the spectrum is almost an accidental matter, with the essential approximation
phenomenon taking place in the pseudospectra. There is a natural reason for
this. If a discrete operator approximates a continuous operator to a certain
accuracy T, it is reasonable to expect that the c--pseudospectra with ?> O(T)
may be approximated nieaningfully, at least in regions near eigenvalues whose
eigenfunctions are well resolved by the grid, but there is little reason to expect
21
accurate approximations 4'below the level of truncation error," i.e., for ? <? r
and in particular for c = 0.
A fine point that comes up in computations such as these is the question of
norms. The quantity II(AJ?)--H' is defined with respect to the norm of L2[O, d],
whereas (Al--H L(Iv))--Hl and II(AI--H L8(pN))?l would appear to be defined with
respect to the discrete 2 norms on the grid (6.1) and (6.3), respectively. To
correct for this discrepancy, in our computation L(N) is actually replaced by
w(N)L(N)(?y(N))--H1 auct Ls(pN) by %pN)Ls(pN) (Ws(pN))?1, where w(N) and W?(pN)
are appropriate weighting matrices. As the effect on the plots is in the end
slight, we omit the details.
We can now explain how Figure 1 was computed. Experiments indicate
that for each A, after weighting as described above, both II(AI--H L(N))?1II and
II(AI?Ls(pN))?lII converge to the correct value II(AI--H)--H'II as N?oo. To
produce Figure 1 we used spectral discretizations Ls(pN) and chose values of iV
adaptively for each A until convergence to four relative digits seemed to be
achieved. The values of A used for this picture lay on a grid of 150>c 70 points
in the complex upper balfplane, so O(10?) such computations were involved
altogether. The largesi value of AT required for the range of A in Figure 1
was 72, with a more typical value being 32. The whole computation required
approximately 1.2 x 10? floating-point operations as measured by Matlab (four
hours on a SPARC 2 workstation).
7. Applications
In this final sedion we briefly examine applications of the pseudospectra
of , and relatedly, some of the limitations of spectral analysis for convection-
diffusion problems.
One common application of spectral analysis is to simplify a problem by
expressing it in the basis of eigenmodes, where the equation becomes diagonal.
The limitation of this idea for convection-diffusion problems is apparent from
(5.3): the basis of eigenmodes is exponentially ill-conditioned. If a function of
norm 0(1) is expanded in this basis, then as a rule the expansion coefficients
22
will be of order edi2, and the time evolution will be largely determined by the
evolving patterns of caiicellation among these coefficients. This is an indica-
tion that the eigenfunction basis is a poor one. It may be simpler and more
natural to work with a basis chosen for its approximation properties, unrelated
to eigenmodes. If edi2 is greater than the inverse of machine precision on the
computer, one will be forced to do so, or else the computation will fail because
of rounding errors.
Another application of spectra is the estimation of norms of evolution pro-
cesses. Suppose one wishes to study energy growth or decay for the time-
dependent linear problem
tit=ti,			?I(O)=?tu.			(7.1)
The solution is u(t) --H 61flu0, with norm governed by the quantity Iletnil. From
spectral analysis we exped that 11e1C11 should be related to a(), the specirat
abscissa of , whose vaiue we know from Theorem 1:
1 7r2
= sup ReA =			4?j2
AEA()
(7.2)
However, a precise connection between IIe?II and o?() again requires the intro-
duction of the exponentially large factor ?(?Ni) --H edI2 of (5.3):
ta(fl) < etC11 <
(7.3)
or by combining (7.2) and (7.3),
e?t/4?t1t2Id2 < 1etC1 < ?d/2--Ht/4--Ht7r2/d2 (7.4)
Figure 9 plots the three quantities of (7.4) for the case d = 40. ?rhat we see is
a great difference between the transient and the asymptotic behavior of IIetCII.
In the transient phase, of length 0(d), IIetCII is nearly constant. Only after
that does the curve bend around to match the slC)pe predicted by the spectrum.
Physically, it is easy to see the explanation for this behavior. A broad pulse
that begins at the upstream end of the in?erval may convect downstream with
23
loIu
10--H10
120			160
0			40			80
Figure 9. Comparison of IIe??II with upper and lower bounds based 011
the spectrum. Spectral information is meaningful asymptotically as t 00,
but fails to capture the transient.
relatively little diffusion for a time period 0(d), decaying rapidly in energy only
when it nears the downstream end.
This gap between transient and asymptotic behavior can be captured better
if one considers the pseudospectra of  as well as the spectra. One has the
following relationships:
IIe??I fort HoC H Ac() for? HO,
IIe?L fort--H> 0 H A?() for ?H??
Iletfl for finite t ,` A?() for finite ?.
(7.5a)
(7.5b)
(7.5c)
The first of these statem?nts is the association just discussed between the asymp-
totic behavior of Iletnil as t H Oc and the spectrum of . The second can be
justified by the equation
6tP(fl)--Ho(t) ? IIe?II <			(7.6)
as t H O? where Q() is the the numerical abseissa of , i.e., the largest real
part of the numerical range, which is very close to 0. It can be shown that the
24
numencal range is determined by the behavior of A?() in the limit ? H 00, and
this implies (7.5b). Finatly, (7.5c) can be made precise in various ways with the
aid of the contour integral
e1			2ir1i F etA(AI --H )--H` dA,			(7.7)
where F is any contour enclosing A(L). From (7.7) one can derive bounds on
Iletnil that depend on tiie location of the pseudospectra in the complex plane;
an example is the Kreiss Niatrix Theorem. We shall not give details here, but
refer the reader to [11,35,42,44].
Can one do better than loose bounds and predict Iletnil as a function of I
exactly? Certainly not, based on the spectrum and/or the numerical range
alone. With the knowledge of all the pseudospectra A?(), however, the answer
may be yes. In [19] it has been conjectured that for any two matrices or operators
? and 2 with Ac(i) A?(2) for all ?>O, one has IIf(1)II = IIf(2)II for all
functions f anatytic in a neighborhood of the spectra. As of this writing, this
conjecture remains unresolved.
All of these remarks apply to the pureJy linear, constant-coefficient problem
(7.1). Another limitation of spectral analysis is that as a rule it is not robust
with respect to the introduction of complications such as
o+ nonlinearity,
o+ variable coefficients,
o+ inhomogeneous fordug data,
o+ other perturbations.
If (7.1) is modified in any of these ways, the effect will typically be that the
transient part of Figure 9 remains approximately unchanged, but the long-time
part of the figure precisely the part related to eigenvalues may change en-
tirely. Such modification of a linear constant-coefficient problem is more the
rule than the exception, for many simple problems in the mathematical sci-
ences arise from more complicated problems by the familiar sequence of steps
of linearization, freezing of coefficients, and diagonalization:
nonlinear			linear			linear problem with			collection of
H			H			H
problem			prnblem			constant coeffients			scalar problems
25
In general, the first two of these simplifications are valid for short times and
the third is valid for I 00. When these time scales fail to overlap, spectral
analysis may be of litth value, even though there may be nothing wrong with
the linearized approximation or with analyzing its behavior by pseudospectra.
Among the complications just listed, the easiest to make precise is the
introduction of inhomogeneous forcing data. Suppose (7.1) is changed from an
initial-value problem on [0,00) to a problem on (--H00,00) driven by a periodic
forcing term:
ut=titeA?v,			--H00?t?00			(7.8)
for some A ? A() and some fixed function v in the domain of . (Values of A
on the imaginary axis correspond to forcing the system at a real frequency.) It
is easily verified that the response will be
u(t)=eAtu,			u=(AI--H?)?'v.			(7.9)
Thus the resolvent (Al --H )--H1 transforms "inputs" v in (7.8) to corresponding
"outputs" u, and the degree of amplification that may occur in the process is
given by
lull
sup =ll(AI--H)--H1ll.
v#() llvll
(7.10)
This equation implies that the resonant or "pseudo-resonant response of a
non-normal system to inputs at a frequenc? A may be large when A E A?() for
small E, regardless of whether or not A is close to the spectrum. In particular,
a convection-diffusion problem may have a larg&amplitude response to low-
frequency input, even tbough the spectrum lies to the left of the line ReA =
For further discussion of these matters, see [23], and for their application
to the physical problem of hydrodynamic stability, see [42].
We shall close with a discussion of a different type of application of the
ideas of this paper. Suppose we wish to solve a linear system of equations
L(N)u?f			(7.11)
by a nonsymmetric Chebyshev iteration of the kind investigated by Manteuffel
[29], where L(N) is a discrete approximation to  as considered in Section 6.
Such an iteration depends upon a pair of parameters a, b ? C, which define
26
an interval with respec? to which a family of Chebyshev polynomials [p?J is
implicitly constructed. The lemniscates Ipn(z)I = const. of these polynomials
approximate ellipses with foci a, b.
The conventional view is that a and b should be chosen so that the spectrum
A(L(N)) lies in one of these ellipses that is in an appropriate sense as small as
possible. This view is supported by the estimate analogous to (7.2),
IIr?I ?
II?oII --H IPnIIA(L(N))K(I?I),			(7.12)
where IIr?II/IIroII denotes the ratio of the nth to the initial residual norms,
IIPnIIA(L(N)) denotes the suprenium of pn(z)I for z E A(L(N)), and K(M) is the
condition number of a transformation that makes L(N) diagonal or self-adjoint.
The trouble with (7.12) is that the factor K(I?) may be exponentially large [37].
An alternative bound can be based on the norm of p? on a pseudospectrum,
rather than the spectrurn. F?or any ?> 0 we have
lirII			4
ir?oii?? IIPn IIA?(L(N)) 2r?'			(7.13)
where 4 denotes the arc length of ?Ac(L(N)); the proof is by a contour in-
tegral analogous to (7.7) [32]. NVhen ?(At) is large and ? is not too small, a
comparison of (7.12) and (7.13) suggests that when L(N) is far from normal,
the parameters of a nonsymmetric matrix iteration may be more appropriately
based on the pseudospectra than on the .9pectrum. Experiments confirm this
prediction dramatically for certain problerns [39]. A striking related example
is the phenomenon of direction-dependent convergence considered in [4,13,22]:
a Gauss-Seidel or SOR iteration for a convedion-diffusion problem may con-
verge in far fewer steps if one sweeps with the convection rather than against
it. This effect cannot l?e explained on tlie basis of spectra alone, which are
typically independent of the sweep direction, but the explanation becomes clear
when one looks at the pseudospectra. A similar gap between spectra and pseu-
dospectra appears if one considers the numerical stability of discrete methods
for convection-diffusion equations [20,31 ,3(i].
It is a curious phenomenon that choices of iteration parameters based on
adaptive or automatic estimates of the spectrum A(L(N)), such as those con-
structed by the Manteuffel algorithm, are often more reliable than the above
27
observations seem to suggest. The same is true of choices based on certain kinds
of a priori analysis of A(L(N)), as in a recent paper by Kellogg [25]. The ex-
planation is that most approximate techniques for estimating spectra both of
matrices and of operators?are robust enough that that they actually estimate
pseudospectra instead, whether or not that that was the intent of their au-
thors. This remark applies for example to Gerschgorin's theorem, the simplest
of all eigenvalue estimates; to the Arnoldi process, one of the more sophisticated
[14,38]; and to approximations based on energy methods, which are so robust
that they are typically valid for the numerical range as well as the spectrum and
the pseudospectra. Many numerical methods based on estimates of spectra of
nonsymmetric matrices would become markedly less reliable if those estimates
could magically be replaced by exact spectral information.
Appendix.
Formulas for arbitrary convection and
diffusion coefficients
The results in this paper have been stated for the nondimensionalized oper-
ator (1.1), with convection and diffusion coefficients both equal to 1. However,
this is no essential restriction. Suppose we are given the fully dimensional equa-
tion
au			a2u dU
__ =			= v			t
we have introduced a time variable for motivation. Then the substitutions
v			v
?=--Hx,			T=--H1
c			c2
reduce the problem to
c28u			c2&2ti &du
vdt =?U=vax2+vax?
With the further substitutions
9
Ar =ffi,
v
28
u(O,t) = ti(?c/v, I) = 0.
??ftv
c
we obtain
j2u a?
at			ax2+ax?			u(O,I)--Hu(d,t)=O,			(A.4)
the time-dependent analogue of (1.1). These reductions show that any convec-
tion-diffusion problem (A. 1) is equivalent to a problem that depends on a single
parameter d = ?4v, the P&Iet number.
By means of the transformations just described, the results of this paper
can be extended to arbitrary convection-diffusion operators A? of the form (A.1).
One simply makes the following substitutions indicated by (A.2) and (A.3):
?=;?,			,			d=?&.
1/
(No other changes should be made; in particular ?j and A do not change.) All
of the results in this paj?er now apply to A?.
For example, Theorem 1 describes tbe spectrum of ` i.e. the spectrum of
vA?/c2. Multiplying by c2/v shows that the eigenvalues of A? are the numbers
--H --H C2 --H 7r2fl2V
--H 4v			?2
The "critical parabola" liAr for ?? is likewise obtained by multiplying by c2/v:
li?=fA??:A= C2 ?? +?,Re?=OJ.
v
For another example, replacing  by v??/c2 in (4.2) and multiplying by v/c2
gives
which reduces to
vic2
II(Ac2I/v--HA?)?'II ? dist(A,fl)'
I (?J?Ar)?lII ?
if we set ? = c2A/v. As a final example, Th??rem 5 likewise remains valid, now as
a statement about the limit ?c/v? 00. The exponential growth of II(AJ--H)--H111
at each fixed point A corresponds to exponential growth of II(u1?A?)?lII at
points ?=c2A/v.
29
References
[1] R. E. Bank, J. F. Biirgler, W. Fichtner, and R. K. Smith, Some upwinding tech-
niques for finite element approximations of convection-diffusion equations, Nu-
mer. Math. 58 (1990), 185--H202.
[2] J. W. Barrett and K. W. Morton, Approximate symmetrization and Petrov-Galer-
kin methods for diffusion-convection problems, Comp. Meth. Appi. Mech. Engr.
45 (1984), 97--H122.
[3] G. Birkhoff, E. C. Gartland, Jr., and R. E. Lynch, Difference meThods for solving
convection-diffusion equations, Computers Math. Applic. 19 (1990), 147--H160.
[4] G. Birkhoff and R. E. Lynch, Some current questions on solving linear elliptic
problems, Numer. Math. 57 (1990), 527--H546.
[5]
[6]
C. Canuto, M. Y. Hussaini, A. Quarteroni, and T. A. Zang, Spectral Methods in
Fluid Dynamics, Springer-Verlag, New York, 1988.
M. A. Celia, T. F. Russell, I. Herrera, and R. E. Ewing, An Eulerian-Lagrangian
localized adjoint method for the advection-diffusion equation, Adv. Water Re-
sources 13 (1990), 187--H206.
[7] F. Chatelin and V. Fr;wsse', Analysis of arithmetic algorithms: a statistical study,
in P. Kornerup and I). Matula, eds., 10th IEEE Symposium on Computer Arith-
metic, Grenoble, France, 1991, pp. 10--H16.
[8] R. C. Y. Chin and T. A. Manteuffd, An analysis of block successive overrelaxation
for a class of matrices with complex sp?ctra, SlAM J. Numer. Anal. 25 (1988),
564--H585.
[9] J. W. Demmel, A cou?tterexampTh for two conjectures about stability, IEEE Trans.
Aut. Control AC-32 (1987), 340--H342.
[10] E. P. Doolan, J. J. II. Miller, and `AT. II. A. Schilders, Uniform Numerical Methods
for Problems wtih Initial and Boundary Layers, Boole Press, Dublin, 1980.
[11] J. L. M. van Dorsselaer, J. F. B. M. Kraaijevanger, and M. N. Spijker, Lin-
ear stability analysis in the numerical solution of initial value problems, Acta
Numerica, to appear.
[12] H. C. Elman, Iterative Methods for Large, Sparse, Nonsymmetric Systems of
Linear Equations, PhD diss. and Res. Rep. #229, Dept. Comp. Sci., Yale U.,
1982.
[13] H. C. Elman and M. 1). Chernesky, Ordering effects on relaxation methods applied
to the discrete one-dimensional convection-d'ffusion equation, SlAM J. Numer.
Anal., to appear.
[14] R. W. Freund, G. H. Golub, and N. M. Nachtigai, Iterative solution of linear
systems, in Acta Numerica 1992, Cambridge U. Press, 1992.
[15] R. W. Freund and N. M. Nachtigal, QMR: A quasi-minimal residual method for
non-Hermitian linear systems, Numer. Math. 60 (1991), 315--H339.
30
[16] 5. K. Godunov, A. (;. Antonov, 0. P. Kirityuk, and V. I. Kostin, Guaranteed
Accuracy of Solving Linear Systems of Equations in Euclidean Spaces, Nauka,
Novosibirsk, 1988 (Russian), to appear in English, Kluwer, 1993.
[17] I. C. Gohberg and I. A. Fet'dman, Convolution Equations and Projection Alethods
for their Solution, Amer. Math. Soc., Providence, 1974.
[18] C. I. Goldstein, Preconditioned iterative methods for convection diffusion and
related boundary value prnblems, Computers Math. Applic. 19 (1990), 11--H29.
[19] A. Greenbaum and L. N. Trefethen, Arnoldi/Lanczos and GMRES/GCR as ma-
trix approximation problems, SlAM J. Matnx Anal. Applics., to appear.
[20] D. F. Griffiths, I. Christie, and A. R. Mitchell, Analysis of error growth for
explicit difference schemes in conduction-convection problems, Int. J. Numer.
Meth. Engr. 15 (1980), 1075--H1081.
[21] P. R. Halmos, A Hilbert Space Problem Book, Springer-Verlag, New York, 1974.
[22] H. Han, V. P. ll'in, R. B. Kellogg, and W. Yuan, Analysis of flow directed itera-
tions, J. Comp. Math. (P. R. China) 10 (1992), 57--H76.
[23] D. J. Higham and L. N. Trefethen, Stiffness of ODEs, BIT, to appear.
[24] T. Kato, Perturbation Theory for Linear Operators, Springer-Verlag, New York,
1976.
[25] R. B. Kellogg, Spectral bounds and iterative methods in convection dominated
flow, in J. J. Miller, ed., Computational Methods for Boundary and Interior
Layers in Several Dimensions, Boole Press, Dublin, 1991.
[26] II. J. Landau, On SzegJ's eigenvalue distribution theory and non-Hermitian ker-
nels, J. d'Analyse Math. 28 (1975), 335?-357.
[27] H. J. Landau, Loss in unstable resonators, J. Opt. Soc. Amer. 66 (1976), 525--H529.
[28] A. Malyshev, An Introduction to Numerical Linear Algebra, Nauka, Novosibirsk,
1991 (Russian).
[29] T. A. Manteuffel, Adaptive procedure for estimating parameters for the nonsym-
metric Tchebychev iteration, Nurner. Math. 31(1978), 183--H208.
[30] J. J. H. Miller, ed., BAIL IV--HPrnceedings of the Fourth International Conference
on Boundary and Interior Layers, Boole Press, Dublin, 1986.
[31] K. W. Morton, Stabllity of finite differ??nce approximations to a diffusion-con-
vection equation, Int. J. Numer. Meth. Engr. 15 (1980), 677--H683.
[32] N. M. Nachtigal, 5. (;. Reddy, and L. N. Trefethen, How fast are nonsymmetric
matrix iterations?, SL?M J. Matrix Anal. Applics. 13 (1992), 778--H795.
[33] F. Nataf, Resolution de l'e'quation de convection-diffusion stationnaire par une
factorisation parabolique, C. R. Acad. Sci. Paris 310 (1990), 869-872.
[34] 5. C. Reddy, Pseudospectra of Wiener-Hopf integral operators and constant-
coefficient differential operators, J. 1nt. ?qs. Applics., submitted.
31
[35] 5. C. Reddy, P. J. Schmid and D. 5. Henningson, Pseudospectra of the Orr-
Sommerfeld operator, SlAM J. Appi. Math. 53 (1993), 15-47.
[36] 5. C. Reddy and L. N. Trefethen, Stability of the method of lines, Numer. Math.
62 (1992), 235--H267.
[37] L. Reichel and L. N. Trefethen, Eigenvalues and pseudo-cigenvalues of Toeplitz
matrices, Lin. Aig. Applics. 162--H164 (1992), 153--H185.
[38] Y. Saad, Numerical AIethods for Large Eigenvalue Problems, Manchester Univer-
sity Press, Manchester, UK, 1992.
[39] L. N. Trefethen, Approximation theory and numerical linear algebra, in J. C.
Mason and M. G. Cox, eds., Algorithins for Approximation II, Chapman and
Hall, London, 1990.
[40] L. N. Trefethen, Pseudospectm of Matrices, in D. F. Griffiths and G. A. Watson,
eds., Numerical Analysis 1991, Longman, 1992.
[41] L. N. Trefethen, Spectra and Pseudospectra: the Behavior of A'on-Normal Matri-
ces and Operators, book in preparation.
[42] L. N. Trefethen, A. E. Trefethen, 5. C. Reddy, T. A. Driscoll, A new direction in
hydrodynamic stability: beyond eigenvalues, Science, to appear.
[43] J. M. Varah, On the scparation of two matrices, SlAM J. Numer. Anal. 16 (1979),
216--H222.
[44] E. Wegert and L. N. Trefethen, From the Buffon needle problem to the Kreiss
matrix theorem, Am"r. Math. Monthly, to appear.
[45] J. H. Wilkinson, Sensitivity of Figenvalues H, Utilitas Math. 30 (1986), 243--H286.
32
