A note on the residual type a posteriori error estimates for finite element eigenpairs of nonsymmetric elliptic eigenvalue problems

Applied Numerical Mathematics 82 (2014) 51–67

Contents lists available at ScienceDirect

Applied Numerical Mathematics www.elsevier.com/locate/apnum

A note on the residual type a posteriori error estimates for finite element eigenpairs of nonsymmetric elliptic eigenvalue problems ✩ Yidu Yang ∗ , Lingling Sun, Hai Bi, Hao Li School of Mathematics and Computer Science, Guizhou Normal University, Guiyang 550001, China

a r t i c l e

i n f o

Article history: Received 28 December 2012 Received in revised form 31 October 2013 Accepted 21 February 2014 Available online 1 April 2014 Keywords: Nonsymmetric eigenvalue problems Finite elements A posteriori error estimates

a b s t r a c t In this paper we study the residual type a posteriori error estimates for general elliptic (not necessarily symmetric) eigenvalue problems. We present estimates for approximations of semisimple eigenvalues and associated eigenvectors. In particular, we obtain the following new results: 1) An error representation formula which we use to reduce the analysis of the eigenvalue problem to the analysis of the associated source problem; 2) A local lower bound for the error of an approximate finite element eigenfunction in a neighborhood of a given mesh element T . © 2014 IMACS. Published by Elsevier B.V. All rights reserved.

1. Introduction The ideas to use a posteriori error estimates and adaptive finite element algorithms were first introduced by Babuska and Rheinboldt in 1978 (see [4]). Since 1981, a posteriori error estimates for finite element methods have been developed greatly, and many kinds of a posteriori error estimates, e.g., of residual type or recovery type, have been proposed. The first a posteriori error estimate proposed by Babuska and Rheinboldt [4] was of residual type. This kind of estimate has been studied deeply and applied to many problems widely (cf. [1,5,6,8,12,15,23,28,29,32] and the references cited therein). A posteriori error estimates of residual type have also been applied to conforming finite element methods (see [2,6,17,19, 21,22,24,26,31,33]), nonconforming finite element methods and mixed finite element methods (see, for example, [16,25,27]) for eigenvalue problems. In Heuveline and Rannacher [21,22], second-order nonsymmetric elliptic eigenvalue problems have been treated by employing duality techniques from optimal control theory. In Gedicke and Carstensen [18], for the case of linear P 1 finite elements and globally constant coefficients, the error estimates of the residual and averaging error estimators are refined. In comparison, the studies in [18,21,22] only considered simple eigenvalues (algebraic multiplicity q = 1) and did not give a local lower bound for the error of a finite element eigenfunction. Although Remark 4 in [21] pointed out that Proposition 2 in [21] is also (asymptotically) true for semisimple eigenvalues (ascent α = 1), it has not been proved that  v h∗  in the proposition is uniformly bounded from above with respect to h. Thus, it seems that a complete theoretical analysis is incomplete. We know that, due to symmetry considerations, eigenvalues of nonsymmetric problems are in general not simple and may even fail to be semisimple. In addition, in order to implement an adaptive finite element calculation, it

✩

*

Supported by the National Natural Science Foundation of China (Grant No. 11161012). Corresponding author. E-mail addresses: [email protected] (Y. Yang), [email protected] (L. Sun), [email protected] (H. Bi), [email protected] (H. Li).

http://dx.doi.org/10.1016/j.apnum.2014.02.015 0168-9274/© 2014 IMACS. Published by Elsevier B.V. All rights reserved.

52

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

is necessary to give a local lower bound estimate of a finite element eigenfunction since it tells us the estimates of local element contributions are efficient. Therefore, based on the above discussion, this paper presents results on a posterior error estimates of finite element approximations of a semisimple eigenvalue of finite multiplicity, and obtains the following new results. (1) We establish an error representation formula which we use to reduce the analysis of the eigenvalue problem to the analysis of the associated source problem. (2) We prove that Theorem 3.2 is true for semisimple eigenvalues with multiplicity larger than one, and we explain in Remark 3.2 that the assumption that the ascent of λ is equal to 1 cannot guarantee that u h∗ in Theorem 3.2 is an eigenfunction corresponding to λh∗ . (3) Using the method from [33,34], we obtain a local lower bound for the error of a finite element eigenfunction in a neighborhood of a given mesh element T . 2. Preliminaries Let Ω ⊂ R d (d  2) be a polyhedral bounded domain with boundary ∂Ω = Γ . For S ⊂ Ω , we denote by H 1 ( S ) the usual Sobolev space endowed with the norm

1/2   v 1, S =  v 2L 2 ( S ) + ∇ v 2L 2 ( S ) . For convenience, we denote  · 1 =  · 1,Ω ,  · 0, S =  ·  L 2 ( S ) and  · 0 =  · 0,Ω . Let V = { v: v ∈ H 1 (Ω), v |Γ = 0} where v |Γ = 0 is in the sense of trace. Consider the eigenvalue problem for the following non-selfadjoint elliptic operator.





Lu ≡ −∇ · a(x)∇ u + b(x) · ∇ u + c (x)u = λu , u = 0,

in Ω,

on Γ.

(2.1) (2.2)

This problem has important physical background such as convection–diffusion problems in fluid mechanics, environmental problems and others. Thus, finite element methods for solving this problem become an important topic which has attracted many attentions from mathematics and physics community. Let



a (u , v ) =



a∇ u · ∇ v dx + Ω

 (u , v ) =

(b · ∇ u v + cuv ) dx, Ω

uv dx. Ω

Assume that a = a(x) ∈ L ∞ (Ω) is a bounded and measurable real function on Ω and has a positive lower bound, b = b(x) ∈ W 1,∞ (Ω)d and c = c (x) ∈ L ∞ (Ω) are bounded and measurable real or complex functions on Ω , and there exist positive constants ca and C a such that

Re a( v , v )  ca  v 21 ,

∀v ∈ V ,   a(u , v )  C a u 1  v 1 , ∀u , v ∈ V .

(2.3) (2.4)

The variational form associated with (2.1), (2.2) is given by: Find λ ∈ C (complex plane), u ∈ V , u 0 = 1, satisfying

a(u , v ) = λ(u , v ),

∀v ∈ V .

(2.5)

Let πh = { T } be a d-simplicial mesh of Ω . For each element T ∈ πh , let h T be the diameter of T and ρ T be the diameter of the biggest ball contained in T , h = max{h T : T ∈ πh }. We further assume that πh is a shape-regular mesh (see [13]): there exists a constant γ ∗ such that

hT

ρT

 γ ∗ ∀T ∈



πh .

h

Let εh denote the set of all (d − 1)-faces in πh . There holds the split εh = εh,Ω ∪ εh,Γ where εh,Ω and εh,Γ refer to interior faces and faces on the boundary Γ , respectively. For each T ∈ πh , we denote by ε ( T ) the set of its edges. For each E ∈ εh , let h E denote the diameter of E. For any piecewise continuous function ϕ and any E ∈ εh,Ω , we use the symbol [[ϕ ]] E to denote the jump of ϕ across E in an arbitrary but fixed direction ν E orthogonal to E. For any T ∈ πh and E ∈ εh , we set

ωT =



ε ( T )∩ε ( T ) =∅

T ,

ωE =



T .

E ⊂∂ T

