Multigrid methods based on shifted inverse iteration for the Maxwell eigenvalue problemProject supported by the National Natural Science Foundation of China (Grant No. 11161012).

Multigrid methods based on shifted inverse iteration for the Maxwell eigenvalue problemthanks: Project supported by the National Natural Science Foundation of China (Grant No. 11161012).

JIAYU HAN School of Mathematical Sciences, Guizhou Normal University, GuiYang, 550001, China (hanjiayu126@126.com).
Abstract

In this paper two types of multgrid methods, i.e., the Rayleigh quotient iteration and the inverse iteration with fixed shift, are developed for solving the Maxwell eigenvalue problem with discontinuous relative magnetic permeability and electric permittivity. With the aid of the mixed form of source problem associated with the eigenvalue problem, we prove the uniform convergence of the discrete solution operator to the solution operator in using discrete compactness of edge element space. Then we prove the asymptotically optimal error estimates for both multigrid methods. Numerical experiments confirm our theoretical analysis.

Key words. Maxwell eigenvalue problem, multigrid method, edge element, error analysis

AMS subject classifications. 65N25, 65N30

1 Introduction

The Maxwell eigenvalue problem is of basic importance in designing resonant structures for advanced waveguide. Up to now, the communities from numerical mathematics and computational electromagnetism have developed plenty of numerical methods for solving this problem (see, e.g., [1, 3, 4, 6, 7, 8, 10, 14, 15, 16, 23, 24, 25, 31, 32, 37]).

The difficulty of numerically solving the eigenvalue problem lies in imposing the divergence-free constraint. For this purpose, the nodal finite element methods utilize the filter, parameterized and mixed approaches to find the true eigenvalues [7, 14]. The researchers in electromagnetic field usually adopt edge finite elements due to the property of tangential continuity of electric field [4, 16, 20, 24]. Using edge finite element methods, when one only considers to compute the nonzero eigenvalue, the divergence-free constraint can be dropped from the weak form and satisfied naturally (see [37]). But this will introduce spurious zero eigenvalues. Since the eigenspace corresponding to zero is infinite-dimensional, usually, the finer the mesh is, the more the spurious eigenvalues there are. However, using this form there is no difficulty in computing eigenvalues on a very coarse mesh. So the work of [37] subtly applies this weak form to two grid method for the Maxwell eigenvalue problem. That is, they first solve a Maxwell eigenvalue problem on a coarse mesh and then solve a linear Maxwell equation on a fine mesh. Another approach is the mixed form of saddle point type in which a Lagrange multiplier is introduced to impose the divergence-free constraint (see [25, 1, 20]). A remarkable feature of the mixed form is its equivalence to the weak form in [37] for nonzero eigenvalues. The mixed form is well known as having a good property of no spurious eigenvalues being introduced. However, it is not an easy task to solve it on a fine mesh (see [1, 2]).

The multigrid methods for solving eigenvalue problems originated from the idea of two grid method proposed by [33]. Afterwards, this work was further developed by [22, 34, 35, 37]. Among them the more recent work [34] makes a relatively systematical research on multigrid methods based on shifted inverse iteration especially on its adaptive fashion.

Inspired by the above works, this paper is devoted to developing multigrid methods for solving the Maxwell eigenvalue problem. We first use the mixed form to solve the eigenvalue problem on a coarser mesh and then solve a series of Maxwell equations on finer and finer meshes without using the mixed form. Roughly speaking, we develop the two grid method in [37] into multigrid method where only nonzero eigenvalues are focused on. We prefer to use the mixed form instead of the one in [37] on a coarse mesh to capture the information of the true eigenvalues. One reason is that the mixed form can include the physical zero eigenvalues and rule out the spurious zero ones simultaneously. Using it, the physical zero eigenvalues can be captured on a very coarse mesh, which is necessary when the resonant cavity has disconnected boundaries(see, e.g., [24]). Another reason lies in that the mixed discretization of saddle point type is not difficult to solve on a coarse mesh.

In this paper, we study two types of multigrid methods based on shifted inverse iteration: Rayleigh quotient iteration and inverse iteration with fixed shift. The former is a well-known method for solving matrix eigenvalues but the corresponding coefficient matrix is nearly singular and difficult to solve to some extend. To overcome this difficulty the latter first performs the Rayleigh quotient iteration at previous few steps and then fixes the shift at the following steps as the estimated eigenvalue obtained by the former. Referring to the error analysis framework in [36] and using compactness property of edge element space, we first prove the uniform convergence of the discrete solution operator to the solution operator in and then the error estimates of eigenvalues and eigenfunctions for the mixed discretization; then we adopt the analysis tool in [34] that is different from the one in [37] and prove the asymptotically optimal error estimates for both multigrid methods. In addition, this paper is concerned about the theoretical analysis for the case of the discontinuous electric permittivity and magnetic permeability in complex matrix form, which has important applications for the resonant cavity being filled with different dielectric materials invariably. It is noticed that our multigrid methods and theoretical results are not only valid for the lowest order edge elements but also for high order ones. More importantly, based on the work of this paper, once given an a posteriori error indicator of eigenpair one can further develop the adaptive algorithms of shifted inverse iteration type for the problem. In the last section of this paper, we present several numerical examples to validate the efficiency of our methods in different cases.

