Pseudospectrum enclosures by discretization

A new method to enclose the pseudospectrum via the numerical range of the inverse of a matrix or linear operator is presented. The method is applied to finite-dimensional discretizations of an operator on an infinite-dimensional Hilbert space, and convergence results for different approximation schemes are obtained, including finite element methods. We show that the pseudospectrum of the full operator is contained in an intersection of sets which are expressed in terms of the numerical ranges of shifted inverses of the approximating matrices. The results are illustrated by means of two examples: the advection-diffusion operator and the Hain-L\"ust operator.


Introduction
Traditional stability analysis of linear dynamic models is based on eigenvalues. Thus determining the eigenvalues of a matrix or, more generally, the spectrum of a linear operator is a major task in analysis and numerics. The explicit computation of the whole spectrum of a linear operator by analytical or numerical techniques is only possible in rare cases. Moreover, the spectrum is in general quite sensitive with respect to small perturbations of the operator. This is in particular true for non-normal matrices and operators. Therefore, one is interested in supersets of the spectrum that are easier to compute and that are also robust under perturbations. One suitable superset is the ε-pseudospectrum, a notion which has been independently introduced by Landau [11], Varah [21], Godunov [10], Trefethen [18] and Hinrichsen and Pritchard [7]. The ε-pseudospectrum of a linear operator A on a Hilbert space H consists of the union of the spectra of all operators on H of the form A + P with P < ε. Besides the fact that the pseudospectrum is robust under pertubations, it is also suitable to determine the transient growth behavior of linear dynamic models in finite time, which may be far from the asymptotic behavior. For an overview on the pseudospectrum and its applications we refer the reader to [20] and [4].
Numerical computation of the pseudospectrum of a matrix has been intensively studied in the literature. Most algorithms use simple grid-based methods, where one computes the smallest singular value of A − z at the points z of a grid, or path-following methods, see the survey [19] or the overview at [4]. Both methods face several challenges. The main problem of grid-based methods is first to find a suitable region in the complex plane and then to perform the computation on a usually very large number of grid points. The main difficulty of path-following algorithms is to find a starting point, that is, a point on the boundary of the pseudospectrum. Moreover, as the pseudospectrum may be disconnected it is difficult to find every component. However, there are several speedup techniques available, see [19], which are essential for applications.
In this article we propose a new method to enclose the pseudospectrum via the numerical range of the inverse of the matrix or linear operator. More precisely, for a linear operator A on a Hilbert space and ε > 0 we show see Theorem 2.2. Here σ ε (A) denotes the ε-pseudospectrum of A, W ((A − s) −1 ) is the numerical range of the resolvent operator (A − s) −1 , B δs (U ) is the δ s -neighbourhood of a set U , and S is a suitable subset of the complex plane. This inclusion holds for matrices as well as for linear operators on Hilbert spaces. The idea to study the numerical range of the inverses stems from the fact that the spectrum of a matrix can be expressed in terms of inverses of shifted matrices [8]. From a numerical point of view this new method faces similar challenges as grid-based methods as a suitable set S of points has to be found and then the numerical ranges of a large number of matrices have to be computed. However, this new method has the advantage that it enables us to enclose the pseudospectrum of an infinite-dimensional operator by a set which is expressed by the approximating matrices. The usual procedure to compute the pseudospectrum of a linear operator on an infinite-dimensional Hilbert space is to approximate it by matrices and then to calculate the pseudospectrum of one of the approximating matrices. In [20, Chapter 43] spectral methods are used for the approximation, but no convergence properties of the pseudospectrum under discretization are proved. So far only few results are available concerning the relations between the pseudospectra of the discretized operator and those of the infinite-dimensional operator. Convergence properties of the pseudospectrum under discretization have been studied for the linearized Navier-Stokes equation [5], for band-dominated bounded operators [14] and for Toeplitz operators [2]. Further, Wolff [22] proves some abstract convergence results for the approximate point spectrum of a linear operator using the pseudospectra of the approximations.
In this article we refine the enclosure (1) of the pseudospectrum of linear operators further and show that it is sufficient to calculate the numerical ranges of approximating matrices. More precisely, we show in Theorem 3.5 that if n is sufficiently large. Here A n is a sequence of matrices which approximates the operator A strongly. We refer to Section 3 for the precise definition of strong approximation. If we even have a uniform approximation of the operator A, then we are able to prove an estimate for the index n such that (2) holds in intersections with compact subsets of the complex plane, see Section 4. In Section 5 we show that finite element discretizations of elliptic partial differential operators yield uniform approximations. Further, as an example of strong approximation we study in Section 6 a class of structured block operator matrices. In the final section we apply our obtained results to the advection-diffusion operator and the Hain-Lüst operator. We conclude this introduction with some remarks on the notation used. Let H be a Hilbert space. Throughout this article we assume that A : D(A) ⊂ H → H is a closed, densely defined, linear operator. We denote the range of A by R(A) and the spectrum by σ(A). The resolvent set is (A) = C\σ(A). Let L(H 1 , H 2 ) denote the set of linear, bounded operators from the Hilbert space H 1 to the Hilbert space H 2 . The operator norm of T ∈ L(H 1 , H 2 ) will be denoted by T L(H1,H2) . To shorten notation, we write L(H) = L(H, H) and denote the operator norm of T ∈ L(H) by T . The identity operator is denoted by I. For every λ ∈ (A), the resolvent (A − λ) −1 := (A − λI) −1 satisfies (A − λ) −1 ∈ L(H). For a set of complex numbers S ⊂ C we denote the δ-neighborhood by B δ (S), i.e., B δ (S) = {z ∈ C | dist(z, S) < δ}, and we also use the notation S −1 = z −1 z ∈ S \ {0} . Further, we use the notation C * := C\{0}.