Let V h ⊂ V be finite element spaces consisting of continuous piecewise polynomials over πh of degree less than or equal to k. The finite element approximation of (2.1) is given by: Find λh ∈ C , u h ∈ V h , u h 0 = 1, such that

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

a(u h , v ) = λh (u h , v ),

∀v ∈ V h .

53

(2.6)

It is shown in Section 8 in [3] that (2.3), (2.4) imply that there are linear bounded operators A : L 2 (Ω) → V and A h : L 2 (Ω) → V h satisfying

a( A f , v ) = ( f , v ),

∀v ∈ V ,

a( A h f , v ) = ( f , v ),

(2.7)

∀v ∈ V h .

(2.8)

Also from [3] we know that A : L 2 (Ω) → L 2 (Ω) is a compact operator, and (2.5) and (2.6) have the following equivalent operator form (2.9) and (2.10), respectively.

Au = λ−1 u ,

(2.9)

−1

(2.10)

A h u h = λh u h . Thanks to [3], we know the adjoint problem of (2.1), (2.2) is as follows:









L ∗ u ∗ ≡ −∇ · a(x)∇ u ∗ − ∇ · b(x)u ∗ + c (x)u ∗ = λ∗ u ∗ , ∗

u = 0,

in Ω,

on Γ.

(2.11) (2.12)

The corresponding variational form and discrete variational form of (2.11), (2.12) are given by: Find λ∗ ∈ C , u ∗ ∈ V , u ∗ 0 = 1, satisfying









a v , u ∗ = λ∗ v , u ∗ ,

∀v ∈ V .

(2.13)

Find λh∗ ∈ C , u h∗ ∈ V h , u h∗ 0 = 1, satisfying









a v , u h∗ = λh∗ v , u h∗ ,

∀v ∈ V h .

(2.14)

Note that the primal and dual eigenvalues are connected via λ = λ∗ and λh = λh∗ . For (2.11), (2.12), from [3] we know that (2.3), (2.4) imply that there are linear bounded operators A ∗ : L 2 (Ω) → V and A h∗ : L 2 (Ω) → V h satisfying









a v , A ∗ f = ( v , f ),

∀v ∈ V ,

(2.15)

a v , A h∗ f = ( v , f ),

∀v ∈ V h .

(2.16)

And (2.13) and (2.14) have the following equivalent operator form (2.17) and (2.18), respectively.

A ∗ u ∗ = λ∗−1 u ∗ , ∗ ∗

∗−1 ∗

A h u h = λh

uh .

(2.17) (2.18)

It can be proved that A ∗ is the adjoint operator of A in the sense of inner product (·, ·). In fact,

    ( Au , v ) = a Au , A ∗ v = u , A ∗ v , ∀u , v ∈ L 2 (Ω),     ( A h u , v ) = a A h u , A h∗ v = u , A h∗ v , ∀u , v ∈ L 2 (Ω). Let λ be an eigenvalue of (2.5). There exists a smallest integer α , called the ascent of λ, such that N ((λ−1 − A )α ) = N ((λ−1 − A )α +1 ), where N denotes the null space. The number q = dim N ((λ−1 − A )α ) is called the algebraic multiplicity of λ. The vectors in N ((μ − A )α ) are called generalized eigenvectors of A corresponding to λ. Similarly, the ascent, algebraic multiplicity and generalized eigenvectors of λh , λ∗ and λh∗ can be defined. In this paper, let λ be an eigenvalue of (2.5) with the algebraic multiplicity q and the ascent α  1. Let λh be the eigenvalue of (2.6) which converges to λ. Let M (λ) be the space spanned by all generalized eigenfunctions corresponding to λ of A, and let M h (λ) be the space spanned by all generalized eigenfunctions corresponding to all eigenvalues of A h that converge to λ. As for the adjoint problems (2.13) and (2.14), the definitions of M ∗ (λ∗ ) and M h∗ (λ∗ ) are analogous to M (λ) and M h (λ). Notice that λ∗ (λ∗ = λ) and λ have the same ascent α . When α = 1, i.e., λ is a semisimple eigenvalue, M (λ) and M ∗ (λ∗ ) are all eigenfunction spaces. In this paper, the letter C (with or without subscripts) denotes a positive constant independent of h, which may not be the same constant in different places. For simplicity, we use the symbol a  b to mean that a  Cb and the symbol R = O (ε ) to mean that | R |  C ε . Regarding the a priori error estimates of finite element approximations, (2.6) and (2.14), we refer to [3,7].

54

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

Lemma 2.1. Assume that M (λ) ⊂ H t1 (Ω), and M ∗ (λ∗ ) ⊂ H t2 (Ω). Then

 |λh − λ|  hτ1 +τ2 −2

sup

v ∈ M (λ), v 1 =1

 v t1

sup

v ∗ ∈ M ∗ (λ∗ ), v ∗ 1 =1

 ∗  α1 v  . t 2

(2.19)

Let u h be a unit vector satisfying (λh−1 − A h )l1 u h = 0 for some positive integer l1  α . Then for any integer l2 with l1  l2  α , there

exists a vector u such that (λ−1 − A )l2 u = 0 and

 uh − u 1  hτ1 −1

sup

v ∈ M (λ), v 1 =1

 v t1

(l2 −l1 +1)/α

.

(2.20)

Let u h∗ be a unit vector satisfying (λh∗−1 − A h∗ )l1 u h∗ = 0 for some positive integer l1  α . Then for any integer l2 with l1  l2  α , there

is a vector u ∗ such that (λ∗−1 − A ∗ )l2 u ∗ = 0 and

  ∗   u − u ∗   h τ 2 −1 h 1

sup

 ∗  (l2 −l1 +1)/α v  .

v ∗ ∈ M ∗ (λ∗ ), v ∗ 1 =1

t2

(2.21)

In the above, τ1 = min(k + 1, t 1 ), τ2 = min(k + 1, t 2 ). Proof. See [3].

2

The following regularity estimates can be found in [20]: For any f ∈ L 2 (Ω), A f ∈ H 01 (Ω) ∩ H 1+γ1 (Ω) and A ∗ f ∈ H 01 (Ω) ∩ H (Ω) satisfying 1+γ2

 A f 1+γ1  C Ω  f 0 ,  ∗  A f   C Ω  f 0 1 +γ

(2.22) (2.23)

2

for some γ1 , γ2 ∈ (0, 1] depending on Ω and the coefficients of L. Thus, by Aubin–Nitsche technique we can get the L 2 -norm error estimate of finite element eigenfunction. 3. A posteriori error identity relation of conforming finite elements for nonsymmetric eigenvalue problems In order to estimate the error of approximate eigenpairs (λh , u h ) we consider the associated source and approximate source problems: Find w ∈ V such that

a( w , v ) = λh (u h , v ),

∀v ∈ V .

(3.1)

Find w h ∈ V h such that

a( w h , v ) = λh (u h , v ),

∀v ∈ V h .

(3.2)

Note that the right hand sides of (3.1) and (3.2) are λh u h . According to the definitions of A and A h , we know that w = λh Au h and w h = λh A h u h = u h are the solutions of (3.1) and (3.2), respectively. Let