Throughout this paper, we use the symbol to mean that , where denotes a positive constant independent of mesh parameters and iterative times and may not be the same in different places.

2 Preliminaries

Consider the Maxwell eigenvalue problem in electric field

\hb@xt@.01(2.1)
\hb@xt@.01(2.2)
\hb@xt@.01(2.3)

where is a bounded Lipschitz polyhedron domain in , is the tangential trace of . The coefficient is the electric permittivity, and is the magnetic permeability and piecewise smooth. In this paper, with being the angular frequency, is defined as the eigenvalue of this problem. We assume that are two positive definite Hermite matrices such that and there exist two positive numbers satifying

\hb@xt@.01(2.4)

2.1 Some weak forms

Let

equipped with the norm . Throughout this paper, and denote the norms in induced by the inner products and respectively. Define the divergence-free space:

The standard weak form of the Maxwell eigenvalue problem (LABEL:s1)-(LABEL:s3) is as follows: Find and such that

\hb@xt@.01(2.5)

where . Denote .

As the divergence-free space in (LABEL:2.5) is difficult to discretize, alternatively, we would like to solve the eigenvalue problem (LABEL:s1)-(LABEL:s3) in the larger space , that is: Find and such that

\hb@xt@.01(2.6)

Note that when (LABEL:s5) and (LABEL:2.5) are equivalent, since (LABEL:s5) implies the divergence-free condition holds for (e.g., see [37]).

According to (LABEL:2.4), we have

\hb@xt@.01(2.7)

In order to study the eigenvalue problem in we need the auxiliary bilinear form

which defines an equivalent norm in .

By Lax-Milgram Theorem we can define the solution operator as

\hb@xt@.01(2.8)

Then the eigenvalue problem (LABEL:2.5) has the operator form

The following mixed weak form of saddle point type can be found in [25, 1, 2, 6, 24]: Find with such that

\hb@xt@.01(2.9)
\hb@xt@.01(2.10)

where for any .

We introduce the corresponding mixed equation: Find and for such that

\hb@xt@.01(2.11)
\hb@xt@.01(2.12)

the following LBB condition can be verified by taking

This yields the existence and uniqueness of linear bounded operators and (see [9]). Due to Helmholtz decomposition , it is easy to see and , for any . Hence and share the same eigenpairs. More importantly, the operator is self-adjoint. In fact, ,

\hb@xt@.01(2.13)

Note that is compact as a operator from to and from to since compactly (see Corollary 4.3 in [20]).

2.2 Edge element discretizations and error estimates

We will consider the edge element approximations based on the weak forms (LABEL:2.5), (LABEL:s5) and (LABEL:3.1ss)-(LABEL:3.2ss). Let be a shaped-regular triangulation of composed of the elements . Here we restrict our attention to edge elements on tetrahedra because the argument for edge elements on hexahedra is the same. The () order edge element of the first family [30] generates the space

where is the polynomial space of degree less than or equal to on , is the homogeneous polynomial space of degree on , and . We also introduce the discrete divergence-free space

where is the standard Lagrangian finite element space vanishing on of total degree less than or equal to and

The standard finite element discretization of (LABEL:2.5) is stated as: Find and such that

\hb@xt@.01(2.14)

It is also equivalent to the following form for nonzero (see [37]): Find and such that

\hb@xt@.01(2.15)

In order to investigate the convergence of edge element discretization (LABEL:2.13), we have to study the convergence of edge element discretization for the associated Maxwell source problem. Then by Lax-Milgram Theorem we can define the solution operator as

\hb@xt@.01(2.16)

Then the eigenvalue problem (LABEL:2.13) has the operator form

Introduce the discrete form of (LABEL:3.1ss)-(LABEL:3.2ss): Find , such that

\hb@xt@.01(2.17)
\hb@xt@.01(2.18)

Introduce the corresponding operators: Find and for any

\hb@xt@.01(2.19)
\hb@xt@.01(2.20)

Due to discrete Helmholtz decomposition , it is easy to know and , for any . Hence and share the same eigenpairs.