Pseudospectrum enclosures using the numerical range
In this section we present the basic idea of considering numerical ranges of shifted inverses of an operator in order to obtain an enclosure of its pseudospectrum. We start by recalling the notions of the numerical range and the ε-pseudospectrum.
The numerical range of an operator A is defined as the set see e.g. [9]. It is always a convex set and, if A is additionally bounded, then W (A) is bounded too. Moreover the numerical range satisfies the inclusions where σ p (A) is the point spectrum of A, i.e., the set of all eigenvalues and σ app (A) is the so-called approximate point spectrum defined by The spectrum, point spectrum and approximate point spectrum are related by . If A has a compact resolvent, then the spectrum consists of eigenvalues only and hence we have equality. For ε > 0 the ε-pseudospectrum of A is given by If we understand (A − λ) −1 to be infinity for λ ∈ σ(A), then this can be shortend to The central idea of this article is the following: If λ ∈ C is such that 1/λ has a certain positive distance δ to the numerical range of the inverse operator A −1 , then this yields an estimate of the form with some constant ε > 0, which will in turn be used to show λ ∈ (A) with (A−λ) −1 ≤ 1 ε , i.e., λ ∈ σ ε (A). This is made explicit with the next proposition: Proposition 2.1. Suppose that 0 ∈ (A). Then for every 0 < ε ≤ Proof. Let us denote U = B δ (W (A −1 )) −1 . As a first step we show that So let λ ∈ C \ U . We consider two cases. First suppose that |λ| > In the other case if |λ| ≤ We have thus shown (3). In particular, λ ∈ C \ U implies λ ∈ σ app (A), i.e., Since B δ (W (A −1 )) is convex and bounded, the set C * \B δ (W (A −1 )) is connected and hence also the image under the homeomorphism C * → C * , z → z −1 . On the other hand, the boundedness of B δ (W (A −1 )) implies that a neighborhood around 0 belongs to C \ U = (C * \ U ) ∪ {0}. Consequently, the set C \ U is connected and satisfies 0 ∈ (A) ∩ C \ U . Using (4) and the fact that ∂σ(A) ⊂ σ app (A), we conclude that C \ U ⊂ (A). Here ∂σ(A) denotes the boundary of the spectrum of A. Now (3) implies that if λ ∈ C \ U then (A − λ) −1 ≤ 1 ε and therefore we obtain λ ∈ σ ε (A).
Applying the last result to the shifted operator A − s and then taking the intersection over a suitable set of shifts, we obtain our first main result on an enclosure of the pseudospectrum: Then for 0 < ε ≤ 1 2M we get the inclusion Proof. For every s ∈ S we can apply Proposition 2.1 to the operator A − s and obtain The following simple example demonstrates that the δ-neighborhood around the numerical range is actually needed to obtain an enclosure of the pseudospectrum.
Since A −1 is normal, its numerical range is simply the convex hull of its eigenvalues. Thus W (A −1 ) is the following square: Then, using the fact that z → 1 z is a Möbius transformation, we obtain for W (A −1 ) −1 the following curve plus its exterior: We see that W (A −1 ) −1 touches the spectrum of A. This is of course clear: if an eigenvalue 1/λ of A −1 is on the boundary of W (A −1 ), then the eigenvalue λ of A is on the boundary of W (A −1 ) −1 . In particular in this example we do not have σ ε (A) ⊂ W (A −1 ) −1 for any ε > 0 since σ ε (A) contains discs with radius ε around the eigenvalues.