e = λh Au h − λh A h u h . When b = 0 and c (x) is a real function, (2.1), (2.2) is a symmetric elliptic eigenvalue problem, and in this case [14,35] discussed the relation between the error of conforming finite element eigenfunction u h and the error of conforming finite element solution of the above boundary value problem, and proved that the a posteriori error indicator of the finite element solution of the boundary value problem is the a posteriori error indicator of the finite element eigenpair. Based on [14,35], we prove the following Theorem 3.1 for (2.1), (2.2). Let D ⊂ Ω . Theorem 3.1. Assume that (λh , u h ) is an eigenpair of (2.6) with u h 0 = 1. Then there exists an eigenpair (λ, u ) of (2.5) such that

    u − uh 1, D = ( A − A h )(λh uh )1, D + O λu − λh uh 0 .

(3.3)

Assume that (λh∗ , u h∗ ) is an eigenpair of (2.14) with u h∗ 0 = 1. Then there exists an eigenpair (λ∗ , u ∗ ) of (2.13) such that

 ∗  u − u ∗ 

h 1, D

      =  A ∗ − A h∗ λh∗ uh∗ 1, D + O λ∗ u ∗ − λh∗ uh∗ 0 .

(3.4)

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

55

Proof. According to (2.9) and (2.10), we have

u − u h = λ Au − λh A h u h

= λ Au − λh Auh + λh Auh − λh A h uh = A (λu − λh uh ) + ( A − A h )(λh uh ).

(3.5)

Denote





R = u − u h 1, D − ( A − A h )(λh u h )1, D . From (3.5) we obtain

u − u h − ( A − A h )(λh u h ) = A (λu − λh u h ). Then we have

    | R | = u − uh 1, D − ( A − A h )(λh uh )1, D      A (λu − λh uh )1, D  λu − λh uh 0 , from which we know that (3.3) is true. Using the similar argument, we can prove that (3.4) is valid.

2

Remark 3.1. Theorem 3.1 attributes the error estimate of finite element eigenfunction u h to the error of finite element solution of the corresponding source problem (with right hand side f = λh u h ) in the sense of local energy norm. And the finite element solution of the source problem is just u h . So the a posteriori error indicator of finite element for the source problem is also the a posteriori error indicator of the corresponding eigenvalue problem. We denote E h∗ as the spectral projection associated with A h∗ and the eigenvalues of A h∗ which converge to λ∗ , and Ran( E h∗ )

as the image space of E h∗ . Let l be the ascent of λh∗ , and M h∗ (λh∗ ) = N ((λh∗−1 − A h∗ )l ) be the generalized eigenfunction space corresponding to the eigenvalue λh∗ . Then





 

Ran E h∗ = M h∗ λ∗ ,

 

 

M h∗ λh∗ ⊂ M h∗ λ∗ .

Theorem 3.2. Assume that (λh , u h ) and (λ, u ) satisfy Theorem 3.1, and λ is a semisimple eigenvalue. Set λh∗ = λh , λ∗ = λ. Then there exists u h∗ ∈ M h∗ (λh∗ ) with u h∗ 0 = 1 such that |(u h , u h∗ )| has a positive lower bound uniformly with respect to h, and there exists u ∗ ∈ M ∗ (λ∗ ) such that

 ∗   u − λ∗ A ∗ u ∗ 

      =  A ∗ − A h∗ λh∗ uh∗ 1, D + O λ∗ u ∗ − λh∗ uh∗ 0 ,      1 ( A − A h )(λh uh )2 +  A ∗ − A ∗ λ∗ u ∗ 2 |λh − λ|  h h h 1 ∗ 1 |(uh , uh )|  2  + λu − λh uh 20 + λ∗ u ∗ − λh∗ uh∗ 0 . h

h h 1, D

(3.6)

(3.7)

Proof. Let u h− be the orthogonal projection from u h to N ((λh∗−1 − A h∗ )l ) in the sense of inner product (·, ·), where l is the ascent of λh∗ . We choose u h∗ = u h− /u h− 0 . On the one hand, when λh and λh are two eigenvalues of A h and λh = λh , from [10] we know that the generalized eigenvector of A h corresponding to λh and any generalized eigenvector of A h∗ corresponding to λh ∗ are perpendicular, and

N ((λh∗−1 − A h∗ )l ) ⊂ Ran( E h∗ ) holds. So the orthogonal projection u h− from u h to N ((λh∗−1 − A h∗ )l ) is the orthogonal projection from u h to Ran( E h∗ ). On the other hand, according to p. 689 in [3], there exists w h ∈ N ((λ∗−1 − A ∗ )α ) (where α = 1) such that (u h , w h ) = 1 and w h is uniformly bounded in h. Let u h = E h∗ w h /(u h , E h∗ w h ), then u h ∈ Ran( E h∗ ),





u h , u h = 1

and (u h , E h∗ w h ) → (u h , w h ) = 1. So u h 0  C  E h∗ w h 0  C 0 is uniformly bounded in h. Noticing that u h− is the orthogonal projection from u h to Ran( E h∗ ) and



uh ,

u h−

(uh , uh− )

we know that



= 1,

56

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

   uh−         (u , u − )   u h 0  C 0 h

0

h

are uniformly bounded with respect to h, and thus



     u h−  uh , u ∗  =  uh ,  1 , h  − uh 0  C 0 i.e., |(u h , u h∗ )| has a positive lower bound uniformly with respect to h. Select u ∗ ∈ M ∗ (λ∗ ) satisfying Lemma 2.1. Since the ascent of λ is equal to 1, M ∗ (λ∗ ) is the eigenfunction space. So, from (2.17) we deduce

u ∗ − λh∗ A h∗ u h∗ = λ∗ A ∗ u ∗ − λh∗ A h∗ u h∗

= λ∗ A ∗ u ∗ − λh∗ A ∗ uh∗ + λh∗ A ∗ uh∗ − λh∗ A h∗ uh∗      = A ∗ λ∗ u ∗ − λh∗ uh∗ + A ∗ − A h∗ λh∗ uh∗ .

(3.8)

Let

     R = u ∗ − λh∗ A h∗ u h∗ 1, D −  A ∗ − A h∗ λh∗ u h∗ 1, D . According to (3.8), we conclude

     | R | =  A ∗ λ∗ u ∗ − λh∗ uh∗ 1, D  λ∗ u ∗ − λh∗ uh∗ 0 , which implies that (3.6) is valid. From (2.16) we get

















a u h , λh∗ A h∗ u h∗ = u h , λh∗ u h∗ = λh∗ u h , u h∗ = λh u h , u h∗ , and from (2.6) and (2.16) we obtain





u h , λh∗ A h∗ u h∗ =

    1  1  a u h , λh∗ A h∗ u h∗ = u h , λh∗ u h∗ = u h , u h∗ .

λh

λh

Thus, by calculation we derive







a u h − u , λh∗ A h∗ u h∗ − u ∗ − λ u h − u , λh∗ A h∗ u h∗ − u ∗



        = a uh , λh∗ A h∗ uh∗ − a uh , u ∗ − a u , λh∗ A h∗ uh∗ + a u , u ∗         − λ uh , λh∗ A h∗ uh∗ − uh , u ∗ − u , λh∗ A h∗ uh∗ + u , u ∗         = λh uh , uh∗ − λ uh , u ∗ − λ u , λh∗ A h∗ uh∗ + λ u , u ∗         − λ uh , uh∗ − uh , u ∗ − u , λh∗ A h∗ uh∗ + u , u ∗   = (λh − λ) uh , uh∗ ,

i.e.,

λh − λ =

1

(uh , uh∗ )







a u h − u , λh∗ A h∗ u h∗ − u ∗ − λ u h − u , λh∗ A h∗ u h∗ − u ∗

 ,