Similar to (LABEL:s2.3)-(LABEL:s2.4), one can verify the corresponding LBB condition for the discrete mixed form (LABEL:s2.1)-(LABEL:s2.2). According to the theory of mixed finite elements (see [9]), we get for all ,

\hb@xt@.01(2.21)
\hb@xt@.01(2.22)

Similar to (LABEL:s2.6) we can prove is self-adjoint in the sense of . In fact, ,

The discrete compactness is a very interesting and important property in edge elements because it is intimately related to the property of the collective compactness. Kikuchi [26] first successfully applied this property to numerical analysis of electromagnetic problems, and more recently it was further developed by [5, 6, 12, 27, 29] and so on. The following lemma, which states the discrete compactness of into , is a direct citation of Theorem 4.9 in [20].

Lemma 2.1

(Discrete compactness property) Any sequence with that is uniformly bounded in contains a subsequence that converges strongly in .

In the remainder of this subsection, we will prove the error estimates for the discrete forms (LABEL:3.1s)-(LABEL:3.2s), (LABEL:2.13) or (LABEL:2.14) with . The authors in [36] have built a general analysis framework for the a priori error estimates of mixed form (see Theorem 2.2 and Lemma 2.3 therein). Although we cannot directly apply their theoretical results to the mixed discretization (LABEL:3.1s)-(LABEL:3.2s), we can use its proof idea to derive the following Lemma 2.2 and Theorem 2.3. The following uniform convergence provides us with the possibility to use the spectral approximation theory in [11].

Lemma 2.2

There holds the uniform convergence

Proof. Since and are dense in and , respectively, we deduce from (LABEL:s2.18) for any

That is, converges to pointwisely. Since are linear bounded uniformly with respect to , is a bounded set in where is the unit ball in . From compactly and the discrete compactness property of in Lemma 2.1, we know that is a relatively compact set in , which implies collectively compact convergence . Noting are self-adjoint, due to Proposition 3.7 or Table 3.1 in [17] we get This ends the proof.     

Prior to proving the error estimates for edge element discretizations, we define some notations as follows. Let be the th eigenvalue of (LABEL:2.5) or (LABEL:3.1ss)-(LABEL:3.2ss) of multiplicity . Let be eigenvalues of that converge to the eigenvalue . Here and hereafter we use to denote the space spanned by all eigenfunctions corresponding to the eigenvalue , and to denote the direct sum of all eigenfunctions corresponding to the eigenvalues . For argument convenience, hereafter we denote and . Now we introduce the following small quantity:

Thanks to (LABEL:s2.18) and (LABEL:2.7) we have

\hb@xt@.01(2.23)

The error estimates of edge elements for the Maxwell eigenvalue problem have been obtained in, e.g., [4, 6, 29, 32]. Here we would like to use the quantity to characterize the error for eigenpairs. From the spectral approximation, we actually derive the a priori error estimates for the discrete eigenvalue problem (LABEL:2.14) with , (LABEL:2.13) or (LABEL:3.1s)-(LABEL:3.2s).

Theorem 2.3

Let be the eigenvalue of (LABEL:2.5) or (LABEL:3.1ss)-(LABEL:3.2ss) and let be the discrete eigenvalue of (LABEL:2.13) or (LABEL:3.1s)-(LABEL:3.2s) converging to . There exist such that if then for any eigenfunction corresponding to with there exists such that

\hb@xt@.01(2.24)

and for any with there exists such that

\hb@xt@.01(2.25)

where the positive constant is independent of mesh parameters.

Proof. We take . Suppose is an eigenfunction of (LABEL:3.1s)-(LABEL:3.2s) corresponding to satisfying . Then according to Theorems 7.1 and 7.3 in [11] and Lemma 2.2 there exists satisfying

\hb@xt@.01(2.26)
\hb@xt@.01(2.27)

By a simple calculation, we deduce

Since the equality (LABEL:2.21) implies , this together with (LABEL:2.28)-(LABEL:2.29) yields (LABEL:2.22). Conversely, suppose is an eigenfunction of (LABEL:3.1ss)-(LABEL:3.2ss) corresponding to satisfying . Then according to Theorems 7.1 in [11] and Lemma 2.1 there exists satisfying

\hb@xt@.01(2.28)

Let where is the eigenfunction corresponding to such that constitutes an orthogonal basis of in . Then

Since , this together with (LABEL:2.29)-(LABEL:2.30) yields (LABEL:2.23).     

Remark 2.1. Based on the estimate (LABEL:2.22), one can naturally obtain the optimal convergence order for using the Rayleigh quotient relation (LABEL:s2.12) in the following section. In addition, note that when in Theorem 2.2 the estimate (LABEL:2.22) implies converges to . Here we introduce then and (LABEL:2.22) gives

\hb@xt@.01(2.29)