A strong approximation scheme
In this section we consider finite-dimensional approximations A n to the full operator A. Our aim is to prove a version of Theorem 2.2 which provides a pseudospectrum enclosure for the full operator A in terms of numerical ranges of the approximating matrices A n ; this will allow us to compute the enclosure by numerical methods.
We suppose that 0 ∈ (A) and consider a sequence of approximations A n of the operator A of the following form: (a) U n ⊂ H, n ∈ N, are finite-dimensional subspaces of the Hilbert space H.
In this case we say that the family (P n , A n ) n∈N approximates A strongly. Note that (5) implies that n∈N U n is dense in H and that sup n∈N P n < ∞ by the uniform boundedness principle.
Lemma 3.1. Let U n , P n be such that (5) holds and let A n ∈ L(U n ) be invertible. Then the following assertions are equivalent: this shows the first part. For the second, let x ∈ D(A) and set y = Ax and x n = A −1 n P n y. Then x n → A −1 y = x and A n x n = P n y → y = Ax as n → ∞.
(b) ⇒ (a). Let y ∈ H. Set x = A −1 y and choose x n ∈ U n according to (b). Then Since both P n Ax → Ax and A n x n → Ax as n → ∞ and A −1 n is uniformly bounded, we obtain (a).
The following lemma shows that if A is approximated by A n strongly, then A−λ is approximated by A n −λ strongly too, provided (A n −λ) −1 is uniformly bounded in n.
Proof. This follows immediately from Lemma 3.1 since We now prove a convergence result for the numerical range of the inverse operator under strong approximations.
Proof. (a) We set y n = P n x/ P n x . Note that y n is well defined for almost all n since P n x → x = 1. We get y n → x as n → ∞ and which yields the assertion.
For every j we have z j = A −1 x j , x j with some x j ∈ H, x j = 1, and by (a) there exists n j ∈ N such that for all n ≥ n j there is a y j ∈ U n , y j = 1 such that Hence for all n ≥ n 0 = max{n 1 , . . . , n m }.
The previous lemma allows us easily to prove an approximation version of the basic enclosure result Proposition 2.1.
Proof. By Proposition 2.1 we have for n ≥ n 0 and the proof is complete.
Combining the previous proposition with shifts of the operator, we get our second main result. It is analogous to Theorem 2.2, but provides an enclosure of the pseudospectrum of the infinite-dimensional operator in terms of numerical ranges of the approximating matrices.
Proof. In view of Lemma 3.2, Proposition 3.4 can be applied to every A − s j . Hence there exists n j ∈ N such that Since σ ε (A) = σ ε (A−s j )+s j , the claim follows with n 0 = max{n 1 , . . . , n m }.