(3.9)

hence

|λh − λ|  

Ca

|(uh , uh∗ )| 1

|(uh , uh∗ )|

      uh − u 1 λh∗ A h∗ uh∗ − u ∗ 1 + |λ|uh − u 0 λh∗ A h∗ uh∗ − u ∗ 0   uh − u 1 λh∗ A h∗ uh∗ − u ∗ 1 .

(3.10)

It is easy to know that when D = Ω , (3.3) and (3.6) are still valid. Combining (3.10), (3.3) and (3.6) with D = Ω , we obtain (3.7). The proof is completed. 2

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

57

Remark 3.2. The assumption that the ascent of λ is equal to 1 can guarantee that the ascent of λ∗ (λ∗ = λ) is 1 but cannot guarantee that the ascent of λh∗ is 1 (see Example 5.5 in [10]). Therefore, in general, it fails to ensure that u h∗ is an eigenfunction corresponding to λh∗ such that (u h , u h∗ ) has a positive lower bound uniformly with respect to h, but it is certain that u h∗ is a generalized eigenfunction, i.e., u h∗ ∈ M h∗ (λh∗ ). For instance, consider



Ah =

1 0 h 1





A=

,

1 0 0 1



.

It is obvious that



  0 0   = h → 0 (h → 0).   Ah − A  =  h 0 

By calculation we know that λ = 1 is the eigenvalue of A with the algebraic multiplicity 2 and the ascent 1; λh = 1 is the eigenvalue of A h with the algebraic multiplicity 2 and the ascent 2 and the corresponding eigenvector e 2 = (0, 1) T . Let

A h∗ =



1 h 0 1



.

A simple calculation shows that λh∗ = 1 is the eigenvalue of A h∗ with the algebraic multiplicity 2 and the ascent 2 and the corresponding eigenvector e ∗1 = (1, 0) T . We can see that e ∗2 = (0, 1) T is not an eigenvector of A h∗ but



∗



 ∗ 2 ∗

A h − λh

e2 =

0 h 0 0

2

e ∗2 = 0.

Hence e ∗2 = (0, 1) T is a generalized eigenvector of A h∗ corresponding to λh∗ . It can be seen that (e 2 , e ∗1 ) = 0 does not have a positive lower bound while (e 2 , e ∗2 ) = 1 has a positive lower bound. Remark 3.3. Theorem 3.2 is applicable to the case that the ascent of λh is greater than or equal to 1. When the ascent of λh is equal to 1, M h∗ (λh∗ ) is spanned by eigenfunctions of A h∗ which are corresponding to λh∗ and u h∗ is an eigenfunction of A h∗ corresponding to λh∗ , thus we have λh∗ A h∗ u h∗ = u h∗ . 4. The residual type a posteriori error estimate of conforming finite elements for nonsymmetric eigenvalue problems 4.1. The global bounds for the error Consider the associated source and approximate source problems: Find w ∈ V such that

a( w , v ) = ( f , v ),

∀v ∈ V .

(4.1)

Find w h ∈ V h such that

a( w h , v ) = ( f , v ),

∀v ∈ V h .

(4.2)

Let

[[a  w h ]] E · ν E = −a  w h+ · ν + − a  w h− · ν − , where E is the common side of elements T + and T − with unit outward normals ωE = T + ∪ T − . Define the element residual R T ( f , w h ) and the jump residual J E ( w h ):

ν + and ν − , respectively, and ν E = ν − . Let

R T ( f , wh) = f − L wh

= f +  · (a  w h ) − b · ∇ w h − c w h in T ∈ πh ,  [[a  w h ]] E · ν E if E ∈ εh,Ω , J E (wh) = 0 if E ∈ εh,Γ .

(4.3) (4.4)

For T ∈ πh , G ⊂ Ω , define





ηG ( f , w h ) =

 T ∈πh , T ⊂G



2



ηT ( f , w h ) = h2T  R T ( f , w h )0,T +

E ∈εh,Ω , E ⊂∂ T



2

h E  J E ( w h )0, E

1/2 ,

(4.5)

1/2

η

2 T ( f , wh)

.

(4.6)

58

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

Consider the associated source problem and approximate source problems: Find w ∗ ∈ V such that





a v , w ∗ = ( f , v ),

∀v ∈ V .

(4.7)

∀v ∈ V h .

(4.8)

Find w h∗ ∈ V h such that





a v , w h∗ = ( f , v ),

Define the element residual R ∗T ( f , w h∗ ) and the jump residual J ∗E ( w h∗ ):













R ∗T f , w h∗ = f − L ∗ w h∗ = f +  · a(x)  w h∗ + ∇ · bw h∗ − c w h∗

  = f +  · a(x)  w h∗ + b · ∇ w h∗ + (∇ · b − c ) w h∗ in T ∈ πh ,    [[a  w h∗ ]] E · ν E if E ∈ εh,Ω , J ∗E w h∗ = 0 if E ∈ εh,Γ .

(4.9) (4.10)

For T ∈ πh , G ⊂ Ω , define













2 η∗T f , w h∗ = h2T  R ∗T f , w h∗ 0,T +

∗



∗







ηG f , w h =

∗



 ∗ 2

ηT f , w h





E ∈εh,Ω , E ⊂∂ T





2 h E  J ∗E w h∗ 0, E

1/2 ,

(4.11)

1/2

(4.12)

.

T ∈πh , T ⊂G

Denote

a S = min a(x), x∈ S

S ⊂ Ω;



1/2 αT = min h T a− , 1 , T ∈ πh ; T



1/2 αω E = min h E a− ωE , 1 ,

E ∈ εh .

T (·, ·) denotes a polyIn the following, we use  F to denote a polynomial approximation of a function F , for example, R nomial approximation of R T (·, ·) on an element T , and  J E (·) is a polynomial approximation of J E (·) on E. We define a norm





| v |2S =



a(x)∇ v · ∇ v + v v dx.

S

Denote | v |Ω = | v | for simplicity. It is obvious that the norm | · | is equivalent to the norm  · 1 . For the problems (4.1), (4.2) and (4.7), (4.8), we need the following results. Lemma 4.1. The following a posteriori error estimates for the source problems are valid:

 w − w h 1  ηΩ ( f , w h ),  ∗     w − w ∗   η∗ f , w ∗ , h 1

Ω

(4.13) (4.14)

h

for any T ∈ πh ,



1/2 ηT ( f , w h )  1 + c L ∞ (ωT ) + a− bL ∞ (ωT ) αT | w − w h |ωT T

+ 







T ⊂ω T





T ( f , w h ) , hT R T ( f , wh ) − R 0, T

(4.15)





1/2 bL ∞ (ωT ) αT  w ∗ − w h∗ ω η∗T f , w h∗  1 + c − ∇ · bL ∞ (ωT ) + a− T T

+



T ⊂ω T

Proof. See [28,33,34].











∗ f , w ∗  . h T  R ∗T f , w h∗ − R h 0, T T

2

Let (λh , u h ) and (λh∗ , u h∗ ) be an eigenpair of (2.6) and (2.14), respectively.

(4.16)

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

59

Selecting f = λh u h in (4.1) and (4.2), then w = A (λh u h ), w h = A h (λh u h ) and the a posterior error indicator of w h is

