aa r X i v : . [ m a t h . NA ] J a n LOCAL MULTIGRID IN H ( curl ) RALF HIPTMAIR ∗ AND
WEIYING ZHENG † Abstract.
We consider H ( curl , Ω)-elliptic variational problems on bounded Lipschitz polyhedraand their finite element Galerkin discretization by means of lowest order edge elements. We assumethat the underlying tetrahedral mesh has been created by successive local mesh refinement, eitherby local uniform refinement with hanging nodes or bisection refinement. In this setting we developa convergence theory for the the so-called local multigrid correction scheme with hybrid smoothing.We establish that its convergence rate is uniform with respect to the number of refinement steps. Theproof relies on corresponding results for local multigrid in a H (Ω)-context along with local discreteHelmholtz-type decompositions of the edge element space. Key words.
Edge elements, local multigrid, stable multilevel splittings, subspace correctiontheory, regular decompositions of H ( curl , Ω), Helmholtz-type decompositions, local mesh refinement
AMS subject classifications.
1. Introduction.
On a polyhedron Ω ⊂ R , scaled such that diam(Ω) = 1, weconsider the variational problem: seek u ∈ H Γ D ( curl , Ω) such that( curl u , curl v ) L (Ω) + ( u , v ) L (Ω) | {z } =: a ( u , v ) = ( f , v ) L (Ω) ∀ v ∈ H Γ D ( curl , Ω) . (1.1)For the Hilbert space of square integrable vector fields with square integrable curl and vanishing tangential components on Γ D we use the symbol H Γ D ( curl , Ω), see [24,Ch. 1] for details. The source term f in (1.1) is a vector field in ( L (Ω)) . The left handside of (1.1) agrees with the inner product of H Γ D ( curl , Ω) and will be abbreviatedby a ( u , v ) (“energy inner product”).Further, Γ D denotes the part of the boundary ∂ Ω on which homogeneous Dirichletboundary conditions in the form of vanishing tangential traces of u are imposed. Thegeometry of the Dirichlet boundary part Γ D is supposed to be simple in the followingsense: for each connected component Γ i of Γ D we can find an open Lipschitz domainΩ i ⊂ R such that Ω i ∩ Ω = Γ i , Ω i ∩ Ω = ∅ , (1.2)and Ω i and Ω j have positive distance for i = j . Further, the interior of Ω ∪ Ω ∪ Ω . . . is expected to be a Lipschitz-domain, too (see Fig. 5.2). This is not a severerestriction, because variational problems related to (1.1) usually arise in quasi-staticelectromagnetic modelling, where simple geometries are common. Of course, Γ D = ∅ is admitted.Lowest order H Γ D ( curl , Ω)-conforming edge elements are widely used for thefinite element Galerkin discretization of variational problems like (1.1). Then, for asolution u ∈ ( H (Ω)) with curl u ∈ ( H (Ω)) we can expect the optimal asymptoticconvergence rate k u − u h k H ( curl , Ω) ≤ CN − / h , (1.3) ∗ SAM, ETH Z¨urich, CH-8092 Z¨urich, [email protected] † LSEC, Institute of Computational Mathematics, Academy of Mathematics and System Sciences,Chinese Academy of Sciences, Beijing, 100080, People’s Republic of China. This author was supportedin part by China NSF under the grant 10401040 ([email protected]).1
R. Hiptmair and W.-Y. Zheng H ( curl , Ω) : Sobolev space of square integrable vector fields on Ω ⊂ R withsquare integrable curl H Γ D ( curl , Ω) : vector fields in H ( curl , Ω) with vanishing tangential compo-nents on Γ D ⊂ ∂ Ω M , T : tetrahedral finite element meshes, may contain hanging nodes N ( M ) : set of vertices (nodes) of a mesh ME ( M ) : set of edges of a mesh M ρ K , ρ M : shape regularity measures h : – local meshwidth function for a finite element mesh– (as subscript) tag for finite element functions U ( M ) : lowest order edge element space on M b E : nodal basis function of U ( M ) associated with edge EV ( M ) : space of continuous piecewise linear functions on M V ( M ) : quadratic Lagrangian finite element space on M e V ( M ) : quadratic surplus space, see (2.19) b p : nodal basis function of V ( M ) (“tent function”) associated withvertex p B X ( M ) : set of nodal basis functions for finite element space X on mesh M Π h : nodal edge interpolation operator onto U ( M ), see (2.7) I h : vertex based piecewise linar interpolation onto V ( M ) P p : space of 3-variate polynomials of total degree ≤ p U ( M ), V ( M ): finite element spaces oblivious of zero boundary conditions ≺ : nesting of finite element meshes ℓ ( K ) : level of element K in hierarchy of refined meshes ω l : refinement zone, see (4.1)Σ l : refinement strip, see (5.57) B lV , B l U : sets of basis functions supported inside refinement zones, see(4.9) Q h : quasi-interpolation operator, Def. 5.1 Table 0.1
Important notation used in this paper on families of finite element meshes arising from global refinement. Here, u h is thefinite element solution, N h the dimension of the finite element space, and C > N h . However, often u will fail to possess the required regularitydue to singularities arising at edges/corners of ∂ Ω and material interfaces [22, 23].Fortunately, it seems to be possible to retain (1.3) by the use of adaptive local meshrefinement based on a posteriori error estimates, see [10, 55] for theory in H -setting,[7, 17] for numerical evidence in the case of edge element discretization, and [8, 34, 52]for related theoretical investigations.We also need ways to compute the asymptotically optimal finite element solutionwith optimal computational effort, that is, with a number of operations proportionalto N h . This can only be achieved by means of iterative solvers, whose convergenceremains fast regardless of the depth of refinement. Multigrid methods are the mostprominent class of iterative solvers that achieve this goal. By now, geometric multigridmethods for discrete H ( curl , Ω)-elliptic variational problems like (1.1) have becomewell established [20, 29, 54, 58]. Their asymptotic theory on sequencies of regularly ocal Multigrid in H ( curl ) asymptotic optimality :the speed of convergence is uniformly fast regardless of the number of refinementlevels involved. In addition, the costs of one step of the iteration scale linearly withthe number of unknowns.Yet, the latter property is lost when the standard multigrid correction scheme isapplied to meshes generated by pronounced local refinement. Optimal computationalcosts can only be maintained, if one adopts the local multigrid policy, which waspioneered by A. Brandt et al. in [5], see also [41]. Crudely speaking, its gist is toconfine relaxations to “new” degrees of freedom located in zones where refinementhas changed the mesh. Thus an exponential increase of computational costs withthe number of refinement levels can be avoided: the total costs of a V-cycle remainproportional to the number of unknowns. An algorithm blending the local multigrididea with the geometric multigrid correction scheme of [29] is described in [54]. Onthe other hand, a proof of uniform asymptotic convergence has remained elusive sofar. It is the objective of this paper to provide it, see Theorem 4.2.We recall the key insight that (1.1) is one member of a family of variationalproblems. Its kin is obtained by replacing curl with grad or div, respectively. Allthese differential operators turn out to be incarnations of the fundamental exteriorderivative of differential geometry, cf. [29, Sect. 2]. They are closely connected in thedeRham complex [3] and, thus, it is hardly surprising that results about the related H D (Ω)-elliptic variational problem, which seeks u ∈ H D (Ω) such that( grad u, grad v ) L (Ω) + ( u, v ) L (Ω) = ( f, v ) L (Ω) ∀ v ∈ H D (Ω) , (1.4)prove instrumental in the multigrid analysis for discretized versions of (1.1). Here H D (Ω) is the subspace of H (Ω) whose functions have vanishing traces on Γ D .Thus, when tackling (1.1), we take the cue from the local multigrid theory for(1.4) discretized by means of linear continuous finite elements. This theory has beendeveloped in various settings, cf. [5,11,14,15,62]. In [1] local refinement with hangingnodes is treated. Recently, H. Wu and Z. Chen [60] proved the uniform convergence oflocal multigrid V-cycles on adaptively refined meshes in two dimensions. Their meshrefinements are controlled by a posteriori error estimators and carried out accordingto the “newest vertex bisection” strategy introduced, independently, in [6, 40].As in the case of global multigrid, the essential new aspect of local multigridtheory for (1.1) compared to (1.4) is the need to deal with the kernel of the curl -operator, cf. [29, Sect. 3]. In this context, the availability of discrete scalar potentialrepresentations for irrotational edge element vector fields is pivotal. Therefore, wedevote the entire Sect. 2 to the discussion of edge elements and their relationshipwith conventional Lagrangian finite elements. Meshes with hanging nodes will receiveparticular attention. Next, in Sect. 3 we present details about local mesh refinement,because some parts of the proofs rest on the subtleties of how elements are split. Thefollowing Sect. 4 introduces the local multigrid method from the abstract perspectiveof successive subspace correction.The proof of uniform convergence (Theorem 4.2) is tackled in Sects. 5 and 6, whichform the core of the article. In particular, the investigation of the stability of the localmultilevel splitting requires several steps, the first of which addresses the issue for thebilinear form from (1.4) and linear finite elements. These results are already availablein the literature, but are re-derived to make the presentation self-contained. Thisalso applies to the continuous and discrete Helmholtz-type decompositions covered inSect. 5.3. Many developments are rather technical and to aid the reader important R. Hiptmair and W.-Y. Zheng notations are listed in Table 0.1. Eventually, in Sect. 7, we report two numericalexperiments to show the competitive performance of the local multigrid method andthe relevance of the convergence theory.
Remark 1.1.
In this article we forgo generality and do not discuss the moregeneral bi-linear form a ( u , v ) := ( α curl u , curl v ) L (Ω) + ( β u , v ) L (Ω) , ∀ u , v ∈ H Γ D ( curl , Ω) , (1.5) with uniformly positive coefficient functions α, β ∈ L ∞ (Ω) . We do this partly for thesake of lucidity and partly, because the current theory cannot provide estimates thatare robust with respect to large variations of α and β , cf. [33]. We refer to [63] forfurther information and references.
2. Finite element spaces.
Whenever we refer to a finite element mesh in thisarticle, we have in mind a tetrahedral triangulation of Ω, see [19, Ch. 3]. In certainsettings, it may feature hanging nodes, that is, the face of one tetrahedron can coincidewith the union of faces of other tetrahedra. Further, the mesh is supposed to resolvethe Dirichlet boundary in the sense that Γ D is the union of faces of tetrahedra. Thesymbol M with optional subscripts is reserved for finite element meshes and the setsof their elements alike.We write h ∈ L ∞ (Ω) for the piecewise constant function, which assumes value h K := diam( K ) in each element K ∈ M . The ratio of diam( K ) to the radius of thelargest ball contained in K is called the shape regularity measure ρ K [19, Ch. 3, § ρ M of M is the maximum of all ρ K , K ∈ M . Provisionally, we consider only finite elementmeshes M that are conforming, that is, each face of a tetrahedron is either con-tained in ∂ Ω or a face of another tetrahedron, see [19, Ch. 2, § H Γ D ( curl , Ω)-conforming edge finite elements, also known as Whitney-1-forms [59], U ( M ) := { v h ∈ H Γ D ( curl , Ω) : ∀ K ∈ M : ∃ a , b ∈ R :(2.1) v h ( x ) = a + b × x , x ∈ K } . (2.2)For a detailed derivation and description please consult [30, Sect. 3] or the monographs[13, 42]. Notice that curl U ( M ) is a space of piecewise constant vector fields. We alsoremark that appropriate global degrees of freedom (d.o.f.) for U ( M ) are given by (cid:26) U ( M ) R v h R E v h · d ~s , E ∈ E ( M ) , (2.3)where E ( M ) is the set of edges of M not contained in Γ D . We write B U ( M ) for thenodal basis of U ( M ) dual to the global d.o.f. (2.3). Basis functions are associatedwith active edges. Hence, we can write B U ( M ) = { b E } E ∈E ( M ) . The support of thebasis function b E is the union of tetrahedra sharing the edge E . We recall the simpleformula for local shape functions b E | K = λ i grad λ j − λ j grad λ i E = [ a i , a j ] ⊂ K (2.4)for any tetrahedron K ∈ M with vertices a i , i = 1 , , ,
4, and associated barycentriccoordinate functions λ i . ocal Multigrid in H ( curl ) U ( M ) with basis B U ( M ) is perfectly suited for the finiteelement Galerkin discretization of (1.1). The discrete problem based on U ( M ) reads:seek u h ∈ U ( M ) such that( curl u h , curl v h ) L (Ω) + ( u h , v h ) L (Ω) = ( f , v h ) L (Ω) ∀ v h ∈ U ( M ) . (2.5)The properties of U ( M ) will be key to constructing and analyzing the local multigridmethod for the large sparse linear system of equations resulting from (2.5). Next, wecollect important facts.The basis B U ( M ) enjoys uniform L -stability, meaning the existence of a con-stant C = C ( ρ M ) > v h = P E ∈E ( M ) α E b E ∈ U ( M ), α E ∈ R , C − k v h k L (Ω) ≤ X E ∈E ( M ) α E k b E k L (Ω) ≤ C k v h k L (Ω) . (2.6)The global d.o.f. induce a nodal edge interpolation operator Π h : dom( Π h ) ⊂ H Γ D ( curl , Ω) U ( M ) v P E ∈E ( M ) (cid:16)R E v · d ~s (cid:17) · b E . (2.7)Obviously, Π h provides a local projection, but it turns out to be unbounded evenon ( H (Ω)) . Only for vector fields with discrete rotation the following interpolationerror estimate is available, see [30, Lemma 4.6]: Lemma 2.1.
The interpolation operator Π h is bounded on { Ψ ∈ ( H (Ω)) , curl Ψ ∈ curl U ( M ) }⊂ ( H (Ω)) , and for any conforming mesh thereis C = C ( ρ M ) > such that (cid:13)(cid:13) h − ( Id − Π h ) Ψ (cid:13)(cid:13) L (Ω) ≤ C | Ψ | H (Ω) ∀ Ψ ∈ ( H (Ω)) , curl Ψ ∈ curl U ( M ) . If Ω is homeomorphic to a ball, then grad H (Ω) = H ( curl , Ω) := { v ∈ H ( curl , Ω) , curl v = 0 } , that is, H (Ω) provides scalar potentials for H ( curl , Ω).To state a discrete analogue of this relationship we need the Lagrangian finite ele-ment space of piecewise linear continuous functions on M V ( M ) := { u h ∈ H D (Ω) : u h | K ∈ P ( K ) ∀ K ∈ M} , (2.8)where P p ( K ) is the space of 3-variate polynomials of degree ≤ p on K . The globaldegrees of freedom for V ( M ) boil down to point evaluations at the vertices of M away from Γ D (set N ( M )). The dual basis of “tent functions” will be denoted by B V ( M ) = { b p } p ∈N ( M ) . Its unconditional L -stability is well known: with a universalconstant C > u h = P p ∈N ( M ) α p b p ∈ V ( M ), α p ∈ R , C − k u h k L (Ω) ≤ X p ∈N ( M ) α p k b p k L (Ω) ≤ C k u h k L (Ω) . (2.9) The symbol C will stand for generic positive constants throughout this article. Its value mayvary between different occurrences. We will always specify on which quantities these constants maydepend. R. Hiptmair and W.-Y. Zheng
For the nodal interpolation operator related to B V we write I h : dom( I h ) ⊂ H D (Ω) V ( M ). Recall the standard estimate for linear interpolation on conforming meshes (i.e., no hanging nodes allowed), [19, Thm. 3.2.1], that asserts the existenceof C = C ( k, ρ M ) > (cid:13)(cid:13) h k − ( Id − I h ) u (cid:13)(cid:13) H k (Ω) ≤ C | u | H (Ω) ∀ u ∈ H (Ω) ∩ H D (Ω) , k ∈ { , } . (2.10)Obviously, grad V ( M ) ⊂ U ( M ), and immediate from Stokes theorem is the crucial commuting diagram property Π h ◦ grad = grad ◦I h on dom( I h ) . (2.11)This enables us to give an elementary proof of Lemma 2.1. Proof . [of Lemma 2.1] Pick one K ∈ M and, without loss of generality, assume0 ∈ K . Then define the lifting operator, cf. the “Koszul lifting” [3, Sect. 3.2], w
7→ L w , L w ( x ) := w ( x ) × x , x ∈ K . (2.12)Elementary calculations reveal that for any constant vectorfield w ∈ ( P ( K )) curl L w = w , (2.13) kL w k L ( K ) ≤ h K k w k L ( K ) , (2.14) L w ∈ U ( K ) . (2.15)The continuity (2.14) permits us to extend L to ( L ( K )) .Given Ψ ∈ ( H ( K )) with curl Ψ ≡ const , by (2.15) we know L curl Ψ ∈ ( P ( K )) . Thus, an inverse inequality leads to |L curl Ψ | H ( K ) ≤ Ch − K kL curl Ψ k L ( K ) (2.14) ≤ C k curl Ψ k L ( K ) , (2.16)with C = C ( ρ K ) >
0. Next, (2.13) implies curl ( Ψ − L curl Ψ ) = 0 ⇒ ∃ p ∈ H ( K ) : Ψ − L curl Ψ = grad p . (2.17)From (2.16) we conclude that p ∈ H ( K ) and | p | H ( K ) ≤ C | Ψ | H ( K ) . Moreover,thanks to the commuting diagram property we have Ψ − Π h Ψ = L curl Ψ − Π h L curl Ψ | {z } =0 by (2.15) + grad ( p − I h p ) , (2.18)which means, by the standard estimate (2.10) for linear interpolation on K , k Ψ − Π h Ψ k L ( K ) = | p − I h p | H ( K ) ≤ Ch K | p | H ( K ) ≤ Ch K | Ψ | H ( K ) . Summation over all elements finishes the proof.As theoretical tools we need “higher order” counterparts of the above finite ele-ment spaces. We recall the quadratic Lagrangian finite element space V ( M ) := { u h ∈ H D (Ω) : u h | K ∈ P ( K ) ∀ K ∈ M} , (2.19)and its subspace of quadratic surpluses e V ( M ) := { u h ∈ V ( M ) : I h u h = 0 } . (2.20) ocal Multigrid in H ( curl ) V ( M ) = V ( M ) ⊕ e V ( M ) , (2.21)which is unconditionally H -stable: there is a C = C ( ρ M ) > C − | u h | H (Ω) ≤ | ( Id − I h ) u h | H (Ω) + |I h u h | H (Ω) ≤ C | u h | H (Ω) , (2.22)for all u h ∈ V ( M ).Next, we examine the space ( V ( M )) of continuous piecewise linear vector fieldsthat vanish on Γ D . Standard affine equivalence techniques for edge elements, see [30,Sect. 3.6], confirm ∃ C = C ( ρ M ) > k Π h Ψ h k L (Ω) ≤ C k Ψ h k L (Ω) ∀ Ψ h ∈ ( V ( M )) . (2.23) Lemma 2.2.
For all Ψ h ∈ ( V ( M )) we can find e v h ∈ e V ( M ) such that Ψ h = Π h Ψ h + grad e v h , and, with C = C ( ρ M ) > , C − k Ψ h k L (Ω) ≤ k Π h Ψ h k L (Ω) + k grad e v h k L (Ω) ≤ C k Ψ h k L (Ω) . For the proof we rely on a very useful insight, which relieves us from all worriesconcerning the topology of Ω:
Lemma 2.3. If v ∈ H Γ D ( curl , Ω) and Π h v = 0 , then v ∈ grad H D (Ω) .Proof . Since the mesh covers Ω, the relative homology group H (Ω; Γ D ) is gen-erated by a set of edge paths. By definition (2.3) of the d.o.f. of U ( M ), the pathintegrals of v along all these paths vanish. As an irrotational vector field with van-ishing circulation along a complete set of Γ D -relative fundamental cycles, v must bea gradient. Proof . [of Lemma 2.2] Given Ψ h ∈ ( V ( M )) , we decompose it according to Ψ h = Π h Ψ h + ( Id − Π h ) Ψ h | {z } =: grad e v h . (2.24)Note that curl ( Id − Π h ) Ψ h is piecewise constant with vanishing flux through alltriangular faces of M . Then Stokes’ theorem teaches that curl ( Id − Π h ) Ψ h = 0.By the projector property of Π h , ( Id − Π h ) Ψ h satisfies the assumptions ofLemma 2.3. Taking into account that, moreover, the field is piecewise linear, it isclear that ( Id − Π h ) Ψ h = grad ψ with ψ ∈ V ( M ). Along an arbitrary edge path γ in M we have R γ ( Id − Π h ) Ψ h · d ~s = 0 so that ψ attains the same value (w.l.o.g. = 0)on all vertices of M . The stability of the splitting is a consequence of (2.23).By definition, the spaces U ( M ) and V ( M ) accommodate the homogeneousboundary conditions on Γ D . Later, we will also need finite element spaces oblivi-ous of boundary conditions, that is, for the case Γ D = ∅ . These will be tagged by abar on top, e.g., U ( M ), V ( M ), etc. The same convention will be employed for notionsand operators associated with finite element spaces: if they refer to the particular caseΓ D = ∅ , they will be endowed with an overbar, e.g. Π h , I h , B U ( M ), N ( M ), etc. R. Hiptmair and W.-Y. Zheng
Now, general tetrahedral meshes withhanging nodes are admitted. We simply retain the definitions (2.8) and (2.19) of thespaces V ( M ) and V ( M ) of continuous finite element functions. Degrees of freedomfor V ( M ) are point evaluations at active vertices of M . A vertex is called active, ifit is not located in the interior of an edge/face of M or on Γ D . A 2D illustration isgiven in Fig. 2.1. M M M M Fig. 2.1 . Active vertices ( • ) of 2D triangular meshes with hanging nodes, Ω =]0 , , Γ D = ∂ Ω .In M , M , M active edges are marked with green arrows. The values of a finite element function at the remaining (“slave”) vertices aredetermined by recursive affine interpolation. A dual nodal basis B V ( M ) and corre-sponding interpolation operator I h can be defined as above.In principle, the definition (2.1) of the edge element space could be retained onnon-conforming meshes, as well. Yet, for this choice an edge interpolation operator Π h that satisfies the commuting diagram property (2.11) is not available. Thus, weconstruct basis functions directly and rely on the notion of active edges , see Fig. 2.1. Definition 2.4.
An edge of M is active , if it is an edge of some K ∈ M , notcontained in Γ D , and connects two vertices that are either active or located on Γ D . We keep the symbol E ( M ) to designate the set of active edges of M . To each E ∈ E ( M ) we associate a basis function b E , which, locally on the tetrahedra of M , isa polynomial of the form (2.2). In order to fix this basis function completely, it sufficesto speficify its path integrals (2.3) along all edges of M . In the spirit of duality, wedemand Z F b E · d ~s = ( F = E , F ∈ E ( M ) \ { E } . (2.25)For the non-active (“slave”) edges of M the path integrals of b E (subsequently called“weights”) are chosen to fit (2.11), keeping in mind that B U ( M ) := { b E } E ∈E ( M ) ,and that the d.o.f. and Π h are still defined according to (2.3) and (2.7), respectively.Ultimately, we set U ( M ) := Span { B U ( M ) } .Let us explain the policy for setting the weights in the case of the subdividedtetrahedron of Fig. 2.2 with hanging nodes at the midpoints of edges, which will turnout to be the only relevant situation, cf. Sect. 5.3. Weights have to be assigned to the“small edges” of the refined tetrahedron, some of which will be active, and some ofwhich will have “slave” status, see the caption of Fig. 2.2. For ease of visualization, we will often elucidate geometric concepts in two-dimensional settings.Their underlying ideas are the same in 2D and 3D.ocal Multigrid in H ( curl ) q − p = ( p − p ) , q − q = ( p + p ) − ( p + p ) = ( p − p ) , p − q = p − ( p + p ) = p − p + ( p − p ) , q − q = ( p + p ) − ( p + p ) = ( p − p ) + ( p − p ) . In a sence, we express slave edges as “linear combinations” of active edges. In adifferent context, this policy is explained in more detail in [26].PSfrag replacements p p p p p p q q q q E Fig. 2.2 . Subdivided tetrahedron, active vertices ( • ) p , . . . , p , slave vertices ( ◦ ) q , . . . , q ,active edges [ p , p ] , [ p , p ] , [ p , p ] , [ p , p ] , [ p , p ] , [ p , p ] , [ p , p ] , [ p , p ] , [ p , p ] , slaveedges [ p , q ] , [ q , p ] , [ p , q ] , [ q , p ] , [ p , q ] , [ q , p ] , [ p , q ] , [ q , q ] , [ p , q ] , [ q , q ] , [ q , q ] , [ q , q ] , [ q , p ] , [ q , p ] , [ q , p ] The coefficients in the combinations tell us the weights. For example, for the activeedge E = [ p , p ] in Fig. 2.2 they are given in Table 2.1. Using these weights and theformula (2.4), b E can be assembled on the tetrahedron by imposing (see Table 2.1 fornotations) Z S b E · d ~s = ( w S for any contributing slave edge S , , S ∈ { “small edges” } . Firstly, the procedure for the selection of weight guarantees that grad V ( M ) ⊂ U ( M ). For illustration, we single out the gradient w h of the nodal basis func-tion belonging to vertex p in Fig. 2.2. Its path integral equals 1 along the (ori-ented) edges [ p , p ] , [ p , p ] , [ p , p ] , [ q , p ] , [ q , p ] , [ q , p ], and vanishes on all0 R. Hiptmair and W.-Y. Zheng
Slave edge S [ q , p ] [ p , q ] [ q , q ] [ q , p ] [ q , q ]weight w S
12 12 12 12 − Table 2.1
Weights for slave edges in Fig. 2.2 relative to active edge E = [ p , p ] . Only slave edges withnon-zero weights are listed. other edges. Hence we expect w h = b [ p , p ] + b [ p , p ] + b [ p , p ] . (2.26)This can be verified through showing equality of path integrals along slave edges.We take a close look at the slave edge [ q , p ]. By construction the basis functionsbelonging to active edges satisfy Z [ q , p ] b [ p , p ] · d ~s = 1 , Z [ q , p ] b [ p , p ] · d ~s = 0 , Z [ q , p ] b [ p , p ] · d ~s = 0 , Z [ q , p ] b [ p , p ] · d ~s = − , Z [ q , p ] b [ p , p ] · d ~s = 0 Z [ q , p ] b [ p , p ] · d ~s = 0 , Z [ q , p ] b [ p , p ] · d ~s = 0 . Then, evidently,1 = Z [ q , p ] w h · d ~s = Z [ q , p ] b [ p , p ] · d ~s + Z [ q , p ] b [ p , p ] · d ~s + Z [ q , p ] b [ p , p ] · d ~s = 1 + 0 + 0 . The same considerations apply to all other slave edges and (2.26) is established.Secondly, the construction ensures the commuting diagram property (2.11): againappealing to Fig. 2.2 we find, for example, Z [ q , q ] grad I h u · d ~s = I h u ( q ) − I h u ( q ) = ( u ( p ) + u ( p )) − ( u ( p ) + u ( p )) = Z [ p , p ] grad u · d ~s + Z [ p , p ] grad u · d ~s = Z [ p , p ] grad u · d ~s + Z [ p , p ] grad u · d ~s + Z [ p , p ] grad u · d ~s . In words, combining the path integrals of grad u along active edges with the relativeweights of the slave edge [ q , q ] yields the same result as evaluating the path integralof the gradient of the interpolant I h u along [ q , q ].The definitions (2.19) and (2.20) also carry over to meshes with hanging nodes.This remains true for the splitting asserted in Lemma 2.2. However, though the al-gebraic relationships like (2.11) remain valid, the estimates and norm equivalences of ocal Multigrid in H ( curl ) cf. Assumption 3.0.1.
Remark 2.5.
Our presentation is confined to tetrahedral meshes and lowest orderedge elements for the sake of simplicity. Extension of all results to hexahedral meshesand higher order edge elements is possible, but will be technical and tedious.
3. Local mesh refinement.
We study the case where the actual finite elementmesh M h of Ω has been created by successive local refinement of a relatively uniforminitial mesh M . Concerning M h and M the following asumptions will be made:1. Given M and M h we can construct a virtual refinement hierarchy of L + 1 nested tetrahedral meshes, L ∈ N : M ≺ M ≺ M ≺ · · · ≺ M L = M h . (3.1)Please note that the virtual refinement hierarchy may be different from theactual sequence of meshes spawned during adaptive refinement .2. Inductively, we assign to each tetrahedron K ∈ M l a level ℓ ( K ) ∈ N bycounting the number of subdivisions it took to generate it from an elementof M .3. For all 0 ≤ l < L the mesh M l +1 is created by subdividing some or all of thetetrahedra in { K ∈ M l : ℓ ( K ) = l } .4. The shape regularity measures of the meshes M l are uniformly boundedindependently of L .Refinement may be local, but it must be regular in the following sense, cf. [46,Sect. 4.2.2] and [60]: we can find a second sequence of nested tetrahedral meshes of Ω M = c M ≺ c M ≺ c M ≺ · · · ≺ c M L . (3.2)that satisfies1. M l ≺ c M l and { K ∈ M l : ℓ ( K ) = l } ⊂ c M l , l = 0 , . . . , L ,2. that the shape regularity measure ρ c M l is bounded independently of l ,3. and that there exist two constants C > < θ < l and L such that C − θ l ≤ h K ≤ Cθ l ∀ K ∈ c M l , ≤ l ≤ L . (3.3)This means that the family { c M l } l is quasi-uniform. Hence, it makes sense torefer to a mesh width h l := max { h K , K ∈ c M l } of c M l . It decreases geomet-rically for growing l .Our analysis targets two popular tetrahedral refinement schemes that generatesequences of meshes that meet the above requirements. This scheme produces M l +1 by splitting someof the tetrahedra of the current mesh M l into eight smaller ones, possibly creatinghanging nodes in the process [1]. An illustrative 2D example with hanging nodes isdepicted in Figure 3.1. The accompanying sequence { c M l } ≤ l ≤ L is produced by globalregular refinement, which implies (3.3) with θ = . Uniform shape-regularity can alsobe guaranteed for repeated regular refinement of tetrahedra, see [9]. two finite element meshes M and T are nested, M ≺ T , if every element of M is the union of R. Hiptmair and W.-Y. Zheng M M M M = M h Fig. 3.1 . Virtual refinement hierarchy for 2D triangular meshes. The quasi-uniform sequence { c M l } ≤ l ≤ L is sketched in blue. Elements of M l eligible for further subdivision are marked yellow. The meshes occurring in the virtual refinement hierarchy need not agree with themeshes that arise during adaptive refinement in an actual computation. Yet, given M h , the virtual refinement hierarchy can always be found a posteriori . Write M hier for the union of all tetrahedra ever created during the refinement process. Then, for0 < l < L , define M l := (cid:26) K ∈ M hier : ℓ ( K ) ≤ l and K does not contain a K ′ ∈ M hier \ { K } with ℓ ( K ′ ) ≤ l (cid:27) . (3.4)Using the construction of finite element spaces detailed in Sect. 2.2, the localmultigrid algorithms can handle any kind of local regular refinement. Yet, convergencemay degrade unless we curb extreme jumps of local meshwidth. Thus, we assume thefollowing throughout the remainder of this paper. Assumption 3.0.1.
Any edge of M h may contain at most one hanging node. This will automatically be satisfied for all meshes M l of the virtual refinementhierarchy. Consequently, hanging nodes can occur only in a few geometric configu-rations, one of which is depicted in Fig. 2.2. This paves the way for using mappingtechniques and scaling arguments, see [30, Sect. 3.6], which confirm the following gen-eralization of results of Sect. 2.1. Of course, we rely on the constructions of finiteelement spaces and interpolation operators described in Sect. 2.2. Proposition 3.1.
Under Assumption 3.0.1 the L -stability of bases, see (2.6) , (2.9) , carries over uniformly to meshes created by local regular refinement. So doLemmas 2.1, 2.2, and Estimates (2.22) , (2.23) . Summing up, Assumption 3.0.1 makes it possible to use the results obtained inSect. 2.1 in the case of local regular refinement as well. To avoid a proliferation oflabels, we are going to quote the statements from Sect. 2.1 even when we mean theirgeneralization to meshes with hanging nodes.
This procedure involves splitting atetrahedron into two by promoting the midpoint of the so-called refinement edgeto a new vertex. Variants of bisection differ by the selection of refinement edges: Theiterative bisection strategy by B¨ansch [4, 6] needs the intermediate handling of hang-ing nodes. The recursive bisection strategies of [36, 38, 57] do not create such hanging elements of T . For the local multigrid algorithm examined in this article the implementation must provideaccess to the virtual refinement hierarchy. This entails suitable bookkeeping data structures, whichare available in the ALBERTA package used for the numerical experiments in Sect. 7ocal Multigrid in H ( curl ) M , the two recursivealgorithms result in exactly the same tetrahedral meshes as the iterative algorithm.Since our implementation relies on the bisection algorithm of [36], we outline its bi-section policy in the following. For more information on bisection algorithms, we referto [49, 56].For the recursive bisection algorithm of [36], the bisections of tetrahedra are to-tally determined by the local vertex numbering of M , plus a prescribed type forevery element in M . Each tetrahedron K is endowed with the local indices 0, 1, 2,and 3 for its vertices. The refinement edge of each element is always set to be the edgeconnecting vertex 0 and vertex 1. After bisection of K , the “child tetrahedron” of K which contains vertex 0 of K is denoted by Child[0] and the other one is denoted byChild[1]. The types of Child[0] and Child[1] are defined by t ype(Child[0]) = t ype(Child[1]) = ( t ype( K ) + 1) mod 3 . The new vertex at the midpoint of the refinement edge of K is always numberedby 3 in Child[0] and Child[1]. The four vertices of K are numbered in Child[0] andChild[1] as follows (see Fig. 3.2):In Child[0] : (0 , , → (0 , , , In Child[1] : (0 , , → (0 , ,
1) , if t ype( K ) = 0 , In Child[1] : (0 , , → (0 , ,
2) , if t ype( K ) > . This recursive bisection creates only a small number of similarity classes of tetrahedra,see [36, 49, 57].Fig 3.3 shows a 2D example of the recursive bisection refinement (the algorithmfor 2D case is called “the newest vertex bisection” in [41]). Similar to the 3D algorithm,for any element K , its three vertices are locally numbered by 0, 1, and 2, its refinementedge is the edge between vertex 0 and 1. The newly created vertex in the two childrenof K are numbered by 2. In the child element containing vertex 0 of K , vertex 0 and2 of K are renumbered by 1 and 0 respectively. In the other child element, vertex 1and 2 of K are renumbered by 0 and 1 respectively.In order to keep the mesh conforming during refinements, the bisection of an edgeis only allowed when such an edge is the refinement edge for all elements which sharethis edge. If a tetrahedron has to be refined, we have to loop around its refinementedge and collect all elements at this edge to create an refinement patch. Then thispatch is refined by bisecting the common refinement edge. A more detailed discussioncan be found in [36].For any mesh M l an associated “quasi-uniform” mesh c M l according to (3.2), M l ≺ c M l , is obtained as follows: the elements in { K ∈ M l : ℓ ( K ) < l } undergobisection until ℓ ( K ) = l for any K ∈ c M l .We still have to make sure that the recursive bisection allows the definition ofa virtual refinement hierarchy. Thus, let M h = M L be generated from the initialmesh M by the bisection algorithm in [36]. Denote by M hier the set of all tetrahedracreated during the bisection process, i.e., for any K ∈ M hier , there is a K ′ ∈ M h such that either K ′ = K or K ′ is created by refining K . Then, the virtual meshes M l , 0 < l < L can again be defined according to (3.4).In the following, we are going to prove that each M l is a conforming mesh, thatis, no hanging nodes occur in M l , 0 ≤ l ≤ L . The proof depends on some mild4 R. Hiptmair and W.-Y. Zheng
1 0 3 2 1 3 2 0 3 1 0 2 0 0 0 0 1 1 1 2 2 2 2 3 3 3 3
Type 0 Type 1 Type 2 K K =Child[0]of K K =Child[1]of K K =Child[0] of K K =Child[1] of K K =Child[1] of K K =Child[0] of K
1 0 2 3 2 0 1 0 1 3 2 0 3 2 1 0 3 1 2 1 1 3 3 2 0 2 0 3 3 0 1 2
Type 0
Children[0] of K Children[1] of K Children[0] of K Children[1] of K Children[1] of K Children[1] of K Children[0] of K Children[0] of K Fig. 3.2 . Bisection of tetrahedra in the course of recursive bisection. Assignment of types tochildren assumptions on M (see assumptions (A1) and (A2) in [36]), which will be taken forgranted. Lemma 3.2. [36, Lemmas 2,3] Let
T, T ′ ∈ M h be a pair of tetrahedra sharing aface F = K ∩ K ′ . It holds true that1. if T contains the refinement edge of T ′ and vice versa, then they have thesame refinement edge,2. if F contains the refinement edges of both K and K ′ , then ℓ ( K ) = ℓ ( K ′ ) ,3. if F contains the refinement edge of K , but does not contain the refinementedge of K ′ , then ℓ ( K ) = ℓ ( K ′ ) + 1 ,4. if F does not contain the refinement edges of K and K ′ , then ℓ ( K ) = ℓ ( K ′ ) . Lemma 3.3.
The meshes M l , ≤ l ≤ L , according to (3.4) are conformingmeshes.Proof . We are going to prove the lemma by backward induction starting from l = L . Since M L = M h is conforming, for any K ∈ M L satisfying ℓ ( K ) = L , thereexists a brother of K , denoted by K ′ ∈ M L , such that ℓ ( K ′ ) = L and K p := K ∪ K ′ ∈M L − . Here K p is called the parent of K and K ′ with ℓ ( K p ) = L − E be the refinement edge of K p . By the recursive bisection algorithm, E must ocal Multigrid in H ( curl ) M M M M M M Fig. 3.3 . Virtual refinement hierarchy for 2D triangular meshes emerging in the course ofsuccessive local newest vertex bisection refinement of M from Fig. 3.1. Accompanying quasi-uniformmeshes outlined in blue, maximally refined triangles marked yellow. be the common refinement edge of all tetrahedra in the refinement patch: P E = [ { K ′ p : K ′ p ∈ M L − and E ⊂ K ′ p } . By Lemma 3.2, ℓ ( K ′ p ) = L − K ′ p ⊂ P E and the midpoint of E , denoted by A new , is the unique new vertex of M L in P E . We conclude that P E = [ { K : K ∈ M L , ℓ ( K ) = L, and A new is a vertex of K } . Coarsen the sub-mesh M L | P E by removing the vertex A new and all edges related to itand adding E to this patch. Thus a conforming sub-mesh M L − | P E is obtained. Dothe above coarsening process for every element K ∈ M L with ℓ ( K ) = L . This provesthat M L − is conforming.Finally, an induction argument confirms that M l is conforming, l = L − , · · · ,
4. Local multigrid.
To begin with, we introduce nested refinement zones asopen subsets of Ω: ω l := interior (cid:16) [ { K : K ∈ M h , ℓ ( K ) ≥ l } (cid:17) ⊂ Ω , (4.1)see Fig. 4.1 and Fig. 4.2. The notion of refinement zones allows a concise definition ofthe local multilevel decompositions of the finite element spaces V ( M h ) and U ( M h )that underly the local multigrid method.6 R. Hiptmair and W.-Y. Zheng A new
0 0 1 1
E E K K’ K p K’ p Fig. 3.4 . The patch around a refinement edge E with vertex 0 and 1. ℓ ( K ) = ℓ ( K ′ ) = L and ℓ ( K p ) = ℓ ( K ′ p ) = L − . “Refinement strips”: set differ-ences of refinement zones: Σ := ω \ ω : Σ := ω \ ω : Σ := ω \ ω : Σ := ω Fig. 4.1 . Refinement zones for the 2D refinement hierarchy of Figure 3.1.
We introduce local multigrid from the perspective of multilevel successive sub-space correction (SSC) [61, 62, 64]. First, we give an abstract description for a linearvariational problem u ∈ H : a ( u, v ) = f ( v ) ∀ v ∈ H , (4.2)involving a positive definite bilinear form a on a Hilbert space H . The method is ocal Multigrid in H ( curl ) (cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1) (cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1) “Refinement strips”: set differ-ences of refinement zones: Σ := ω \ ω : Σ := ω \ ω : Σ := ω \ ω : Σ := ω \ ω (cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1) : Σ := ω \ ω (cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1) : Σ := ω \ ω (cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1)(cid:1) : Σ := ω Fig. 4.2 . Refinement zones for the 2D refinement hierarchy of Figure 3.3. completely defined after we have provided a finite subspace decomposition H = J X j =0 H j , H j ⊂ H closed subspaces , j = 0 , . . . , J, J ∈ N . (4.3)Then the correction scheme implementation of one step of SSC acting on the iterate u m − reads: for m = 1 , , · · · u m − − = u m − for j = 0 , , · · · , J Let e j ∈ H j solvea ( e j , v j ) = f ( v j ) − a ( u m − j − , v j ) ∀ v j ∈ H j u m − j = u m − j − + e j endfor u m = u m − J endfor This amounts to a stationary linear iterative method with error propagation op-erator E = ( I − P J )( I − P J − ) · · · ( I − P ) , (4.4)where P j : H H j stands for the Galerkin projection defined through(4.5) a ( P j v, v j ) = a ( v, v j ) ∀ v j ∈ H j . The convergence theory of SSC for an inner product a and induced energy norm k·k A rests on two assumptions. The first one concerns the stability of the space de-composition . We assume that there exists a constant C stab independent of J such8 R. Hiptmair and W.-Y. Zheng that(4.6) inf n J X j =0 k v j k A : J X j =0 v j = v o ≤ C stab k v k A ∀ v ∈ H. The second assumption is a strengthened Cauchy-Schwartz inequality , namely, thereexist two constants 0 ≤ q < C orth independent of j and k such that(4.7) a ( v j , v k ) ≤ C orth q | k − j | k v j k A k v k k A ∀ v j ∈ H j , v k ∈ H k . The above inequality states a kind of quasi-orthogonality between the subspaces.From [61, Theorem 4.4] and [66, Theorem 5.1] we cite the following central convergencetheorem:
Theorem 4.1.
Provided that (4.6) and (4.7) hold, the convergence rate of Algo-rithm SSC is bounded by (4.8) k E k A ≤ − C stab (1 + Θ) with Θ = C orth q − q , where the operator norm is defined by k E k A := sup v ∈ H,v =0 k Ev k A k v k A . The bottom line is that the subspace splitting (4.3) already provides a full de-scription of the method. Showing that both constants C stab from (4.6) and C orth from(4.7) can be chosen independently of the number L of refinement levels is the challengein asymptotic multigrid analysis.In concrete terms, the role of the linear variational problem (4.2) is played by(1.1) considered on the edge element space U ( M h ), which replaces the Hilbert space H . To define the local multilevel decomposition of U ( M h ), we define “sets of newbasis functions” on the various refinement levels B V := B V ( M ) , B lV := { b h ∈ B V ( M l ) : supp b h ⊂ ω l } , B U := B U ( M ) , B l U := { b h ∈ B U ( M l ) : supp b h ⊂ ω l } , ≤ l ≤ L . (4.9)A 2D drawing of the sets B lV is given in Fig. 4.3 where Γ D = ∂ Ω. Note that wealso have to deal with V ( M h ), because, as suggested by the reasoning in [29], a localmultilevel decomposition of U ( M h ) has to incorporate an appropriate local multileveldecomposition of V ( M h ).Then, a possible local multigrid iteration for the linear system of equations arisingfrom a finite element Galerkin discretization of a H D (Ω)-elliptic variational problemboils down to a successive subspace correction method based on the local multileveldecomposition V ( M h ) = V ( M ) + L X l =1 X b h ∈ B lV Span { b h } . (4.10)Similarly, the local multilevel splitting of U ( M h ) is based on the multilevel de-composition U ( M h ) = U ( M ) + L X l =1 X b h ∈ B lV Span { grad b h } + L X l =1 X b h ∈ B l U Span { b h } . (4.11) ocal Multigrid in H ( curl ) l = 0 l = 1 l = 2 l = 3 Fig. 4.3 . Active vertices (red) carrying “tent functions” in B lV , Γ D = ∂ Ω , refinement hierarchyof Fig. 3.1 These splittings induce SSC iterations that can be implemented as non-symmetricmultigrid V-cycles with only one (hybrid) Gauss-Seidel post-smoothing step, see [29,Sect. 6]. Duplicating components of (4.11) results in more general multigrid cycleswith various numbers of pre- and post-smoothing steps.The splitting (4.11) is motivated both by the design of multigrid methods for(1.1) and U ( M ) in the case of uniform refinement and local multigrid approachesto H D (Ω)-elliptic variational problems after discretization by means of linear finiteelements [41,60]. The occurrence of gradients of “tent functions” b h in (4.11) is relatedto the hybrid local relaxation , which is essential for the performance of multigridin H ( curl , Ω), see [29] for a rationale. A rigorous justification will emerge duringthe theoretical analysis in the following sections. It will establish the following maintheorem.
Theorem 4.2 (
Asymptotic convergence of local multigrid for edge elements ). Under the assumptions on the meshes made above and allowing at most one hangingnode per edge, the decomposition (4.11) leads to an SSC iteration whose convergencerate is bounded away from uniformly in the number L of refinement steps.
5. Stability.
First we tackle the stability estimate (4.6) for the local multileveldecomposition (4.10), which is implicitly contained in (4.11). V ( M ) . Quasi-interpolation operatorsare projectors onto finite element spaces that have been devised to accommodate twoconflicting goals: locality and boundedness in weak norms [18, 21, 46, 50, 53]. As keytool they will be used in Sect. 5.2 and the proof of Lemma 5.10. As in [46, Sect. 2.1.1],we resort to a construction employing local linear L -dual basis functions. We followthe analysis of [50] that permits us to take into account Dirichlet boundary conditions.For a generic tetrahedron K define ψ Kj , j = 1 , , ,
4, by L ( K )-duality to thebarycentric coordinate functions λ i , i = 1 , , ,
4, of K : ψ Kj ∈ P ( K ) : Z K ψ Kj ( x ) λ i ( x ) d x = δ ij , i, j ∈ { , . . . , } . (5.1)Computing an explicit representation of the ψ Kj we find C − ≤ | K | (cid:13)(cid:13) ψ Kj (cid:13)(cid:13) L ( K ) ≤ C , C − ≤ (cid:13)(cid:13) ψ Kj (cid:13)(cid:13) L ( K ) ≤ C , (5.2)with an absolute constant
C >
0. We can regard ψ Kj as belonging to the j -th vertexof K . Thus, we will also write ψ K p , p ∈ N ( K ), N ( K ) the set of vetices of K .0 R. Hiptmair and W.-Y. Zheng
Let M be one of the tetrahedral meshes M l or c M l of Ω. In order to introducequasi-interpolation operators we take for granted some “node → cell”–assignment, amapping N ( M )
7→ M , p ∈ N ( M ) K p ∈ M . Definition 5.1.
Writing { b p } p ∈N ( M ) := B V ( M ) , define the local quasi-interpolation operator (5.3) Q h : ( L (Ω) V ( M ) u P p ∈N ( M ) R K p ψ K p p ( x ) u ( x ) d x · b p . Analoguously, we introduce the local quasi-interpolation Q h : L (Ω) V ( M ) . We point out that Q h respects u = 0 on Γ D , because the sum does not cover basisfunctions attached to vertices on Γ D . From (5.1) it is also evident that both Q h and Q h are projections, for instance,(5.4) Q h u h = u h ∀ u h ∈ V ( M ) . Moreover, they satisfy the following strong continuity and approximation properties:
Lemma 5.2.
The quasi-interpolation operators from Def. 5.1 allow the estimates(set Γ D = ∅ for Q M ) ∃ C = C ( ρ M ) : k Q h u k L (Ω) ≤ C k u k L (Ω) ∀ u ∈ L (Ω) , (5.5) ∃ C = C ( ρ M , Ω , Γ D ) : | Q h u | H (Ω) ≤ C | u | H (Ω) ∀ u ∈ H D (Ω) , (5.6) ∃ C = C ( ρ M , k ) : (cid:13)(cid:13) h − k ( u − Q h u ) (cid:13)(cid:13) L (Ω) ≤ C | u | H k (Ω) ∀ u ∈ H k (Ω) ∩ H D (Ω) , (5.7) and k = 1 , .Proof . [Part I] Continuity in L (Ω) is a simple consequence of the stability (2.9)of the nodal bases B V ( M ) and of the Cauchy-Schwarz inequality: k Q h u k L (Ω) ≤ C X p ∈N ( M ) | Q h u ( p ) | k b p k L (Ω) = C X p ∈N ( M ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)Z K p ψ K p p ( x ) u ( x ) d x (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) k b p k L (Ω) ≤ C X p ∈N ( M ) (cid:13)(cid:13)(cid:13) ψ K p p (cid:13)(cid:13)(cid:13) L ( K p ) k b p k L (Ω) k u k L ( K p ) ≤ C k u k L (Ω) , with C = C ( ρ M ) >
0, because (cid:13)(cid:13)(cid:13) ψ K p p (cid:13)(cid:13)(cid:13) L ( K p ) k b p k L (Ω) ≤ C , too.The following estimate is instrumental in establishing continuity of Q h in H D (Ω): Theorem 5.3 (
Generalized Hardy inquality ). ∃ C = C (Ω , Γ D ) > Z Ω (cid:12)(cid:12)(cid:12)(cid:12) u dist( x , Γ D ) (cid:12)(cid:12)(cid:12)(cid:12) d x ≤ C | u | H (Ω) ∀ u ∈ H D (Ω) . ocal Multigrid in H ( curl ) Proof . By density it suffices to consider u ∈ C ∞ (Ω), supp( u ) ∩ Γ D = ∅ . Using apartition of unity, we can confine the estimate to neighborhoods of Γ D , in which ∂ Ω isthe graph of a Lipschitz-continuous function. Thus, after bi-Lipschitz transformations,we need only investigate three canonical situations, see Fig. 5.1:1. Γ D = { z = 0 } , for which the 1D Hardy inequality gives the estimate, see theproof of Thm. 1.4.4.4 in [27].2. Γ D = { z = 0 ∧ x > } , which can be treated using polar coordinates in the( x, z )-plane and then integrating in y -direction: ∞ Z π Z (cid:12)(cid:12)(cid:12)(cid:12) u ( r, ϕ ) r (cid:12)(cid:12)(cid:12)(cid:12) d ϕ r d r ≤ ∞ Z π Z (cid:12)(cid:12)(cid:12)(cid:12) πr ∂u∂ϕ ( r, ϕ ) (cid:12)(cid:12)(cid:12)(cid:12) d ϕ r d r ≤ π Z z> | grad x,z u | d x d z .
3. Γ D = { z = 0 ∧ x > ∧ y > } , for which we obtain a similar estimate usingspherical coordinates.This ends the proof.PSfrag replacements xyz Γ D ∂ Ω \ Γ D Γ D = { z = 0 } PSfrag replacements xyz Γ D ∂ Ω \ Γ D ϕr Γ D = { z = 0 ∧ x > } PSfrag replacements xyz Γ D ∂ Ω \ Γ D ∂ Ω \ Γ D Γ D = { z = 0 ∧ x > ∧ y > } Fig. 5.1 . Canonical situations to be examined in the proof of Thm. 5.3
Proof . [of Lemma 5.2, part II] In order to tackle the H (Ω)-continuity of Q h , weuse that grad V ( M ) ⊂ U ( M ) along with the stability estimate (2.6) k grad Q h u k L (Ω) ≤ C X E =[ p , q ] ∈E ( M ) ( Q h u ( p ) − Q h u ( q )) k b E k L (Ω) , (5.8)with the notation { b E } E ∈E ( M ) := B U ( M ).(i) for the case E = [ p , q ] ∈ E ( M ), p , q Γ D , we adapt arguments from [50].For any u ∈ H D (Ω), by (5.1), we have the identity | ( Q h u )( p ) − ( Q h u )( q ) | = (cid:12)(cid:12)(cid:12)Z K p Z K q ψ K p p ( x ) ψ K q q ( y )( u ( x ) − u ( y )) d y d x (cid:12)(cid:12)(cid:12) = (cid:12)(cid:12)(cid:12)Z K p Z K q ψ K p p ( x ) ψ K q q ( y ) Z grad u ( y + τ ( x − y )) · ( x − y ) d τ d y d x (cid:12)(cid:12)(cid:12) . Then split the innermost integral and transform Z f ( y + τ ( x − y )) d τ = Z f ( y + τ ( x − y )) d τ + Z f ( x + τ ( y − x )) d τ . R. Hiptmair and W.-Y. Zheng
We infer | ( Q h u )( p ) − ( Q h u )( q ) |≤ Z Z K p Z K q | ψ K p p ( x ) || ψ K q q ( y ) | | grad u ( y + τ ( x − y )) || x − y | d y d x d τ + Z Z K p Z K q | ψ K p p ( x ) || ψ K q q ( y ) | | grad u ( x + τ ( y − x )) | | x − y | d y d x d τ The transformation formula for integrals reveals Z K f ( x + τ ( y − x )) d y = τ − Z K ′ f ( z ) d z , K ′ := x + τ ( K − x ) . Appealing to the bounds for (cid:13)(cid:13) ψ Kj (cid:13)(cid:13) L ( K ) , (cid:13)(cid:13) ψ Kj (cid:13)(cid:13) L ( K ) , K ∈ M , from (5.2), the Cauchy-Schwarz inequality yields | ( Q h u )( p ) − ( Q h u )( q ) | k b E k L (Ω) (5.9) ≤ C " | p − q | k b E k L (Ω) min {| K q | , | K p | } ≤ C = C ( ρ M ) Z τ − / d τ · | u | H ( h Ω E i ) . Here h Ω E i stands for the convex hull of all tetrahedra adjacent to the edge E .(ii) now consider E = [ p , q ] ∈ E ( M ), p ∈ Γ D . Then, for any u ∈ H D (Ω) | ( Q h u )( q ) − ( Q h u )( p ) | = | ( Q h u )( q ) | = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)Z K q ψ K q q ( x ) u ( x ) d x (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) = (cid:12)(cid:12)(cid:12)Z K q dist( x , Γ D ) | {z } ≤ C | p − q | ψ K q q ( x ) u ( x )dist( x , Γ D ) d x (cid:12)(cid:12)(cid:12) ≤ C | p − q | | K q | · Z K q (cid:12)(cid:12)(cid:12)(cid:12) u ( x )dist( x , Γ D ) (cid:12)(cid:12)(cid:12)(cid:12) d x ≤ C | p − q | Z Ω E (cid:12)(cid:12)(cid:12)(cid:12) u ( x )dist( x , Γ D ) (cid:12)(cid:12)(cid:12)(cid:12) d x , with (different) constants C = C ( ρ M ) > M in the form ∃ C = C ( ρ M ) : ♯ { E ∈ E ( M ) : x ∈ h Ω E i} ≤ C ∀ x ∈ Ω , and appealing to Thm. 5.3 confirm | Q h u | H (Ω) ≤ C | u | H (Ω) . Observe that the Hardyinequality makes the constant depend on Ω and Γ D in addition.The quasi-interpolation error estimate (5.7) results from scaling arguments. Pick K ∈ M , u ∈ H (Ω) ∩ H D (Ω), and write I K u ∈ P ( K ) for the linear interpolant of u on K . Thanks to the projection property, we deduce as in Part I of the proof that,with C = C ( ρ M ), k ( Id − Q h ) u k L ( K ) = k ( Id − Q h )( u − I K u ) k L ( K ) ≤ C k u − I K u k L (Ω K ) ≤ Ch K | u | H (Ω K ) . ocal Multigrid in H ( curl ) K := S { K ′ : K ′ ∩ K = ∅} , and the final estimate can be shownby a simple scaling argument, cf. (2.10). Estimate (5.7) for k = 1 follows by scalingarguments and interpolation between the Sobolev spaces H (Ω K ) and L (Ω K ). V ( M L ) . In this section we first revisit the well-known [45, 61, 67] uniform stability of multilevel splittings of H (Ω)-conforming La-grangian finite element functions in the case of mesh hierarchies generated by uniform, i.e. non-local, regular refinement.We take for granted a virtual refinement hierarchy (3.1) of tetrahedral meshes asintroduced in Sect. 3 and its accompanying quasi-uniform family of meshes (3.2).Owing to the inf in (4.6), it is enough to find a concrete family of admissible“candidate” decompositions that enjoys the desired L -uniform stability. We aim forcandidates that fit the locally refined mesh hierarchy.The principal idea, borrowed from [46, Sect. 4.2.2], is to use a sequenceof quasi-interpolation operators Q l : L (Ω) V ( c M l ) based on a judiciouslychosen node → element–assignments. For c M l we introduce a “coarsest neighbornode → element–assignment”: First, for any p ∈ N ( c M l ), l = 1 , . . . , L , we pick K ∈ M l such that ℓ ( K ) = min { ℓ ( K ) : p ∈ K, K ∈ M l } . Secondly, we select a “coarsest neighbor” K p ∈ c M l among those elements of c M l thatare contained in K . This defines a mapping N ( c M l ) d M l , p K p . We write Q l : L (Ω) V ( c M l ) for the induced quasi-interpolation operator according to Def. 5.1.Next, we examine the candidate multilevel splitting u h = Q u h + L X l =1 ( Q l − Q l − ) u h , u h ∈ V ( c M L ) . (5.10) Lemma 5.4.
There holds, with a constant
C > depending only on Ω and theuniform bound for the shape regularity measures ρ c M l , ≤ l ≤ L , (5.11) (cid:12)(cid:12) Q u h (cid:12)(cid:12) H (Ω) + L X l =1 h − l (cid:13)(cid:13) ( Q l − Q l − ) u h (cid:13)(cid:13) L (Ω) ≤ C | u h | H (Ω) ∀ u h ∈ V ( c M L ) . Proof . We take the cue from the elegant approach of Bornemann and Yserentantin [11], who discovered how to bring techniques of real interpolation theory of Sobolevspaces [37], [39, Appendix B] to bear on (5.10). The main tools are the so-called K -functionals given by K ( t, u ) := inf w ∈ H (Ω) n k u − w k L (Ω) + t | w | H (Ω) o ,K R ( t, u ) := inf w ∈ H ( R ) n k u − w k L ( R ) + t | w | H ( R ) o . The estimates (5.6) and (5.7) of Lemma 5.2 create a link between the terms in (5.11)and K ( t, u ): owing to (5.5) and (5.7) there holds for any u ∈ L (Ω) (cid:13)(cid:13) ( Q l − Q l − ) u (cid:13)(cid:13) L (Ω) ≤ (cid:13)(cid:13) ( Q l − Q l − )( u − w ) (cid:13)(cid:13) L (Ω) + (cid:13)(cid:13) ( Q l − Q l − ) w (cid:13)(cid:13) L (Ω) ≤ C (cid:0) k u − w k L (Ω) + h l | w | H (Ω) (cid:1) ∀ w ∈ H (Ω) . R. Hiptmair and W.-Y. Zheng
Here and below the generic constants C may depend on shape regularity max ≤ l ≤ L ρ c M l and the (quasi-uniformity) constants in (3.3). We conclude (cid:13)(cid:13) ( Q l − Q l − ) u (cid:13)(cid:13) L (Ω) ≤ C K ( h l , u ) ∀ u ∈ L (Ω) , (5.12)which implies(5.13) (cid:12)(cid:12) Q u (cid:12)(cid:12) H (Ω) + L X l =1 h − l (cid:13)(cid:13) ( Q l − Q l − ) u (cid:13)(cid:13) L (Ω) ≤ C n | u | H (Ω) + L X l =1 h − l K ( h l , u ) o . Let e u ∈ H ( R ) be the Sobolev extension of u such that, with C = C (Ω) > e u | Ω = u and | u | H ( R ) ≤ C | u | H (Ω) . Define the Fourier Transform of e u by b u ( ξ ) = 1(2 π ) / Z R e u ( x ) e − ı x · ξ d x . By the equivalent definition of Sobolev-norms on R | e u | H i ( R ) ≈ Z R | ξ | i | b u ( ξ ) | d ξ , i = 0 , , we have K R ( t, e u ) ≤ C inf w ∈ H ( R ) Z R n | b u ( ξ ) − b w ( ξ ) | + t | ξ | | b w | o d ξ = C Z R t | ξ | t | ξ | | b u ( ξ ) | d ξ , because the infimum is attained for [39, Thm. B7] b w ( ξ ) = b u ( ξ ) / (1 + t | ξ | ). Since K ( t, u ) = inf w ∈ H (Ω) n k u − w k L (Ω) + t | w | H (Ω) o , = inf w ∈ H ( R ) n k u − w k L (Ω) + t | w | H (Ω) o ≤ K R ( t, e u ) , we deduce that L X l =1 h − l K ( h l , u ) ≤ C L X l =1 Z R h l | ξ | h l | ξ | | b u ( ξ ) | d ξ (5.14) ≤ C sup ξ ∈ R n L X l =1 θ l | ξ | θ l | ξ | o Z R | ξ | | b u ( ξ ) | d ξ ≤ C | b u | H ( R ) ≤ C | u | H (Ω) , where we have used assumption (3.3). The proof is finished by combining (5.13) and(5.14).Now we restrict ourselves to u h ∈ V ( M h ). Then, thanks to the particular design ofthe node → element–assignment underlying Q l , the terms in the decomposition (5.10)turn out to be localized. ocal Multigrid in H ( curl ) Lemma 5.5.
For all u h ∈ V ( M h ) and ≤ l ≤ j ≤ L , Q j u h = u h in Ω \ ω l +1 . (5.15) Proof . If p ∈ N ( c M j ) and p ω l +1 (open set !), then K p ω l +1 ( K p ∈ c M j ).Recall that K p was deliberately chosen such that there is K ∈ M l with K p ⊂ K .Since u h is linear on K , the same holds for K p and (5.1) guarantees( Q j u h )( p ) = u h ( p ) . When restricted to Ω \ ω l +1 , the mesh c M j is a refinement of M h . Hence, agreementof the M h -piecewise linear function u h with Q j u h in all nodes of c M j outside ω l +1 implies Q j u h | Ω \ ω l +1 = u h | Ω \ ω l +1 .Consequently, for any u h ∈ V ( M h ), outside ω l both Q l u h and Q l − u h agree with u h . Corollary 5.6.
For any u h ∈ V ( M h ) and ≤ l ≤ L , supp(( Q l − Q l − ) u h ) ⊂ ω l . In other words, the components of (5.10) are localized inside refined regions of Ω.In light of the definition (4.1) of the refinement zones, we also find( Q l − Q l − ) u h ∈ V ( M l ) !(5.16)However, having used Q l we cannot expect the splitting to match potentialhomogeneous Dirichlet boundary conditions. This can be remedied using Oswald’strick [44, Cor. 30]. We fix u h ∈ V ( M h ) and abbreviate u = Q u h ∈ V ( M ), u l := ( Q l − Q l − ) u h ∈ V ( M l ), l ≥
1. Then, we consider the partial sums s l := l X j =0 u j ∈ V ( M l ) l ≥ . (5.17)Dropping those basis functions in B V ( M l ) that belong to vertices in Γ D in the rep-resentation of s l we arrive at s l ∈ V ( M l ) ∈ H D (Ω).Due to Cor. 5.6, we observe that s l and s l − agree on Ω \ ω l .(5.18)Hence, away from ω l ∩ Γ D the same basis contribution are removed from both functionswhen building s l and s l − , respectively. This permits us to conclude s l and s l − agree on Ω \ ω l .(5.19)Putting it differently, supp( s l − s l − ) ⊂ ω l . (5.20)6 R. Hiptmair and W.-Y. Zheng
Hence, for all 1 ≤ l ≤ L we can estimate k s l − s l − k L (Ω) = k s l − s l − k L ( ω l ) ≤ k s l − s l k L ( ω l ) + k s l − − s l − k L ( ω l ) + k u l k L ( ω l ) . (5.21)The benefit of zeroing in on ω l is that on this subdomain s l has the same “uniformscale” h l as u l . Thus, repeated application of uniform L -stability estimates (2.9) forbasis representations and elementary Cauchy-Schwarz inequalities make possible theestimates (for arbitrary 0 < ǫ < ) k s l − s l k L ( ω l ) ≤ Ch l X p ∈N (Γ l ) s l ( p ) ≤ Ch l (cid:13)(cid:13) s l | ∂ Ω (cid:13)(cid:13) L (Γ l ) = Ch l (cid:13)(cid:13)(cid:13)(cid:13)(cid:13)(cid:13) L X j =0 u j − l X j =0 u j (cid:13)(cid:13)(cid:13)(cid:13)(cid:13)(cid:13) L (Γ l ) ≤ Ch l (cid:16) L X j = l +1 k u j k L (Γ l ) (cid:17) ≤ Ch l (cid:16) L X j = l +1 h − j k u j k L ( ω l ) (cid:17) ≤ Ch l · L X j = l +1 h − ǫj · L X j = l +1 h ǫ − j k u j k L ( ω l ) ≤ Ch − ǫl · L X j = l +1 h ǫ − j k u j k L ( ω l ) . Here the set N (Γ l ) comprises the nodes of N ( c M l ) that lie on ω l ∩ Γ D and we makeheavy use of the geometric decay of h l . The latter also yields L X l =1 h − l k s l − s l k L ( ω l ) ≤ C L X l =1 h − ǫl L X j = l +1 h ǫ − j k u j k L ( ω l ) = C L X j =2 j X l =1 h − ǫl ! h − ǫj k u j k L ( ω l ) ≤ C L X j =2 h − j k u j k L (Ω) ≤ C | u h | H (Ω) , by virtue of Lemma 5.4. Except for the last line, all constants depend only on ρ c M l and the constants in (3.3). Merging the last estimate with (5.21) gives us L X l =1 h − l k s l − s l − k L (Ω) ≤ C | u h | H (Ω) . (5.22)Thus, in light of (5.20) and the following identity s + L X l =1 ( s l − s l − ) = s L = s L = u h , we have accomplished the proof of the following theorem: Theorem 5.7.
For any u h ∈ V ( M h ) we can find u l ∈ V ( M l ) such that u h = L X l =0 u l , supp( u l ) ⊂ ω l , (5.23) ocal Multigrid in H ( curl ) and | u | H (Ω) + L X l =1 h − l k u l k L (Ω) ≤ C | u h | H (Ω) , with C > independent of L . Notice that in combination with the L -stability (2.9) of nodal bases and inverseinequalities, this theorem asserts an L -uniform estimate of the form (4.6) for thesplitting (4.10) w.r.t. the energy norm |·| H (Ω) . From (5.23) it is clear that the basisfunctions admitted in (4.10) can represent the functions u l of Thm. 5.7. Remark 5.8.
It is interesting to note that, in contrast to other analyses [1, 11],the above proof does not hinge on Assumption 3.0.1. Thm. 5.7 remains valid for anarbitrary number of hanging nodes on an active edge. Howver, this does not translateinto asymptotically optimal convergence of local H -multigrid in this case, because, inorder to infer it from Thm. 5.7, we also need uniform L -stability of the bases. Helmholtz-type decompositions, alsocalled regular decompositions , have emerged as a powerful tool for answering questionsconnected with H ( curl , Ω). In particular, they have paved the way for a rigorousmultigrid theory for H ( curl , Ω)-elliptic problems [17, 25, 29, 31–33, 35, 47]. We referto [30, Sect. 2.4] for more information.We will need a very general version provided by the following theorem.
Theorem 5.9.
Let Ω meet the requirements stated in Sect. 1. Then, for any v ∈ H Γ D ( curl , Ω) , there exists a p ∈ H D (Ω) and Ψ ∈ ( H D (Ω)) such that v = ∇ p + Ψ , (5.24) | p | H (Ω) ≤ C k v k H ( curl , Ω) , k Ψ k H (Ω) ≤ C k curl v k L (Ω) , (5.25) where the constant C depends only on Ω . PSfrag replacements ΩΩ Ω Ω Fig. 5.2 . Buffer zones attached to connected components of (red) Dirichlet boundary part Γ D R. Hiptmair and W.-Y. Zheng
Proof . Given u ∈ H Γ D ( curl , Ω), we define e u ∈ H ( curl , e Ω), e Ω := interior(Ω ∪ Ω ∪ Ω ∪ . . . ) (see Sect. 1 and Fig. 5.2 for the meaning of Ω i ), by e u ( x ) = ( u ( x ) for x ∈ Ω , x ∈ Ω i for some i . (5.26)Notice that the tangential components of e u are continuous across ∂ Ω, which ensures e u ∈ H ( curl , e Ω). Then extend e u to u ∈ H ( curl , R ), see [16].Since curl u ∈ H (div 0 , R ), Fourier techniques [24, Sect. 3.3] yield a Φ ∈ ( H ( R )) that fulfills curl Φ = curl u , k Φ k H ( R ) ≤ C k curl u k L ( R ) , (5.27)with C = C (Ω) >
0. As a consequence curl ( u − Φ ) = 0 ⇒ u − Φ = grad q in R . (5.28)On every Ω i , by definition u = 0, which implies q | Ω i ∈ H (Ω i ). As the attacheddomains Ω i are well separated Lipschitz domains, see Fig. 5.2, the H -extension of q | S i Ω i to q ∈ H ( R ) is possible. Moreover, it satisfies k q k H ( R ) ≤ C k q k H ( S i Ω i ) ≤ C k Φ k H ( R ) ≤ k curl u k L (Ω) . (5.29) u = Φ − grad q + grad ( q + q ) . (5.30)Finally, set Ψ := ( Φ − grad q ) | Ω , p := q + q , and observe k Ψ k H (Ω) ≤ k Φ k H ( R ) + k q k H ( R ) ≤ C k curl u k L (Ω) . (5.31)The constants may depend on Ω, Γ D , and the chosen Ω i .The stable Helmholtz-type decomposition (5.24) immediately suggests the follow-ing idea: when given v h ∈ U ( M h ), first split it according to (5.24) and then attackboth components by the uniformly H -stable local multilevel decompositions exploredin the previous section. Alas, the idea is flawed, because neither of the terms in (5.24)is guaranteed to be a finite element function, even if this holds for v h .Fortunately, the idea can be mended by building a purely discrete counterpart of(5.24) as in [33, Lemma 5.1] (called there “discrete regular decomposition”). For thesake of completeness we elaborate the proof below. Lemma 5.10.
For any v h ∈ U ( M h ) , there is Ψ h ∈ ( V ( M h )) , p h ∈ V ( M h ) ,and e v h ∈ U ( M h ) such that v h = e v h + Π h Ψ h + ∇ p h , (5.32) k p h k H (Ω) ≤ C k v h k H ( curl , Ω) , (5.33) (cid:13)(cid:13) h − e v h (cid:13)(cid:13) L (Ω) + k Ψ h k H (Ω) ≤ C k curl v h k L (Ω) , (5.34) where the constant C depends only on Ω , Γ D , and the shape regularity of M h .Proof . ( cf. [33, Lemma 5.1]) We fix a v h ∈ U ( M h ) and use the stable regulardecomposition of Thm. 5.9 to split it according to v h = Ψ + grad p , Ψ ∈ ( H D (Ω)) , p ∈ H D (Ω) . (5.35) ocal Multigrid in H ( curl ) Ψ and p satisfy k Ψ k H (Ω) ≤ C k curl v h k L (Ω) , k grad p k L (Ω) ≤ C k v h k H ( curl , Ω) , (5.36)with constants depending only on Ω and Γ D .Next, note that in (5.35) curl Ψ = curl v h ∈ curl U ( M h ), and, owing toLemma 2.1, Π h Ψ is well defined. Further, a commuting diagram property togetherwith Lemma 2.3 implies curl ( Id − Π h ) Ψ = 0 ⇒ ∃ q ∈ H D (Ω) : ( Id − Π h ) Ψ = grad q . (5.37)The estimate of Lemma 2.1 together with (5.36) yields (cid:13)(cid:13) h − grad q (cid:13)(cid:13) L (Ω) = (cid:13)(cid:13) h − ( Id − Π h ) Ψ (cid:13)(cid:13) L (Ω) ≤ C | Ψ | H (Ω) ≤ C k curl v h k L (Ω) . (5.38)In order to push Ψ into a finite element space, a quasi-interpolation operator Q h : ( L (Ω)) ( V ( M h )) is the right tool. We simply get it from componentwiseapplication of an operator according to Def. 5.1 where any node → element–assignmentwill do. Thus, we can define the terms in the decomposition (5.32) as e v h := Π h ( Ψ − Q h Ψ ) ∈ U ( M h ) , (5.39) Ψ h := Q h Ψ ∈ ( V ( M h )) , (5.40) grad p h := grad ( p + q ) , p h ∈ V ( M h ) . (5.41)Indeed, grad ( p + q ) ∈ U ( M h ) such that p + q ∈ V ( M h ). The stability of the decom-position (5.32) can be established as follows: first, make use of Lemma 2.1 and (5.7)to obtain, with C = C ( ρ M h ) > (cid:13)(cid:13) h − e v h (cid:13)(cid:13) L (Ω) ≤ (cid:13)(cid:13) h − ( Id − Π h )( Ψ − Q h Ψ ) (cid:13)(cid:13) L (Ω) + (cid:13)(cid:13) h − ( Id − Q h ) Ψ (cid:13)(cid:13) L (Ω) ≤ C | ( Id − Q h ) Ψ | H (Ω) + | Ψ | H (Ω) ≤ C | Ψ | H (Ω) ≤ C k curl v h k L (Ω) . Due to the definition (5.40), the next estimate is a simple consequence of (5.6) andThm. 5.9 k Ψ h k H (Ω) ≤ C k Ψ k H (Ω) ≤ C k curl v h k L (Ω) . (5.42)Finally, the estimates established so far plus the triangle inequality yield k grad p h k L (Ω) ≤ C k v h k H ( curl , Ω) . (5.43) ( M h ) . With the discrete Helmholz-typedecomposition of Lemma 5.10 at our disposal, we can now tackle its piecewise linearand continuous components with Thm. 5.7.
Lemma 5.11.
For any v h ∈ U ( M h ) , there exists a constant C only depending onthe domain, the Dirichlet boundary part Γ D , the shape regularity of the meshes M l , c M l , ≤ l ≤ L , and the constants in (3.3) , such that v h = L X l =0 (cid:16) v l + ∇ p l (cid:17) , v l ∈ Span (cid:8) B l U (cid:9) , p l ∈ Span (cid:8) B lV (cid:9) , (5.44)0 R. Hiptmair and W.-Y. Zheng and (5.45) k v k H ( curl , Ω) + | p | H (Ω) + L X l =1 h − l (cid:16) k v l k L (Ω) + k p l k L (Ω) (cid:17) ≤ C k v h k H ( curl , Ω) , where B lV and B l U are defined in (4.9) .Proof . We start from the discrete Helmholtz-type decomposition of v h in (5.32): v h = e v h + Π h Ψ h + ∇ p h , Ψ h ∈ ( V ( M h )) , p h ∈ V ( M h ) , e v h ∈ U ( M h ) . We apply the result of Thm. 5.7 about the existence of stable local multilevel splittingsof V ( M h ) componentwise to Ψ h : this gives Ψ h = L X l =0 Ψ l , Ψ l ∈ Span (cid:8) B lV (cid:9) , (5.46) | Ψ | H (Ω) + L X l =1 h − l k Ψ l k L (Ω) ≤ C | Ψ h | H (Ω) . (5.47)Observe that the functions Ψ l do not belong to U ( M l ). Thus, we target them withedge element interpolation operators Π l onto U ( M l ), see (2.7), and obtain the split-ting described in Lemma 2.2:(5.48) Ψ l = Π l Ψ l + ∇ w l , w l ∈ e V ( M l ) . The gradient terms introduced by (5.48) are well under control: writing s h := P Ll =0 w l ,the L -stability of (5.48), see Lemma 2.2, yields k Π l Ψ l k L (Ω) ≤ C k Ψ l k L (Ω) , | s h | H (Ω) ≤ C (cid:16) L X l =0 k Ψ l k L (Ω) (cid:17) ≤ C L X l =0 h l · L X l =0 h − l k Ψ l k L (Ω) (5.47) ≤ C | Ψ h | H (Ω) . Because of curl Π Ψ = curl Ψ , we infer from (5.47) k Π Ψ k H ( curl , Ω) + L X l =1 h − l k Π l Ψ l k L (Ω) ≤ C | Ψ h | H (Ω) . (5.49)Above and throughout the remainder of the proof, constants are independent of L .By the projector property Π h ◦ Π l = Π l , l = 0 , . . . , L , and the commuting diagramproperty (2.11), we arrive at v h = e v h + L X l =0 Π l Ψ l + grad ( I h s h + p h ) , (5.50)where I h is the nodal linear interpolation operator onto V ( M h ). Recall (2.22) to seethat |I h s h + p h | H (Ω) ≤ C | s h | H (Ω) + | p h | H (Ω) ≤ C k v h k H ( curl , Ω) . ocal Multigrid in H ( curl ) I h s h + p h according to Thm. 5.7 gives I h s h + p h = L X l =0 p l , p l ∈ Span (cid:8) B lV (cid:9) , (5.51) | p | H (Ω) + L X l =1 h − l k p l k L (Ω) ≤ C |I h s h + p h | H (Ω) ≤ C k v h k H ( curl , Ω) . (5.52)Still, the contribution e v h does not yet match (4.11). The idea is to distribute e v h to the terms Π l Ψ l by scale separation. To that end, we assign a level to each activeedge of M h ℓ ( E ) := min { ℓ ( K ) : K ∈ M h , E is edge of K } , E ∈ E ( M h ) . (5.53)Thus, we distinguish parts of e v h on different levels: given the basis representation e v h = X E ∈E ( M h ) α E b E , { b E } E ∈E ( M h ) = B U ( M h ) , (5.54)we split e v h = L X l =0 e v l , e v l := X E ∈E ( M h ) ℓ ( E )= l α E b E , supp( e v l ) ⊂ ω l . (5.55)The estimate (cid:13)(cid:13) h − e v h (cid:13)(cid:13) L (Ω) ≤ C k curl v h k L (Ω) from Lemma 5.10 means that e v h is“small on fine scales”. Thanks to the L -stability (2.6) of the edge bases, this carriesover to e v l : L X l =0 h − l k e v l k L (Ω) ≤ C L X l =0 h − l X E ∈E ( M h ) ,ℓ ( E )= l α E k b E k L (Ω) ≤ C L X l =0 h − l X E ∈E ( M h ) ,ℓ ( E )= l α E k b E k L ( T E ) (5.56) ≤ C L X l =0 h − l k e v h k L (Σ l ) ≤ C (cid:13)(cid:13) h − e v h (cid:13)(cid:13) L (Ω) , where T E ∈ M h is coarsest element adjacent to E , cf. (5.53), and refinement stripsare defined by Σ l := ω l \ ω l +1 , ≤ l < L, Σ L := ω L , (5.57)see Figs. 4.1 and 4.2.Yet, in the case of bisection refinement, e v l may not be spanned by basis functionsin B l U , because the basis functions of U ( M h ) attached to each edge on Σ l T ω l +1 ,0 ≤ l < L do not belong to any B l U !Take any E ⊂ Σ l T ω l +1 . Let b E , b lE , and b l +1 E be the basis functions of U ( M h ), U ( M l ), and U ( M l +1 ) associated with E , see Fig. 5.3 for a 2D illustration. Denoteby K , . . . , K n all elements in ω l +1 and M l which contain E , and by E , . . . , E m their2 R. Hiptmair and W.-Y. Zheng
PSfrag replacements E Edge E , support of basisfunction b E PSfrag replacements E Support of b lE PSfrag replacements E E Edges supporting b l +1 E , b l +1 E Fig. 5.3 . Basis function with which b E can be represented new edges connecting E but not contained in the refinement edges of K , . . . , K n (see Fig. 5.4). Supposing the orientations of each E i and E point to their commonendpoint, we have b E = b lE + 12 m X i =1 b l +1 E i . (5.58)This decomposition is L -stable with constants merely depending on shape regularity. E K i E j Fig. 5.4 . Situation at an edge E lying on the interface between Σ l and ω l +1 . Since P mi =1 b l +1 E i ∈ B l +1 U , we may move the component of e v l associated with thisterm to e v l +1 for any E . Then the decomposition (5.55) and the stability estimate(5.56) remain valid.Summing up, the stability estimate (5.49) is preserved after replacing Π l Ψ h with Π l Ψ h + e v l ∈ U ( M l ).Eventually, the proof of Thm. 4.2 is readily accomplished. With Lemma 5.11at our disposal, we merely appeal to the L -stabilities expressed in (2.6) and (2.9) ocal Multigrid in H ( curl ) B l U and B lV , respectively.
6. Quasi-orthogonality.
The strengthened Cauchy-Schwartz inequality (4.7)has been established in [61, 65] for H -conforming linear Lagrangian finite elementspaces, in [28, Sect. 6] for H (div)-elliptic variational problems and so-called face ele-ments. It is discussed in [31, Sect. 4] for (1.1), edge elements, and geometric multigridwith global refinement. The considerations for locally refined meshes are fairly similar,but will be elaborated for the sake of completeness.The trick is, not to consider the one-dimensional spaces spanned by individualbasis functions as building blocks of the splitting (4.3), but larger aggregates. Thus,we put the nodal basis functions in B lV and B l U into a small number of classes, suchthat the supports of any two basis functions in the same class do not overlap. Sincethese basis functions are attached to vertices and edges respectively, the definition ofthose classes can be based on a partitioning the vertices/edges of M l | ω l into disjointsets such that any two vertices/edges of the same set do not belong to the sametetrahedron. The is formally stated in the following “colouring lemma”: Lemma 6.1.
There exist P N , P E ∈ N depending only on shape regularity such thatthe sets N ( M l ) | ω l and E ( M l ) ω l of vertices and edges of M l inside the refinementzone ω l can be partitioned into subsets N ( M l ) | ω l := { p ∈ V ( M l ) , supp b p ⊂ ω l } = N l ∪ · · · ∪ N P N l , E ( M l ) | ω l := { E ∈ E ( M l ) , supp b E ⊂ ω l } = E l ∪ · · · ∪ E P E l , and for any K ∈ M l , K ⊂ ω l , two of its vertices/edges will belong to different subsets. Here, b p is the nodal basis function of V ( c M l ) attached to the vertex p , and b E is thenodal basis function of U ( c M l ) associated with the edge E , see (2.4). Proof . A crude argument cites the fact that each vertex and each edge belongs toonly a finite number of elements. A bound for this number can be deduced from theshape regularity measure. The rest is elementary combinatorial arguments.Next, define subspaces of V ( c M l ) and U ( c M l ) by V il := Span (cid:8) b p , p ∈ N il (cid:9) ⊂ Span (cid:8) B lV (cid:9) , i = 1 , . . . , P N , U il := Span (cid:8) b E , E ∈ E il (cid:9) ⊂ Span (cid:8) B l U (cid:9) , i = 1 , . . . , P E . Note that the basis functions spanning both V il and U il are mutually orthogonal(w.r.t. a and the H (Ω)-inner product). Thus, it suffices to establish the strength-ened Cauchy-Schwarz inequality (4.7) for the family of subspaces { H j } j = { U il } l,i ∪{ grad V il } l,i of U ( M L ). This will yield the relevant constants in (4.8). In other words,we analyze the quasi-orthogonality property of the multilevel decomposition U ( M h ) = U ( M ) + L X l =1 P N X i =1 grad V il + L X l =1 P E X i =1 U il . (6.1)Note that (6.1) gives rise to a multigrid algorithm, for which Thm. 4.1 gives exactlythe same convergence estimate as for the method induced by (4.11)! Lemma 6.2.
For all v m ∈ U ( M m ) and u il ∈ U il , ≤ m ≤ l ≤ L , i = 1 , . . . , P E , itholds that, with C > depending only on the bound for the shape regularity measures R. Hiptmair and W.-Y. Zheng of the meshes M l , (cid:0) curl v m , curl u il (cid:1) L (Ω) ≤ Ch l h − m k curl v m k L (Ω) (cid:13)(cid:13) curl u il (cid:13)(cid:13) L (Ω) , (6.2) (cid:0) v m , u il (cid:1) L (Ω) ≤ Ch l k v m k L (Ω) (cid:13)(cid:13) curl u il (cid:13)(cid:13) L (Ω) . (6.3) Proof . Pick any (open) tetrahedron K ∈ c M m . Use the basis representation of u il to isolate “interior” and “boundary” parts u il | K = X E ⊂ ¯ K,E ∈E il (cid:16)Z E u il · d ~s (cid:17) · b E = u il,bd + u il,int , where u il,bd := X E ⊂ ∂K,E ∈E il (cid:16)Z E u il · d ~s (cid:17) · b E and u il,int = X E ⊂ K,E ∈E il (cid:16)Z E u il · d ~s (cid:17) · b E . Since curl v m is a constant vector in K and u il,int × n = on ∂K , by Green’s formula,it is easy to see Z K curl u il · curl v m d x = Z K curl u il,bd · curl v m d x = Z Σ curl u il,bd · curl v m d x , where Σ := [ { supp b E : E ⊂ ∂K, E ∈ E il } is contained in a narrow strip along the boundary of K of width ≈ h l . Hence, wearrive at the area ratio | Σ | ≤ Ch l h − m | K | . Here and throughout the remainder of the proof,
C > U il are mutually orthogonal, we have Z K curl u il · curl v m d x ≤ (cid:13)(cid:13) curl u il,bd (cid:13)(cid:13) L (Σ) | Σ | / | curl v m | (6.4) ≤ C r h l h m (cid:13)(cid:13) curl u il (cid:13)(cid:13) L ( K ) | K | / | curl v m | = C r h l h m (cid:13)(cid:13) curl u il (cid:13)(cid:13) L ( K ) k curl v m k L ( K ) . To estimate the L -inner product, we recall the following simple fact about the normsof edge basis functions on level l : k b k L ( K ) ≤ Ch l k curl b k L ( K ) ∀ b ∈ B l U . Since the basis functions of U il do not interact, we have Z K u il · v m d x ≤ (cid:13)(cid:13) u il (cid:13)(cid:13) L ( K ) k v m k L ( K ) ≤ Ch l (cid:13)(cid:13) curl u il (cid:13)(cid:13) L ( K ) k v m k L ( K ) . (6.5) ocal Multigrid in H ( curl ) M m and anotherCauchy-Schwarz inequality.After replacing U il with V il in the proof of Lemma 6.2, similar arguments establishthe following estimate: Lemma 6.3.
For all v m ∈ U ( M m ) and u il ∈ V il , ≤ m ≤ l ≤ L , i = 1 , . . . , P N , itholds that, with C > depending only on the bound for the shape regularity measuresof the meshes M l , (cid:0) v m , grad u il (cid:1) L (Ω) ≤ Ch l h − m k v m k L (Ω) (cid:13)(cid:13) grad u il (cid:13)(cid:13) L (Ω) . (6.6) Proof . Again, pick K ∈ M m . By separating interior and boundary parts of u il asabove and noting div v m | K = 0 on K , we find by Green’s formula Z K v m · grad u il d x = Z Σ v m · grad u il d x . As above, we infer Z K v m · grad u il d x ≤ C r h l h m k v m k L ( K ) (cid:13)(cid:13) grad u il (cid:13)(cid:13) L ( K ) . Summation over all K and a Cauchy-Schwarz inequality finish the proof.Because of the geometric decay of the meshwidths h l of the (uniformly refined)meshes c M l , these estimates clearly imply the desired quasi-orthogonality for (6.1). Theorem 6.4. (Strengthened Cauchy-Schwartz inequality) For any u il ∈ U il or u il ∈ grad V il and any v jl ∈ U jl or v jl ∈ grad V jl , ≤ i, j ≤ P E or ≤ i, j ≤ P N ,resp., the estimate (6.7) a ( u il , v jm ) ≤ Cθ | l − m | / (cid:13)(cid:13) u il (cid:13)(cid:13) A (cid:13)(cid:13) v jm (cid:13)(cid:13) A ≤ i, j ≤ N l , ≤ l, m ≤ L holds, where C = C (max l ρ l ) > and < θ < is the decrease rate of the meshwidthsdefined in (3.3) .
7. Numerical experiments.
In the reported numerical experiments the im-plementation of adaptive mesh refinement was based on the adaptive finite elementpackage ALBERTA [48], which uses the bisection strategy of [36], see Sect. 3.Let M be an initial mesh satisfying the two assumptions (A1) and (A2) in [36,P. 282], the adaptive mesh refinements are governed by a residual based a posteriorierror estimator. In the experiments we assume the current density f ∈ H (div , Ω) anduse the estimator given by [17, § u h ∈ U ( M h ),for any T ∈ M h η T := h T k f − u h k H (div ,T ) + h T X F ⊂ ∂T n k [ u h ] F k ,F + k [ curl u h × ν ] F k ,F o , where F is a face of T , ν is the unit normal of F , and [ u h ] F is the jump of u h across F . The global a posteriori error estimate and the maximal estimated element erroron M h are defined by(7.1) η h := X T ∈M h η T ! / , η max = max T ∈M h η T . R. Hiptmair and W.-Y. Zheng
Using η h and η max , we use [17, Algorithm 5.1] to mark and refine M h adaptively.In the following, we report two numerical experiments to demonstrate the compet-itive behavior of the local multigrid method and to validate our convergence theory. example 7.1. We consider the Maxwell equation on the three-dimensional “L-shaped” domain
Ω = ( − , \ { (0 , × ( − , × ( − , } . The Dirichlet boundarycondition and the righthand side f are chosen so that the exact solution is u := ∇ n r / sin( φ/ o in cylindrical coordinates ( r, φ, z ) . Table 7.1 shows the numbers of multigrid iterations required to reduce the initialresidual by a factor 10 − on different levels. We observe that the multigrid algorithmconverges in almost the same small number of steps, though the number of elementsvaries from 156 to 100,420. Table 7.1
The number of adaptive iterations N it , the number of elements N el , the number of multigriditerations N itrs required to reduce the initial residual by a factor − , the relative error between thetrue solution u and the discrete solution u h : E rel = k u − u h k H ( curl , Ω) / k u k H ( curl , Ω) (Example7.1). N it N el
156 388 1,900 4,356 9608 19,424 48,088 100,420 E rel N itrs
11 21 19 19 19 19 19 19Fig. 7.1 (left) plots the CPU time versus the number of degrees of freedom ondifferent adaptive meshes. It shows that the CPU time of solving the algebraic systemincreases roughly linearly with respect to the number of elements. Fig. 7.1 (rught) de-picts a locally refined mesh of 100,420 elements created by the adaptive finite elementalgorithm. S e c ond Fig. 7.1 . Example 7.1, left: execution time for local multigrid method, right: instance of a locallyrefined mesh (100,420 elements) example 7.2.
This example uses the same solution as Example 7.1 u := ∇ n r / sin( φ/ o ocal Multigrid in H ( curl ) in cylindrical coordinates ( r, φ, z ) . But the computational domain is changed to a three-dimensional non-Lipschitz domain with an inner crack-type boundary, which is definedby Ω = ( − , \ { ( x, , z ) : 0 ≤ x < , − < z < } . The Dirichlet boundary condition and the source function f are the same as above. Table 7.2 records the numbers of multigrid iterations required to reduce the initialresidual by a factor 10 − on different levels. We observe that the multigrid algorithmconverges in less than 30 steps, with the number of elements soaring from 128 to135,876. Table 7.2
The number of adaptive iterations N it , the number of elements N el , the number of multigriditerations N itrs required to reduce the initial residual by a factor − , the relative error between thetrue solution u and the discrete solution u h : E rel = k u − u h k H ( curl , Ω) / k u k H ( curl , Ω) (Example7.2). N it N el
128 404 1,236 3,416 12,420 29,428 81,508 135,876 E rel N itrs
14 30 25 26 26 27 27 27Fig. 7.2 (left) shows the CPU time versus the number of degrees of freedom ondifferent adaptive meshes. Obviously, the CPU time for solving the algebraic systemincreases nearly linearly with respect to the number of elements.Fig. 7.2 (right) displays a locally refined mesh of 135,876 elements using adaptivefinite element algorithm. In addition, the restriction of the mesh to the cross-section { y = 0 } , which contains the inner boundary, is drawn. This reveals strong localrefinement.This experiment bears out that the local multigrid is also efficient for the problemsin non-Lipschitz doamins, which are outside the scope of our theory. S e c ond Fig. 7.2 . Example 7.2, left: CPU time for solving the algebraic system by multigrid method,right: a locally refined mesh (135,876 elements) R. Hiptmair and W.-Y. Zheng
Acknowledgement.
The authors would like to thank Dr. L. Wang of ComputerNetwork Information Center, Prof. Z. Chen and Prof. L. Zhang of the Institute ofComputational Mathematics, Chinese Academy of Sciences, for their support in theimplementation of the local multigrid method. They are grateful to one referees whodetected an error in an earlier version of the manuscript.
REFERENCES[1]
M. Ainsworth and W. McLean , Multilevel diagonal scaling preconditioners for boundaryelement equations on locally refined meshes , Numer. Math., 93 (2003), pp. 387–413.[2]
D. Arnold, R. Falk, and R. Winther , Multigrid in H (div) and H ( curl ), Numer. Math., 85(2000), pp. 175–195.[3] , Finite element exterior calculus, homological techniques, and applications , Acta Nu-merica, 15 (2006), pp. 1–155.[4]
D. Arnold, A. Mukherjee, and L. Pouly , Locally adapted tetrahedral meshes using bisection ,SIAM Journal on Scientific Computing, 22 (2000), pp. 431–448.[5]
D. Bai and A. Brandt , Local mesh refinement multilevel techniques , SIAM J. Sci. Stat. Com-put., 8 (1987), pp. 109–134.[6]
E. B¨ansch , Local mesh refinement in 2 and 3 dimensions , IMPACT Comput. Sci. Engrg., 3(1991), pp. 181–191.[7]
R. Beck, P. Deuflhard, R. Hiptmair, R. Hoppe, and B. Wohlmuth , Adaptive multilevelmethods for edge element discretizations of Maxwell’s equations , Surveys on Mathematicsfor Industry, 8 (1998), pp. 271–312.[8]
R. Beck, R. Hiptmair, R. Hoppe, and B. Wohlmuth , Residual based a-posteriori errorestimators for eddy current computation , M AN, 34 (2000), pp. 159–182.[9]
J. Bey , Tetrahedral grid refinement , Computing, 55 (1995), pp. 355–378.[10]
P. Binev, W. Dahmen, and R. DeVore , Adaptive finite element methods with convergencerates , Numerische Mathematik, 97 (2004), pp. 219–268.[11]
F. Bornemann and H. Yserentant , A basic norm equivalence for the theory of multilevelmethods , Numer. Math., 64 (1993), pp. 455–476.[12]
A. Bossavit , Whitney forms: A class of finite elements for three-dimensional computations inelectromagnetism , IEE Proc. A, 135 (1988), pp. 493–500.[13] ,
Computational Electromagnetism. Variational Formulation, Complementarity, EdgeElements , vol. 2 of Electromagnetism Series, Academic Press, San Diego, CA, 1998.[14]
J. Bramble and J. Pasciak , New estimates for multilevel methods including the V–cycle ,Math. Comp., 60 (1993), pp. 447–471.[15]
J. Bramble, J. Pasciak, J. Wang, and J. Xu , Convergence estimates for product iterativemethods with applications to domain decomposition , Math. Comp., 57 (1991), pp. 1–21.[16]
A. Buffa, M. Costabel, and D. Sheen , On traces for H ( curl , Ω) in Lipschitz domains , J.Math. Anal. Appl., 276 (2002), pp. 845–867.[17] Z.-M. Chen, L. Wang, and W.-Y. Zheng , An adaptive multilevel method for time-harmonicMaxwell equations with singularities , SIAM J. Sci. Comp., (2006). To appear.[18]
S. Christiansen and R. Winther , Smoothed projections in finite element exterior calcu-lus
P. Ciarlet , The Finite Element Method for Elliptic Problems , vol. 4 of Studies in Mathematicsand its Applications, North-Holland, Amsterdam, 1978.[20]
M. Clemens, S. Feigh, and T. Weiland , Geometric multigrid algorithms using the conformalfinite integration technique , IEEE Trans. Magnetics, 40 (2004), pp. 1065–1068.[21]
P. Cl´ement , Approximation by finite element functions using local regularization , RAIROAnal. Num´er., 2 (1975), pp. 77–84.[22]
M. Costabel and M. Dauge , Singularities of electromagnetic fields in polyhedral domains ,Arch. Rational Mech. Anal., 151 (2000), pp. 221–276.[23]
M. Costabel, M. Dauge, and S. Nicaise , Singularities of eddy current problems , ESAIM:Mathematical Modelling and Numerical Analysis, 37 (2003), pp. 807–831.[24]
V. Girault and P. Raviart , Finite element methods for Navier–Stokes equations , Springer,Berlin, 1986.[25]
J. Gopalakrishnan, J. Pasciak, and L. Demkowicz , Analysis of a multigrid algorithm fortime harmonic Maxwell equations , SIAM J. Numer. Anal., 42 (2003), pp. 90–108.[26]
V. Gradinaru and R. Hiptmair , Whitney elements on pyramids , Electron. Trans. Numer.ocal Multigrid in H ( curl ) Anal., 8 (1999), pp. 154–168.[27]
P. Grisvard , Elliptic Problems in Nonsmooth Domains , Pitman, Boston, 1985.[28]
R. Hiptmair , Multigrid method for Maxwell’s equations , Tech. Rep. 374, Institut f¨ur Mathe-matik, Universit¨at Augsburg, 1997. USE HIP 99.[29] ,
Multigrid method for Maxwell’s equations , SIAM J. Numer. Anal., 36 (1999), pp. 204–225.[30] ,
Finite elements in computational electromagnetism , Acta Numerica, 11 (2002), pp. 237–339.[31] ,
Analysis of multilevel methods for eddy current problems , Math. Comp., 72 (2003),pp. 1281–1303.[32]
R. Hiptmair, G. Widmer, and J. Zou , Auxiliary space preconditioning in H ( curl , Ω), Numer.Math., 103 (2006), pp. 435–459.[33]
R. Hiptmair and J. Xu , Nodal auxiliary space preconditioning in H(curl) and H(div) spaces ,SIAM J. Numer. Anal., 45 (2007), pp. 2483–2509.[34]
R. Hoppe and J. Sch¨oberl , Convergence of adaptive edge element methods for the 3D eddycurrents equation , J. Comp. Math., (2008).[35]
T. Kolev and P. Vassilevski , Parallel auxiliary space AMG for H ( curl ) problems , J. Comp.Math., (2008).[36] I. Kossaczk´y , A recursive approach to local mesh refinement in two and three dimensions , J.Comput. Appl. Math., 55 (1994), pp. 275–288.[37]
J. Lions and F. Magenes , Nonhomogeneous boundary value problems and applications ,Springer–Verlag, Berlin, 1972.[38]
J. Maubach , Local bisection refinement for n –simplicial grids generated by reflection , SIAMJ. Sci. Stat. Comp., 16 (1995), pp. 210–227.[39] W. McLean , Strongly Elliptic Systems and Boundary Integral Equations , Cambridge Univer-sity Press, Cambridge, UK, 2000.[40]
W. Mitchell , A comparison of adaptive refinement techniques for elliptic problems , ACMTrans. Mathematical Software, 15 (1989), pp. 326–347.[41] ,
Optimal multilevel iterative methods for adaptive grids , SIAM J. Sci. Stat. Comput, 13(1992), pp. 146–167.[42]
P. Monk , Finite Element Methods for Maxwell’s Equations , Clarendon Press, Oxford, UK,2003.[43]
J. N´ed´elec , Mixed finite elements in R , Numer. Math., 35 (1980), pp. 315–341.[44] P. Oswald , On function spaces related to the finite element approximation theory , Z. Anal.Anwendungen, 9 (1990), pp. 43–64.[45] ,
On discrete norm estimates related to multilevel preconditioners in the finite elementmethod , in Constructive Theory of Functions, Proc. Int. Conf. Varna 1991, K. Ivanov,P. Petrushev, and B. Sendov, eds., Bulg. Acad. Sci., 1992, pp. 203–214.[46] ,
Multilevel finite element approximation , Teubner Skripten zur Numerik, B.G. Teubner,Stuttgart, 1994.[47]
J. Pasciak and J. Zhao , Overlapping Schwarz methods in H(curl) on polyhedral domains , J.Numer. Math., 10 (2002), pp. 221–234.[48]
A. Schmidt and K. Siebert , ALBERTA – An adaptive hierarchical finite element toolbox
A. Schmidt and K. Siebert , Design of Adaptive Finite Element Software: The Finite ElementToolbox ALBERTA , Lecture Notes in Computational Science and Engineering, Springer,Heidelberg, 2005.[50]
J. Sch¨oberl , Commuting quasi-interpolation operators for mixed finite elements , PreprintISC-01-10-MATH, Texas A&M University, College Station, TX, 2001.[51] ,
A multilevel decomposition result in H ( curl ), in Proceedings of the 8th European Multi-grid Conference 2005, Scheveningen, P. H. P. Wesseling, C.W. Oosterlee, ed., 2006.[52] , A posteriori error estimates for Maxwell equations , Math. Comp., 77 (2008), pp. 633–649.[53]
L. R. Scott and Z. Zhang , Finite element interpolation of nonsmooth functions satisfyingboundary conditions , Math. Comp., 54 (1990), pp. 483–493.[54]
O. Sterz, A. Hauser, and G. Wittum , Adaptive local multigrid methods for solving time-harmonic eddy current problems , IEEE Trans. Magnetics, 42 (2006), pp. 309–318.[55]
R. Stevenson , Optimality of a standard adaptive finite element method , Foundations of Com-putational Mathematics, 7 (2007), pp. 245–269.[56] ,
The completion of locally refined simplicial partitions created by bisection , Math. Comp.,77 (2008), pp. 227–241.[57]
C. Traxler , An algorithm for adaptive mesh refinement in n dimensions , Computing, 59 R. Hiptmair and W.-Y. Zheng(1997), pp. 115–137.[58]
B. Weiss and O. Biro , Multigrid for time-harmonic 3-d eddy-current analysis with edge ele-ments , IEEE Trans. Magnetics, 41 (2005), pp. 1712–1715.[59]
H. Whitney , Geometric Integration Theory , Princeton University Press, Princeton, 1957.[60]
H.-J. Wu and Z.-M. Chen , Uniform convergence of multigrid V -cycle on adaptively refinedfinite element meshes for second order elliptic problems , Science in China: Series A Math-ematics, 49 (2006), p. 1C28.[61] J. Xu , Iterative methods by space decomposition and subspace correction , SIAM Review, 34(1992), pp. 581–613.[62] ,
An introduction to multilevel methods , in Wavelets, Multilevel Methods and EllipticPDEs, M. Ainsworth, K. Levesley, M. Marletta, and W. Light, eds., Numerical Mathematicsand Scientific Computation, Clarendon Press, Oxford, 1997, pp. 213–301.[63]
J. Xu and Y.-R. Zhu , Uniformly convergent multigrid methods for elliptic problems withstrongly discontinuous coefficients , Math. Models Methods Appl. Sci., 18 (2008), pp. 77–105.[64]
J. Xu and L. Zikatanov , The method of alternating projections and the method of subspacecorrections in Hilbert space , J. Am. Math. Soc., 15 (2002), pp. 573–597.[65]
H. Yserentant , On the multi–level splitting of finite element spaces , Numer. Math., 58 (1986),pp. 379–412.[66] ,
Old and new convergence proofs for multigrid methods , Acta Numerica, (1993), pp. 285–326.[67]