A uniform approximation scheme
In this section we pose additional assumptions on the approximations A n of the infinite-dimensional operator A, that will allow us to estimate the starting index n 0 for which the pseudospectrum enclosures from Proposition 3.4 and Theorem 3.5 hold on bounded sets. Throughout this section we assume that A has a compact resolvent, 0 ∈ (A) and that D(A) ⊂ W ⊂ H where the Hilbert space W is continuously and densely embedded into H. The closed graph theorem then implies A −1 ∈ L(H, W ). Further, we suppose that there is a sequence of approximations of the operator A in the following sense: (a) U n ⊂ H, n ∈ N, are finite-dimensional subspaces of H.
(b) There exist projections P n ∈ L(H) onto U n , n ∈ N, not necessarily orthogonal, with sup n∈N P n < ∞ and (I − P n )| W L(W,H) → 0 as n → ∞.
(c) There exist invertible operators A n ∈ L(U n ), n ∈ N, such that A −1 − A −1 n P n → 0 as n → ∞. We say that (P n , A n ) n∈N approximates A uniformly. For (I − P n )| W L(W,H) we will write abbreviatory I − P n L(W,H) . indeed A −1 is the uniform limit of the finite rank operators A −1 n P n and hence compact.
(b) If (P n , A n ) n∈N approximates A uniformly, then also strongly. Note here that from (b) we first obtain P n x → x for x ∈ W , which can then be extended to all x ∈ H by the density of W in H and the uniform boundedness of the P n . One particular consequence of the strong approximation is sup In order to obtain improved enclosures of the pseudospectrum under a uniform approximation scheme, that is, addtional estimates of the starting index n 0 for which the pseudospectrum enclosures from Proposition 3.4 and Theorem 3.5 hold on bounded sets, we refine the results from Section 2 in terms of certain subsets of the full numerical range of A −1 . For d > 0 we define Clearly Proof. (a) Let λ ∈ σ(A) with |λ| ≤ L. Then there exists x ∈ D(A) with x = 1 and Ax = λx. This implies and thus we obtain Consequently and first show Let λ ∈ B L−ε (0) \ U , x ∈ D(A), x = 1. We consider three cases. Suppose first that λ > In the second case assume x W ≥ d. Then Finally if λ ≤ 1 2 A −1 , the same reasoning as in the proof of Proposition 2.1 yields once again that (A − λ)x ≥ ε, and therefore (9) is proved. Now, since A has a compact resolvent (9) implies that From Proposition 4.2 we get again a shifted version: For 0 < ε ≤ 1 2M0 , L > ε, d = LM 1 and δ s = 2 (A−s) −1 2 ε we get the inclusion Proof. Apply Proposition 4.2(b) to A − s for all s ∈ S and note that For a uniform approximation scheme, the numerical range of A −1 can now be approximated with explicit control on the starting index n 0 : Lemma 4.5. Suppose that (P n , A n ) n∈N approximates A uniformly. Let (a) If d > 0, 0 < δ ≤ C0 2 and n 0 ∈ N are such that for every n ≥ n 0 A −1 − A −1 n P n + dC 0 I − P n L(W,H) < δ, then (b) If δ > 0 and n 0 ∈ N are such that for every n ≥ n 0 we have Proof. Let x ∈ W with x = 1 and x W ≤ d. Then we obtain as well as Let n ≥ n 0 . Then |1 − P n x | ≤ d I − P n L(W,H) < δ C 0 ≤ 1 2 and hence P n x ≥ 1 2 . Let x n = Pnx Pnx . Then x n = 1 and

This implies
and thus for n ≥ n 0 we arrive at