ηG (λh uh , w h ). Selecting f = λh∗ uh∗ in (4.7) and (4.8), then w ∗ = A ∗ (λh∗ uh∗ ), w h∗ = A h∗ (λh∗ uh∗ ) and the a posterior error indicator of w h∗ is ηG∗ (λh∗ u h∗ , w h∗ ).

For the eigenvalue problem, we use R T (λh u h , u h ) and J E (u h ) as the element residual and the jump residual for u h , respectively, and take η T (λh u h , u h ) and ηG (λh u h , u h ) as the a posteriori error indicators for u h . We use R ∗T (λh∗ u h∗ , u h∗ ) and J ∗E (u h∗ ) as the element residual and the jump residual for u h∗ , respectively, and take η∗T (λh∗ u h∗ , u h∗ ) and ηG∗ (λh∗ u h∗ , u h∗ ) as the a posteriori error indicators for u h∗ . Theorem 4.1. Assume that (λh , u h ) is an eigenpair of (2.6) with u h 0 = 1. Then there exists an eigenpair (λ, u ) of (2.5) such that

u − uh 1  ηΩ (λh uh , uh ) + λu − λh uh 0 ,   1/2  ηΩ (λh uh , uh )  1 + c L ∞ (Ω) + max a− bL ∞ (ωT ) αT u − uh 1 T +



T ∈πh

α



2 T R T (λh u h , u h ) −

 T (λh uh , uh )2 R

(4.17)

1/2

0, T

T ∈πh

+ λu − λh uh 0 .

(4.18)

Assume that the conditions of Theorem 3.2 hold. Then there exists u ∗ ∈ M ∗ (λ∗ ), satisfying

 ∗ ∗ ∗ ∗ ∗ 2 ∗ λh uh , λh A h uh |λh − λ|  ηΩ (λh uh , uh )2 + ηΩ 2  ∗ ∗ 2 + λu − λh uh 0 + λ u − λh∗ uh∗ 0 .

(4.19)

Proof. Selecting f = λh u h in (4.1) and (4.2), then w = λh Au h , w h = λh A h u h = u h . Therefore, from (4.13) we obtain

λh Auh − λh A h uh 1  ηΩ (λh uh , uh ).

(4.20)

Substituting (4.20) into (3.3) (taking D = Ω ), we get (4.17). According to (4.15) we derive





1/2 ηΩ (λh uh , uh )  1 + c L ∞ (Ω) + max a− bL ∞ (ωT ) αT T

+

T ∈πh



α



2 T R T (λh u h , u h ) −



|λh Auh − λh A h uh |Ω

 T (λh uh , uh )2 R

1/2 (4.21)

.

0, T

T ∈πh

Substituting (4.21) into (3.3) (taking D = Ω ), we obtain (4.18). In (4.7) and (4.8), selecting f = λh∗ u h∗ , then w ∗ = λh∗ A ∗ u h∗ , w h∗ = λh∗ A h∗ u h∗ . Therefore, from (4.14) we obtain

 ∗ ∗ ∗     λ A u − λ ∗ A ∗ u ∗   η ∗ λ ∗ u ∗ , λ∗ A ∗ u ∗ . Ω h h h h h h h h h h 1

(4.22)

Substituting (4.20) and (4.22) into (3.7), it is easy to get (4.19). The proof is completed.

2

Theorem 4.2. Let (λh∗ , u h∗ ) be an eigenpair of (2.14) with u h∗ 0 = 1. Then there exists an eigenpair (λ∗ , u ∗ ) of (2.13) such that

 ∗       u − u ∗   η ∗ λ ∗ u ∗ , u ∗ +  λ∗ u ∗ − λ ∗ u ∗  , Ω h h h h a h h 0   ∗ ∗ ∗   −1/2   ∗ λh uh , uh  1 + c − ∇ · bL ∞ (Ω) + max a T bL ∞ (ωT ) αT u ∗ − uh∗ 1 ηΩ +



T ∈πh

    2 α 2  R ∗ λ∗ u ∗ , u ∗ − R∗ λ∗ u ∗ , u ∗  T

T ∈πh

T

h h

h

T

h h

h

(4.23)

1/2

0, T

  + λ∗ u ∗ − λh∗ uh∗ 0 . Proof. Using the same argument of Theorem 4.1, we can prove (4.23) and (4.24).

(4.24)

2

Corollary 4.1. Assume that the coefficients a(x), b(x) and c (x) are properly smooth, (λh , u h ) is an eigenpair of (2.6) with u h 0 = 1, and u − u h 0 is an infinitesimal of higher order compared with u − u h 1 . Then ηΩ (λh u h , u h ) is a reliable and efficient a posteriori error indicator for u h . Assume that (λh , u h ) and (λh∗ , u h∗ ) satisfy Theorem 4.1, and the ascent of λ is equal to 1. Then ηΩ (λh u h , u h )2 + ∗ (λ∗ u ∗ , λ∗ A ∗ u ∗ )2 is a reliable a posteriori error indicator for λ . ηΩ h h h h h h

60

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

Proof. We use Theorem 4.1 to complete the proof. Note that the first items of the right hand sides of (4.17)–(4.19) are the dominant terms, then we conclude that Corollary 4.1 is true. 2 Remark 4.1. It is well known that (4.19) holds for simple eigenvalues (λh∗ A h∗ u h∗ = u h∗ ), see, for example, p. 4 in [9] or (6) in [21]. And in this paper we give and prove that (4.19) is valid for semisimple eigenvalues, which is a development of the existed work. Remark 4.2. From the proof in this paper, it is easy to know that Theorem 3.1, Theorem 3.2, (4.17), (4.19) and (4.23) are still valid for finite element spaces consisting of continuous piecewise polynomials over πh (d-simplices mesh, d-rectangles mesh, etc.). 4.2. Local lower bound for the error in a neighborhood of T Next, we will prove that the local error estimator η T (λh u h , u h ) provides a local lower bound for the error in a neighborhood of an element T . Given T ∈ πh . Let E be a (d − 1)-face of element T , we refer to p. 647 in [34] about the definition of ψ T , ψ E ,θ E and P E . Note that ψ T ∈ V , ψ T = 0 on Ω \ T , and ψ E ,θ ∈ V , ψ E ,θ = 0 on Ω \ ω E . We denote by Pk the set of complex-valued polynomials of degree at most k. The following Lemma 4.2 and Lemma 4.3 come from [34]. Lemma 4.2. The following estimates hold for all v ∈ Pk and T ∈ πh :

 v 20,T  ( v , ψT v )T ,  v ψT 0,T   v 0,T , −1/2 −1  v 0,T . | v ψT |T  min h T a T , 1 1/ 2 1 For E ∈ εh , set θ E = min{aω E h− E , 1}. Then the following estimates hold for all σ ∈ Pk | E and E ∈ εh :

σ 20, E  (σ , ψ E ,θ E P E σ ) E ,

1/4 −1/2 1/2 σ 0, E , ψ E ,θ E P E σ 0,ω E  aω E min h E aω E , 1 1/4 −1/2 −1/2 σ 0, E . |ψ E ,θ E P E σ |ω E  aω E min h E aω E , 1 Proof. For any x ∈ S, a simple calculation shows that

a(x)  max a( y )  C min a( y ) = Ca S . y∈ S

y∈ S

It is easy to know that Lemma 3.3 in [34] is also valid for complex-valued polynomials. Thus, combining the above inequality and Lemma 3.3 in [34], we obtain Lemma 4.2. 2 Lemma 4.3. The following trace inequality holds for all T ∈ πh , E ⊂ ∂ T , and v ∈ H 1 ( T ): −1/2

 v 0, E  Ch T