For simplicity of notation, we still use the same and in the above estimate as in (LABEL:2.22)-(LABEL:2.23).

Remark 2.2. When is a Lipschitz polyhedron and are properly smooth, it is known that () (see [19, 5, 28]) and . In particular, if () then (see Theorem 5.41 in [28]).

3 Multigrid schemes based on shifted inverse iteration

3.1 Multigrid Schemes

In practical computation, the information on the physical zero eigenvalues can be easily captured on a coarse mesh using the mixed discretization (LABEL:3.1s)-(LABEL:3.2s). In this section we shall present our multigrid methods for solving nonzero Maxwell eigenvalue. The following schemes are proposed by [34, 35]. Note that we assume in the following schemes the numerical eigenvalue approximates the nonzero eigenvalue .
Scheme 3.1.   Rayleigh quotient iteration.
Given the maximum number of iterative times .
Step 1. Solve the eigenvalue problem (LABEL:s1)-(LABEL:s3) on coarse finite element space : find , such that

Step 2. .
Step 3. Solve an equation on : find such that

Set
Step 4. Compute the Rayleigh quotient

Step 5. If , then output , stop; else, , and return to step 3.

In Step 3 of the above Scheme, when the shift is close to the exact eigenvalue enough, the coefficient matrix of linear equation is nearly singular. Hence the following algorithm gives a natural way of handling this problem.

Scheme 3.2.   Inverse iteration with fixed shift.
Given the maximum number of iterative times and .
Step 1-Step 4. The same as Step 1-Step 4 of Scheme 3.1.
Step 5. If then and return to Step 6; else and return to Step 3.
Step 6. Solve an equation on : find such that

Set .
Step 7. Compute the Rayleigh quotient

Step 8. If , then output , stop; else, , and return to step 6.

Remark 3.1. The mixed discretization (LABEL:3.1s)-(LABEL:3.2s) was adopted by the literatures [25, 1, 24]. As is proved in Theorem 2.2, using this discretization we can compute the Maxwell eigenvalues without introducing spurious eigenvalues. However, it is also a saddle point problem that is difficult to solve on a fine mesh (see [1, 2]). Therefore, the multigrid schemes can properly overcome the difficulty since we only solve (LABEL:3.1s)-(LABEL:3.2s) on a coarse mesh, as shown in step 1 of Schemes 3.1 and 3.2. Moreover, in order to further improve the efficiency of solving the equation in Steps 3 and 6 in Schemes 3.1 and 3.2 the HX preconditioner in [21] is a good choice (see [37]).

3.2 Error Analysis

In this subsection, we aim to prove the error estimates for Schemes 3.1 and 3.2. We shall analyze the constants in the error estimates are independent of mesh parameters and iterative times . First of all, we give two useful lemmas.

Lemma 3.1

For any nonzero there hold

\hb@xt@.01(3.1)

Proof. See [35].     

Lemma 3.2

Let be an eigenpair of (LABEL:2.5) or of (LABEL:s5) with , then for any , the Rayleigh quotient satisfies

\hb@xt@.01(3.2)

Proof. See pp.699 of [11].     

The basic relation in Lemma 3.2 cannot be directly applied to our theoretical analysis, so in the following we shall further simplify the estimate (LABEL:s2.12). Let then according to the definition of ,

If , , and , then by Lemma 3.1 we deduce

which together with yields

Hence, from Lemma LABEL:l2 we get the following estimate

\hb@xt@.01(3.3)

where Define the operators and as

\hb@xt@.01(3.4)
\hb@xt@.01(3.5)

The following lemma turns our attention from the spectrum of and to that of and .

Lemma 3.3

, and share the eigenvalues greater than and the associated eigenfunctions. The same conclusion is valid for , and . Moreover, and .

Proof. The assertions regarding the relations among , , and have been described in section 2. Next we shall prove the relations between and and between and . By the definition of and , the eigenpair () of satisfies and the eigenpair () of satisfies . Note that the above two weak forms are equivalent when (since this implies the eigenfunction of the latter satisfies divergence-free constraint). Hence and share the eigenvalues and the associated eigenfunctions. Similarly one can check and share the eigenvalues and the associated eigenfunctions. Thanks to Helmholtz decomposition and (LABEL:s2.7), we also have for any

This together with (LABEL:3.4) yields . Thanks to discrete Helmholtz decomposition and (LABEL:s2.5), we also have for any

This together with (LABEL:3.5) yields .     

Denote . For better understanding of notations, hereafter we write , and .
The following lemma (see [34]) is valid since and share the same eigenpairs. It provides a crucial tool for analyzing the error of multigrid Schemes 3.1 and 3.2.

Lemma 3.4

Let is an approximate eigenpair of , where is not an eigenvalue of and with . Suppose that