This yields
n )) if n ≥ n 0 and proves (a). In order to show part (b), let x ∈ U n with x = 1. As x = P n x we have Corollary 4.6. If (P n , A n ) n∈N approximates A uniformly, then or, equivalently, Proof. We first show the inclusion "⊃". Let (λ n ) n∈N be a convergent sequence in C with λ n ∈ W (A −1 n ) and define λ = lim n→∞ λ n . Let δ > 0 be arbitrary. Lemma 4.5(b) implies that there exists n 0 ∈ N such that λ n ∈ B δ (W (A −1 )) for every n ≥ n 0 . This implies λ ∈ B δ (W (A −1 )) for every δ > 0, and thus λ ∈ W (A −1 ).
The last result shows that W (A −1 ) can be represented as the pointwise limit of the finite-dimensional numerical ranges W (A −1 n ). Lemma 4.5 even yields a uniform approximation, but this is asymmetric, since one inclusion only holds for the restricted numerical range W (A −1 , d). A more symmetric result is discussed in the next remark: Remark 4.7. If U n ⊂ W for some n ∈ N then, due to the fact that the space U n is finite-dimensional, Using the same reasoning as in the proof of Lemma 4.5(b), we then obtain Note however that for finite element discretization schemes the condition U n ⊂ W will usually not be fulfilled. In our examples for instance U n are piecewise linear finite elements while W ⊂ H 2 (Ω) is a second order Sobolev space, and thus U n ⊂ W .
Under a uniform approximation scheme the pseudspectrum can be approximated as follows.
Proposition 4.8. Suppose that (P n , A n ) n∈N approximates A uniformly. Let If we choose n 0 ∈ N such that for every n ≥ n 0 where C 0 is defined in (10), then we obtain Next note that We can therefore apply Lemma 4.5 with δ replaced by δ − δ and n 0 chosen as stated above and obtain and hence the assertion.
5 Finite element discretization of elliptic partial differential operators As an example for a uniform approximation scheme defined in Section 4 we now consider finite element discretizations. We use the standard textbook approach via form methods, which can be found e.g. in [1,17]. Let V and H be Hilbert spaces with V ⊂ H densely and continuously embedded. In particular there is a constant c > 0 such that Moreover, we consider a bounded and coercive sesquilinear form a : V × V → C, that is, there exists constants M, γ > 0 such that Let A : D(A) ⊂ H → H be the operator associated with a, which is given by Then A is a densely defined, closed operator with 0 ∈ (A) and A −1 ≤ c 2 γ , where c > 0 is the constant from (11).
Let (U n ) n∈N be a sequence of finite-dimensional subspaces of V which are nested, that is U n ⊂ U n+1 . We denote by a n := a| Un the restriction of a from V to U n . The form a n is again bounded and coercive with the same constants M and γ. Let A n ∈ L(U n ) be the operator associated with a n , i.e.
a(x, y) = A n x, y , x, y ∈ U n .
Then again 0 ∈ (A n ) and A −1 n ≤ c 2 γ . Let P n ∈ L(H) be the orthogonal projection onto U n . Thus P n = 1 and A n = P n A n+1 | Un , that is, A n is a compression of A n+1 .
To obtain a uniform approximation scheme, we now consider an additional Hilbert space W which is densely and continuously embedded into H such that D(A) ⊂ W ⊂ V . We assume that there exists a sequence of operators Q n ∈ L(W, V ) with R(Q n ) ⊂ U n and hold. In particular, the family (P n , A n ) n∈N approximates A uniformly.
Proof. For w ∈ W we calculate which shows the first assertion. Moreover, for f ∈ H we set x = A −1 f and x n = A −1 n P n f . Then we obtain a(x, y) = Ax, y = f, y , y ∈ V, a(x n , u) = A n x n , u = P n f, u = f, u , u ∈ U n .
Using the Lemma of Cea [17, Theorem VII.5.A], we find which implies the second assertion.
Theorem 5.2. Let A be the operator associated with the coercive form a and let A n , Q n be as above. Let If n 0 ∈ N is such that for every n ≥ n 0 for all n ≥ n 0 .
Proof. We check that the conditions of Proposition 4.8 are satisfied: Using Lemma 5.1, we estimate for n ≥ n 0 and with C 0 from (10), where a ij ∈ C 0,1 (Ω) and b i , c ∈ L ∞ (Ω). We suppose that a is coercive and uniformly elliptic. Let {T n } n∈N be a family of nested, admissible and quasiuniform triangulations of Ω satisfying sup T ∈Tn diam(T ) ≤ 1 n . Let equipped with the H 2 -norm, and U n = u ∈ C 0 (Ω) u| T ∈ P 1 (T ), T ∈ T n , u| Γ = 0 , n ∈ N.
Here P 1 (T ) denotes the set of polynomials of degree 1 on the triangle T . We get U n ⊂ V . Moreover, the operator A associated with a is given by By the Sobolev embedding theorem we have H 2 (Ω) → C 0 (Ω). For u ∈ W we define Q n u as the unique element of U n satisfying (Q n u)(x) = u(x) for every vertex of the triangulation T n . Then Q n ∈ L(W, V ) with R(Q n ) ⊂ U n . Moreover, [1,Theorem 9.27] implies that there is a constant K > 0 such that We conclude that Theorem 5.2 can be applied in this example with n 0 ∈ N chosen such that Note that in Example 5.3 we can also consider Ω to be an open interval in R. All results continue to hold in an analogous way.