Proof. See Lemma 3.1 in [34]. For S ⊂ Ω , let



a S (u , v ) =

1/2

1/2

 v 0,T +  v 0,T ∇ v 0,T .

(4.25)

2





a(x)∇ u · ∇ v + b · ∇ uv + c (x)u v dx.

S

It is obvious that aΩ (u , v ) = a(u , v ). Lemma 4.4. The following inequality holds

  a S ( v , w )  | v | S | w | S 1 + c L

∞(S)

−1/2

+ aS

bL ∞ ( S ) | v | S  w 0, S ,

∀v , w ∈ V .

(4.26)

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

61

Proof. First we have

        a S ( v , w ) =  a(x)∇ v · ∇ w + c (x) v w + b · ∇ v w dx   S

          a(x)∇ v · ∇ w dx + c L ∞ ( S )  v 0, S  w 0, S +  b · ∇ v w dx S



S

 | v | S | w | S 1 + c L ∞ ( S ) + ∇ v 0, S  w 0, S bL ∞ ( S ) . Then it follows from a calculation that 1 ∇ v 20, S = a− S



1 a S |∇ v |2 dx  a− S

S



1 2 a(x)|∇ v |2 dx  a− S | v | S .

S

Combining the above two inequalities, we obtain (4.26).

2

Theorem 4.3. Assume that (λh , u h ) is an eigenpair of (2.6) with u h 0 = 1. Then there exists an eigenpair (λ, u ) of (2.5) such that for any T ∈ πh ,





1/2 ηT (λh uh , uh )  1 + c L ∞ (ωT ) + a− ωT b L ∞ (ωT ) α T |u − u h |ωT

  T (λh uh , uh ) + αT λu − λh uh 0,ωT + αT  R T (λh uh , uh ) − R 0,ω T  1/2    + αω E J E (uh ) − J E (uh ) 0, E .

(4.27)

Proof. Integration by parts elementwise yields that for v ∈ V

a (u − u h , v ) =



λu + ∇ · (a∇ uh ) − b · ∇ uh − cuh , v

T ∈πh

=



R T (λh u h , u h ), v



T ∈πh



+

T

+

 

 T

J E (u h ), v

E ∈εh,Ω

+

 

a[[∇ u h · ν E ]] E , v

E ∈εh,Ω



 E

E

(λu − λh uh , v )T .

(4.28)

T ∈πh

Given T ∈ πh , we set

T (λh uh , uh ). w T = ψT R

(4.29)

Taking v = w T in (4.28), we get



T (λh uh , uh ), w T a (u − u h , w T ) = R 

 T

T (λh uh , uh ) − R T (λh uh , uh ), w T − R

 T

+ (λu − λh uh , w T )T .

(4.30)

From (4.26) and Lemma 4.2 we deduce

    a(u − uh , w T ) = a T (u − uh , w T )  

−1/2  |u − uh |T | w T |T 1 + c L ∞ (T ) + a T bL ∞ (T )  w T 0,T   −1/2 −1  C |u − uh |T 1 + c L ∞ (T ) min h T a T , 1 

 −1/2 + a T bL ∞ (T )  R T (λh u h , u h )0, T . By Lemma 4.2 we have

 2  R T (λh u h , u h )

0, T

  T (λh uh , uh ), w T ,  R T

thus, from (4.29)–(4.31) we get the lower bound

(4.31)

62

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67



 , 1  R T (λh u h , u h )0, T   −1/2 −1/2

 |u − uh |T 1 + c L ∞ (T ) + a T bL ∞ (T ) min h T a T , 1  −1/2  T (λh uh , uh ) + min h T a T , 1  R T (λh uh , uh ) − R 0, T  + λu − λh uh 0,T .

−1/2

min h T a T

(4.32)

According to the inequality

   R T (λh uh , uh )

0, T

    T (λh uh , uh ) ,   R T (λh u h , u h )0, T +  R T (λh u h , u h ) − R 0, T

we know that (4.32) can be written as



 , 1  R T (λh uh , uh )0,T   −1/2 −1/2

 |u − uh |T 1 + c L ∞ (T ) + a T bL ∞ (T ) min h T a T , 1  −1/2  T (λh uh , uh ) + min h T a T , 1  R T (λh uh , uh ) − R 0, T  + λu − λh uh 0,T .

−1/2

min h T a T

(4.33)

Given E ∈ εh,Ω , we set

ω E = ψ E ,θ E P E  J E (u h ).

(4.34)

Choosing v = ω E in (4.28) we get

a (u − u h , ω E ) =



R T (λh u h , u h ), ω E

T ∈ω E

+



 T

  + J E (uh ), ω E E

(λu − λh uh , ω E )T .

(4.35)

T ∈ω E

From (4.26) and Lemma 4.2 we deduce

|a(u − uh , ω E )| = |aω E (u − uh , ω E )|  

−1/2  |u − uh |ω E |ω E |ω E 1 + c L ∞ (ω E ) + aω E bL ∞ (ω E ) ω E 0,ω E   1/4 −1/2 −1/2  |u − uh |ω E 1 + c L ∞ (ω E ) aω E min h E aω E , 1  −1/2 1/4 −1/2 1/2   + aω E bL ∞ (ω E )aω E min h E aω E , 1 J E (u h )0, E   −1/2 −1/2

 |u − uh |ω E 1 + c L ∞ (ω E ) + aω E bL ∞ (ω E ) min h E aω E , 1  1/4 −1/2 −1/2   J E (u h )0, E . × aω E min h E aω E , 1

(4.36)

From (4.32) and Lemma 4.2 we derive

 T ⊂ω E

R T (λh u h , u h ), ω E

 T



    R T (λh uh , uh ) T ⊂ω E

     R T (λh u h , u h ) T ⊂ω E

 0, T

ω E ω E

0, T

      + R T (λh u h , u h ) − R T (λh u h , u h ) 0, T ω E 0,ω E T ⊂ω E

  −1/2 −1/2

 |u − uh |ω E 1 + c L ∞ (ω E ) + aω E bL ∞ (ω E ) min h T aω E , 1  −1/2  T (λh uh , uh ) + min h T aω E , 1  R T (λh uh , uh ) − R 0,ω E   1/4 −1/2 −1/2   + λu − λh uh 0,ω E × aω E min h E aω E , 1 J E (u h )0, E . A simple calculation shows that

(4.37)

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

63

    −1/2  (λu − λh uh , ω E )T   λu − λh uh 0,ω E min h E aω E , 1  T ∈ω E

 1/4 −1/2 −1/2   × aω E min h E aω E , 1 J E (u h )0, E .

(4.38)

Substituting (4.36), (4.37) and (4.38) into (4.35), we obtain

    J E (uh ), ω E  E   −1/2 −1/2

 |u − uh |ω E 1 + c L ∞ (ω E ) + aω E bL ∞ (ω E ) min h E aω E , 1  −1/2  T (λh uh , uh ) + min h T aω E , 1  R T (λh uh , uh ) − R 0,ω E   1/4 −1/2 −1/2   + λu − λh uh 0,ω E × aω E min h E aω E , 1 J E (u h )0, E .

(4.39)