Discretization of a structured block operator matrix
In this section we investigate discretizations of a certain kind of block operator matrices. We consider block matrices of the form Additionally we assume that 0 ∈ (A), 0 ∈ (D) and that both A and −D are uniformly accretive, i.e., there exist constants γ A , γ D > 0 such that In the next lemma we show that under the above assumptions there is a gap in the spectrum of A along the imaginary axis, and we also prove an estimate for the norm of the resolvent. Similar results were obtained in [12,13] under the additional assumption that A is sectorial and, in [13], without the condition that B and D are bounded. However, no corresponding resolvent estimates were shown. We remark that the boundedness of D is not essential in Lemma 6.1 but will be used thereafter.
Proof. Consider the block operator matrix A simple calculation shows that for λ ∈ U := {λ ∈ C | −γ D < Re λ < γ A } and x ∈ D(A), y ∈ H, where c λ = min{γ A − Re λ, γ D + Re λ}. It follows that and therefore, since Jw = w for all w ∈ H × H, In particular λ ∈ σ app (A), i.e., U ∩ σ app (A) = ∅. The adjoint of A is the block operator matrix which also satisfies the assumptions of this lemma. Indeed, (16) obviously also holds for D * . Moreover, the uniform accretivity (15) of A together with 0 ∈ (A) imply that A−γ A is m-accretive, see [9, §V.3.10]. This in turn yields that A * −γ A is m-accretive too and hence It follows that (17) also holds for A * . In particular ker A * = {0} or, equivalently, R(A) ⊂ H × H is dense. On the other hand, (17) implies that ker A = {0} and that R(A) is closed. Consequently R(A) = H×H and therefore 0 ∈ (A). Using ∂σ(A) ⊂ σ app (A) and the connectedness of the set U , we obtain U ⊂ (A). Now (17) implies (A − λ) −1 ≤ 1/c λ for all λ ∈ U .
We consider approximations A n of A of the form where (a) (P n , A n ) n∈N is a family which approximates A strongly in the sense of Section 3; (b) all projections P n are orthogonal and all A n are uniformly accretive with the same constant γ A as in (15); (c) B n = P n B| Un , D n = P n D| Un where U n = R(P n ) (b) (P n , A n ) n∈N approximates A strongly where P n = diag(P n , P n ).

Proof. (a) From
it follows that −D n is uniformly accretive with constant γ D from (16). Consequently Lemma 6.1 can be applied to A n .
(b) In view of (a) and Lemma 3.1 it suffices to show that for all (x, y) ∈ D(A) × H there exist (x n , y n ) ∈ U n × U n such that Let (x, y) ∈ D(A) × H. From Lemma 3.1 we get x n ∈ U n with x n → x and A n x n → Ax as n → ∞. Set y n = P n y. Then y n → y and D n y n − Dy ≤ P n (Dy n − Dy) + P n Dy − Dy ≤ Dy n − Dy + P n Dy − Dy → 0, n → ∞, i.e., D n y n → Dy. The proof of B n y n → By and B * n x n → B * x is the same after the additional observation B * n = P n B * | Un . Hence we have shown (18). Then there exists n 0 ∈ N such that Proof. Lemma 6.1 and Lemma 6.2 imply and hence the assertion follows from Theorem 3.5.
Remark 6.4. Suppose that A is the operator associated with a coercive sesquilinear form a on V ⊂ H and that U n , W , P n ∈ L(H), A n ∈ L(U n ) are chosen as in Section 5. Then (P n , A n ) approximates A uniformly, and hence also strongly, see Remark 4.1. Moreover, the coercivity of a implies that A and all A n are uniformly accretive with constant γ A = γ from (12). Hence all assumptions of this section are fulfilled in this case.