From Lemmas 4.2 and 4.3 we have

 2    J E (u h )0, E  C  J E (u h ), ω E E     = C J E (uh ), ω E E + C  J E (u h ) − J E (u h ), ω E E      −1/2 1/2  J E (uh ), ω E E +  J E (u h ) − J E (u h )0, E min h E aω E , 1  1/4 −1/2 −1/2  × aω min h E aω , 1 J E (u h ) . E

(4.40)

0, E

E

Thus, from (4.39) and (4.40) we get the lower bound



1/2   J E (u h )0, E  −1/2 −1/2

 |u − uh |ω E 1 + c L ∞ (ω E ) + aω E bL ∞ (ω E ) min h E aω E , 1  −1/2 1/2   J E (u h ) − J E (u h )0, E + min h E aω E , 1  −1/2  T (λh uh , uh ) + min h T aω E , 1  R T (λh uh , uh ) − R 0,ω E  + λu − λh uh 0,ω E .

1/4



−1/2

aω E min h E aω E

,1 

From (4.33) and (4.41) we get (4.27). The proof is completed.

(4.41)

2

Theorem 4.4. Let (λh∗ , u h∗ ) be an eigenpair of (2.14) with u h∗ 0 = 1. Then there exists an eigenpair (λ∗ , u ∗ ) of (2.13) such that for any T ∈ πh ,



η∗T λh∗ uh∗ , uh∗



   −1/2  1 + c − ∇ · bL ∞ (ωT ) + aωT bL ∞ (ωT ) αT u ∗ − uh∗ ω T        ∗ λ∗ u ∗ , u ∗  + αT λ∗ u ∗ − λh∗ uh∗ 0,ω + αT  R ∗T λh∗ uh∗ , uh∗ − R h h h T 0,ω E T  ∗   1/2  ∗ ∗  ∗  + αω E J E uh − J E uh 0, E .

Proof. Using the proof method of Theorem 4.3, we can deduce (4.42).

(4.42)

2

5. Numerical experiments 5.1. Computational method for (λh∗ , u h∗ ) (see [36]) Assume that (λh , u h ) and (λh∗ , u h∗ ) satisfy Theorem 3.2, and let m be the algebraic multiplicity of λh and l be the ascent

of λh . It is obvious that λh∗ = λh . We provide the following method for computing u h∗ (see [36]).

Let u h− be the orthogonal projection of u h to N (( λ1∗ − A h∗ )l ), and u h∗ = u h− /u h− 0 . When u h ∈ N (( λ1∗ − A h∗ )l ) it is clear that h

h

1 ∗ l / N (( λ1∗ − A h∗ )l ), to find uh∗ , first we seek a basis {φi }m u h∗ = u h . When u h ∈ 1 of N (( λ∗ − A h ) ), and solve the following equations h

m  i =1

then let

αi (φi , φ j ) = (uh , φ j ),

h

j = 1, 2, . . . , m ,

(5.1)

64

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

Fig. 1. Ω = (0, 1) × (0, 1), the first eigenvalue.

u h− =

m  i =1

αi φi , 

(5.2)



u h∗ = u h− /u h− 0 .

(5.3)

Obviously, u h∗ satisfies



1

λ∗

l − A h∗ uh∗ = 0,

h

 ∗  u  = 1. h 0

1 ∗ l Thus, to find u h∗ which satisfying (5.2) and (5.3) in V h , the key is to seek a basis {φi }m 1 of N (( λ∗ − A h ) ). When l = 1, it

is actually to solve the following equations to obtain a basis in the solution space.

h

(1 )

uh ∈ V h ,



(1 ) 

a v , uh

 (1 )  − λh b v , uh = 0,

∀v ∈ V h .

(5.4)

(When λh is a simple eigenvalue, l = 1 and N (( λ1∗ − A h∗ )l ) is a one-dimensional space spanned by the eigenfunction u h∗ .) h

When l > 1, choose i = 1, 2, . . . , l to solve the following equations to obtain the basis. (i )

uh ∈ V h ,



(i ) 

a v , uh (0)

    − λh b v , uh(i ) = a v , uh(i −1) , (i −1)

∀v ∈ V h ,

(5.5) (i −1)

(i −1)

= 0 and uh is a vector which runs over the basis of the solution space of “a( v , u h ) − λh b( v , uh ) = (i −2) a( v , u h ), ∀ v ∈ V h ”. These equations have the same coefficient matrix. 1 ∗ l How to seek a basis {φi }m 1 of N (( λ∗ − A h ) ) efficiently is an important issue of linear algebra. Further algorithmic considwhere u h

h

erations are outside the scope of this paper. 5.2. Numerical examples In this section, we will report some numerical experiments by using linear finite elements on triangle meshes. We use MATLAB 2011b under the package of Chen (see [11]) to solve the following Example 1 and Example 2. In the numerical examples we use the following notations: λi ,h : the i-th finite element eigenvalue with smallest real part; λi : the i-th exact eigenvalues with smallest real part; ∗ (λ∗ u ∗ , λ∗ A ∗ u ∗ )2 is the a posteriori error indicator for λ . Φ(λi ,h ) = ηΩ (λi ,h u i ,h , u i ,h )2 + ηΩ i ,h i ,h i ,h i ,h h i ,h

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

65

Fig. 2. Ω = (0, 1) × (0, 1), the second eigenvalue.

Fig. 3. ((−1, 1) × (−1, 1)) \ ([0, 1] × [−1, 0]), the first eigenvalue.

Example 1. Consider the convection–diffusion equation

−u + b · ∇ u = λu ,

in Ω,

u = 0,

on ∂Ω,

(5.6)

where Ω = (0, 1) × (0, 1) and b = (b1 , b2 ) . The eigenvalues of (5.6) are T

λ j ,i =

b21 + b22 4

  + π 2 j2 + i2 ,

where j , i ∈ N+ . The adjoint eigenvalue problem has eigenvalues λ∗j ,i = λ j ,i (see [30]). We see that λ1 =

λ 2 = λ3 =

b21 +b22 4

b21 +b22 4

+ 2π 2 ,

+ 5π 2 , and the multiplicity of λ2 is equal to 2. For (5.6) with b = (1, 0) T , b = (3, 0) T and b = (10, 0) T ,

respectively, the numerical results on uniform isosceles right triangle meshes are shown in Figs. 1 and 2. From Figs. 1 and 2 we can see that the a posteriori error indicators presented in this paper are efficient and reliable, which is consistent with our theoretical analysis. Example 2. Consider the convection–diffusion equation (5.6) where Ω = ((−1, 1) × (−1, 1)) \ ([0, 1] × [0, −1]). The exact eigenvalues, which are unknown, are thereby replaced by approximate eigenvalues with high accuracy. For this problem

66

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

Fig. 4. ((−1, 1) × (−1, 1)) \ ([0, 1] × [−1, 0]), the eighth eigenvalue.