Numerical examples
In order to exemplify the previously developed theory we take a look at the results of numerical computations. We investigate the steps that were involved in the discretisation of a given operator and describe a visualisation of supersets of the pseudospectrum.
Example 7.1. In this example we will examine the Hain-Lüst operator which fits into the framework of section 6. See [15] and [16] for results on the approximation of the quadratic numerical range of such a block operator. The Hain-Lüst operator under consideration here is defined by Let {T 1 n } n∈N be the family of decompositions of the interval (0, 1) where every subinterval T ∈ T 1 n is of length 1 n and let U n = {u ∈ C(0, 1) | u| T ∈ P 1 (T ), T ∈ T 1 n , u(0) = u(1) = 0}, n ∈ N. Here P 1 (T ) denotes the set of polynomials of degree 1 on the subinterval T . The piecewise linear functions for i ∈ {1, . . . , n − 1} form a basis of U n and therefore the functions for i ∈ {1, . . . , 2(n − 1)} form a basis of U n × U n . Evaluating (19) on these basis functions, the finite-element discretisation matrices A n of A are given by  Figure 1 for shifts s 1 , . . . , s m ∈ ρ(A). The choice of the shifts was determined by the expected shape of the pseudospectrum aiming to obtain a relatively small superset thereof. They are located on two circles around −3 with radii greater and smaller than 2 and on lines parallel to the real axis in the right half plane. Here n = 800, δ j = 2.1 (A n − s j ) −1 2 ε and ε = 0.9 2 maxj=1,...,m (A−sj ) −1 . The red dots are the eigenvalues of A n while the black lines correspond to the boundaries of the pseudospectrum of the approximation matrix σ ε (A n ) computed by eigtool, see [3]. Note that according to Theorem 3.5 the intersection of the blue areas form an enclosure of the pseudospectrum of the actual operator σ ε (A), while the black lines only give the information for the discretized operator. Furthermore the spectral gap mentioned in Lemma 6.1 becomes visible.
Example 7.2. Let us consider the the advection-diffusion operator A : D(A) ⊂ L 2 (0, 1) → L 2 (0, 1) defined by with D(A) = {u ∈ H 2 (0, 1) | u(0) = u(1) = 0}, which has also been examined in [20, pp. 115]. For u ∈ D(A) and v ∈ C ∞ (0, 1) we have x ∈ ( i−1 n , i n ), i + 1 − nx, x ∈ ( i n , i+1 n ), 0, else, for i ∈ {1, . . . , n−1} form a basis of U n . Evaluating (20) on these basis functions, the finite-element discretisation matrices A n of A are given by With the choice of η = 0.015, Figure 2 shows the eigenvalues of A n for n = 40 (red) and the sets B δj (W ((A n − s j ) −1 )) −1 + s j (blue) for a number of shifts s 1 , . . . , s m where δ j = 2.1 (A n − s j ) −1 2 ε and ε = 0.9 2 maxj=1,...,m (A−sj ) −1 . The shifts are located at a certain distance to the expected pseudospectrum so as to obtain a relatively small superset thereof. The black line corresponds to the boundary of σ ε (A n ) computed by eigtool, see [3]. This demonstrates the result of Theorem 3.5 which actually yields an enclosure for the pseudospectrum of the operator A while the black line only shows the boundary of the pseudospectrum of the approximation matrix A n .