with b = (1, 0) T , b = (3, 0) T and b = (10, 0) T , respectively, the numerical results on adaptive grids are shown in Figs. 3 and 4. From Figs. 3 and 4 we can see that if the convection parameter b = (1, 0) T , b = (3, 0) T and b = (10, 0) T , the a posteriori error indicators can reflect the total trend of the error of discrete eigenvalues but it is not perfect especially when b = (10, 0) T . This is probably the consequence of the performance of linear algebra routine on a convection dominated problem. To treat such problems to get high accurate approximation much more careful design of the algorithm is needed. 6. Concluding remarks In this paper, we discuss the residual type a posteriori error estimates for finite element eigenpairs of nonsymmetric elliptic eigenvalue problems and obtain some new results, which develop the previous corresponding work. Based on the work in this paper, one may further consider the a posteriori error estimates of eigenspace. A natural idea is to choose a q q q q basis {u i ,h }1 and {u ∗i ,h }1 in M h (λ) and M h∗ (λ∗ ), respectively, and to use the a posteriori error estimates of {u i ,h }1 and {u ∗i ,h }1 ∗ ∗ to get the a posteriori error estimates of M h (λ) and M h (λ ), respectively. Regarding the theoretical analysis and numerical implementation of this approach, it needs a lot of work which will be very valuable. Acknowledgements The authors cordially thank the referees and the editor for their valuable comments and suggestions that led to the great improvement of this paper, especially the idea of analyzing the a posteriori error of eigenspace mentioned in Concluding remarks. And the authors also thank Jiayu Han and Yuanyuan Yu for their help. References [1] M. Ainsworth, J.D. Oden, A Posteriori Error Estimates for the Finite Element Analysis, John Wiley & Sons, New York, 2000. [2] M.G. Armentano, C. Padra, A posteriori error estimates for the Steklov eigenvalue problem, Appl. Numer. Math. 58 (2008) 593–601. [3] I. Babuska, J.E. Osborn, Eigenvalue problems, in: P.G. Ciarlet, J.L. Lions (Eds.), Finite Element Methods (Part 1), in: Handbook of Numerical Analysis, vol. 2, Elsevier Science Publishers, North-Holland, 1991, pp. 640–787. [4] I. Babuska, W.C. Rheinboldt, A posteriori error estimates for the finite element method, Int. J. Numer. Methods Eng. 12 (1978) 1597–1615. [5] R. Becker, R. Rannacher, An optimal control approach to a posteriori error estimation in finite element methods, Acta Numer. 10 (2001) 1–102. [6] R. Becker, S. Mao, Z. Shi, A convergent nonconforming adaptive finite element method with quasi-optimal complexity, SIAM J. Numer. Anal. 47 (6) (2010) 4639–4659. [7] J.H. Bramble, J.E. Osborn, Rate of convergence estimates for nonselfadjoint eigenvalue approximation, Math. Comput. 27 (1973) 525–545. [8] C. Carstensen, J. Hu, A. Orlando, Framework for the a posteriori error analysis of nonconforming finite elements, SIAM J. Numer. Anal. 45 (1) (2007) 68–82. [9] C. Carstensen, J. Gedicke, V. Mehrmann, A. Miedlar, A adaptive homotopy approach for non-selfadjoint eigenvalue problems, Numer. Math. 119 (3) (2011) 557–583. [10] F. Chatelin, Spectral Approximations of Linear Operators, Academic Press, New York, 1983. [11] L. Chen, iFEM: an innovative finite element methods package in MATLAB, 2008. [12] Z. Chen, R.H. Nochetto, Residual type a posteriori error estimates for elliptic obstacle problems, Numer. Math. 84 (2000) 527–548. [13] P.G. Ciarlet, Basic error estimates for elliptic problems, in: P.G. Ciarlet, J.L. Lions (Eds.), Finite Element Methods (Part 1), in: Handbook of Numerical Analysis, vol. 2, Elsevier Science Publishers, North-Holland, 1991, pp. 21–343. [14] X. Dai, J. Xu, A. Zhou, Convergence and optimal complexity of adaptive finite element eigenvalue computations, Numer. Math. 110 (2008) 313–355.

Y. Yang et al. / Applied Numerical Mathematics 82 (2014) 51–67

67

[15] S. Du, X. Xie, Residual-based a posteriori error estimates of non-conforming finite element method for elliptic problems with Dirac delta source terms, Sci. China Math. 51 (8) (2008) 1440–1460. [16] R.G. Durán, L. Gastaldi, C. Padra, A posteriori error estimators for mixed approximations of eigenvalue problems, Math. Models Methods Appl. Sci. 9 (1999) 1165–1178. [17] R.G. Durán, C. Padra, R. Rodriguez, A posteriori error estimators for the finite element approximations of eigenvalue problems, Math. Models Methods Appl. Sci. 13 (2003) 1219–1229. [18] J. Gedicke, C. Carstensen, A posteriori error estimators for non-symmetric eigenvalue problems, Preprint 659, DFG Research Center Matheon, Berlin, 2009. [19] I. Giani, G. Graham, A convergent adaptive method for elliptic eigenvalue problems, SIAM J. Numer. Anal. 47 (2009) 1067–1091. [20] P. Grisvard, Elliptic Problems in Nonsmooth Domains, Pitman, London, MA, 1985. [21] V. Heuveline, R. Rannacher, A posteriori error control for finite element approximations of elliptic eigenvalue problems, Adv. Comput. Math. 15 (2001) 107–138. [22] V. Heuveline, R. Rannacher, Adaptive FE eigenvalue approximation with application to hydrodynamic stability analysis, in: W. Fitzgibbon, et al. (Eds.), Proceedings of the International Conference on Advances in Numerical Mathematics, Moscow, Sept. 16–17, 2005, Institute of Numerical Mathematics RAS, Moscow, 2006, pp. 109–140. [23] J. Hu, Z. Shi, A new a posteriori error estimate for the Morley element, Numer. Math. 112 (2009) 24–40. [24] M.G. Larson, A posteriori and a priori error analysis for finite element approximations of self-adjoint elliptic eigenvalue problems, SIAM J. Numer. Anal. 38 (2) (2000) 608–625. [25] Y. Li, A posteriori error analysis of nonconforming methods for eigenvalue problem, J. Syst. Sci. Complex. 22 (2009) 495–502. [26] H. Li, Y. Yang, The adaptive finite element method based on multi-scale discretizations for eigenvalue problems, Comput. Math. Appl. 65 (2013) 1086–1102. [27] C. Lovadina, M. Lyly, R. Stenberg, A posteriori estimates for the Stokes eigenvalue problem, Numer. Methods Partial Differ. Equ. 25 (1) (2008) 244–257. [28] K. Mekchay, R.H. Nochetto, Convergence of adaptive finite element methods for general second order linear elliptic PDEs, SIAM J. Numer. Anal. 43 (2005) 1803–1827. [29] P. Morin, R.H. Nochetto, K. Siebert, Convergence of adaptive finite element methods, SIAM Rev. 44 (2002) 631–658. [30] A. Naga, Z. Zhang, Function value recovery and its application in eigenvalue problems, SIAM J. Numer. Anal. 50 (1) (2012) 272–286. [31] V.B. Noël, A posteriori error estimator for eigenvalue problem associated to the Schrödinger operator with magnetic field, Numer. Math. 99 (2004) 325–348. [32] Z. Shi, M. Wang, Finite Element Methods, Science Press, Beijing, 2010 (in Chinese). [33] R. Verfürth, A Review of a Posteriori Error Estimates and Adaptive Mesh-Refinement Techniques, Wiley-Teubner, New York, 1996. [34] R. Verfürth, A posteriori error estimates for convection–diffusion equations, Numer. Math. 80 (4) (1998) 641–663. [35] Y. Yang, Finite Element Methods for Eigenvalue Problems, Science Press, Beijing, 2012 (in Chinese). [36] Y. Yang, X. Fan, Generalized Rayleigh quotient and finite element two-grid discretization schemes, Sci. China Ser. A 52 (9) (2009) 1955–1972.