# A unified error analysis of HDG methods for the static Maxwell equations

We propose a framework that allows us to analyze different variants of HDG methods for the static Maxwell equations using one simple analysis. It reduces all the work to the construction of projections that best fit the structures of the approximation spaces. As applications, we analyze four variants of HDG methods (denoted by B, H, B+, H+), where two of them are known (variants H, B+) and the other two are new (variants H+, B). We show that all the four variants are optimally convergent and that variants B+ and H+ achieve superconvergence without post-processing. For the two known variants, we prove their optimal convergence under weaker requirements of the meshes and the stabilization functions thanks to the new analysis techniques being introduced. At the end, we provide numerical experiments to support the analysis.

## Authors

• 2 publications
• 6 publications
04/08/2021

### Analysis of Normal-Form Algorithms for Solving Systems of Polynomial Equations

We examine several of the normal-form multivariate polynomial rootfindin...
12/31/2021

### C^1-VEM for some variants of the Cahn-Hilliard equation: a numerical exploration

We consider the C^1-Virtual Element Method (VEM) for the conforming nume...
02/15/2022

### Stochastic Gradient Descent-Ascent: Unified Theory and New Efficient Methods

Stochastic Gradient Descent-Ascent (SGDA) is one of the most prominent a...
12/30/2021

### DPG methods for a fourth-order div problem

We study a fourth-order div problem and its approximation by the discont...
09/09/2021

### Versions of Gradient Temporal Difference Learning

Sutton, Szepesvári and Maei introduced the first gradient temporal-diffe...
02/18/2020

### A note on devising HDG+ projections on polyhedral elements

In this note, we propose a simple way of constructing HDG+ projections o...
07/02/2014

### How Many Dissimilarity/Kernel Self Organizing Map Variants Do We Need?

In numerous applicative contexts, data are too rich and too complex to b...
##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

Maxwell equations describe the interaction between electric and magnetic fields and play a central role in modern sciences and engineering. To understand the solution of Maxwell equations in various application scenarios, numerical treatments are necessary. The finite element method (FEM) is one of these numerical tools and it has some nice features such as easy handling of complicate geometry, exponential rate of convergence by hp-refinements, etc.

Finite element methods can be divided into two categories – conforming and non-conforming. For Maxwell equations, conforming elements usually refer to -conforming elements since is used as the energy space for the solution of Maxwell equations. -conforming elements (also called edge elements) have been widely studied since they were first proposed by Nédélec in [29, 30]; see, for instance, [19, 21, 26, 27, 28, 35].

For non-conforming elements, one popular choice is the discontinuous Galerkin (DG) finite element method (see [1] for a general introduction and see, for instance, [4, 15, 18, 22, 32, 33]

for DG methods for Maxwell equations). Since DG methods allow the use of independent approximation spaces on each element, they possess certain nice properties such as the flexibility of choosing local spaces, allowance of triangulation with hanging nodes, high parallel efficiency, easiness of implementation, simple treatment of boundary conditions, etc. Despite their advantages, DG methods in general use more degrees of freedom compared to the corresponding conforming methods. To overcome this difficulty, the hybridizable discontinuous Galerkin (HDG) method was proposed

[13]. By introducing a Lagrange multiplier on the skeleton of the mesh and using the hybridization techniques, HDG method allows the solution of a much smaller system only involving the Lagrange multiplier and then to recover locally the rest of the degrees of freedom on each element.

Recently, there has been considerable interest in developing HDG methods for Maxwell equations and many variants [7, 8, 9, 24, 25, 31] of HDG methods have been proposed and analyzed. However, to the best of our knowledge, there is no work that provides a unified point of view of understanding these variants. This leads to a possibility of repeated or unnecessary arguments being generated and a lack of recognition of the connections among these variants. This motivates us to consider a unified analysis. In this paper, we propose a framework that enables us to clearly decouple the error analysis techniques into two groups – those related to the PDE and those related to the HDG variants (namely, the choices of the approximation spaces and stabilization functions). The benefits of doing so include the following:

• Recycling existing error analysis techniques. We demonstrate this by using only one analysis to obtain the error estimates for four variants of HDG methods. In this way, we can avoid introducing repeated arguments for each variants.

• Providing guidelines for systematically discovering new optimal convergent and super convergent HDG methods. We discover two new HDG variants and by using this framework, where variant achieves superconvergence in the sense of the degrees of freedom of the numerical trace (the discrete electric field achieves convergence while its numerical trace only lives in a proper subspace of on each face ; see the end of Section 4 for a detailed discussion about this).

• Simple analysis of mixed type HDG methods where the local spaces and stabilization functions vary from element to element. This is doable since we use local projections to capture the features of the HDG variants (the main part of which is how to choose local approximation spaces and stabilization functions).

Let us mention two inspirations of this work. The first one is [14], where a tailored projection is proposed to analyze a class of HDG methods in a unified way under the setting of elliptic problem (this is inspired by the celebrated Raviart-Thomas (RT) and Brezzi-Douglas-Marini (BDM) projections). This approach to the analysis is often referred to as “projection-based error analysis”. The second one is our previous work about HDG methods for elastic waves [17], in which we show that we can use projection-based error analysis for those HDG methods whose approximation spaces do not admit -decomposition [10]. The work of this paper can be regarded as a generalization of the work in [17] to the setting of Maxwell equations.

To proceed with the discussion, we shall now introduce the model problem. Let be a bounded simply connected polyhedral domain with connected Lipschitz boundary . We consider the following static Maxwell equations in a mixed form:

 w−∇×u =0 in Ω, (1a) ∇×w+∇p =f in Ω, (1b) ∇⋅u =0 in Ω, (1c) n×u =g on Γ, (1d) p =0 on Γ. (1e)

In the above, variables and are the electric and the magnetic fields respectively, and is a Lagrange multiplier introduced to have a better control of (see [2, 3]). Note that when is divergence free, admits trivial solution. We remark that (1) with a different boundary condition can be also derived from the Stokes equations by using vorticity formulations; see, for instance, [11, 12].

The rest of the paper is organized as follows. In Section 2, we propose an HDG framework with unspecified approximation spaces and stabilization functions; we then give an analysis by using a projection satisfying certain criteria. In Section 3, we review some well known projections and construct some new projections that we shall use later. In Section 4, we consider four variants of HDG methods for Maxwell equations (denoted by , , , ). We give a unified analysis to the four variants by using the abstract analysis setting established in Section 2 combined with suitable projections discussed in Section 3. We show that all the variants are optimal and variants , achieve superconvergence. After that, we compare these four variants and discuss their connections. Finally in Section 5, we present some numerical tests to support the analysis.

## 2 The framework

### 2.1 Notation

We begin by introducing some notation that will be used extensively in the paper. Let be a conforming triangulation of , where each element is a star-shaped polyhedron. Let and be the collections of all faces of and , respectively. We use the standard notation as the diameter of and denote by as the mesh size of . For , we denote by the polynomial space of degree supported on , where can be an element in or a face in . Let be a large integer. For any , let and be two subspaces of , and be a subspace of . For any , let be a subspace of

(for a vector field

supported on certain surface , we denote by the tangential component of the vector field), and be a subspace of . Denote by and . Let and be the projections to their range spaces respectively. We assume all the spaces introduced above are non-empty. Define

 Wh:=∏K∈ThW(K),Vh:=∏K∈ThV(K),Qh:=∏K∈ThQ(K), Nh:=∏F∈EhN(F),Mh:=∏F∈EhM(F).

We use the following notation for the discrete inner products on and :

 (∗1,∗2)Th=∑K∈Th(∗1,∗2)K,⟨∗1,∗2⟩∂Th=∑K∈Th⟨∗1,∗2⟩∂K,

where and denote the inner products on and respectively.

### 2.2 HDG methods

Depending on the choices of the approximation spaces and , we obtain different variants of HDG methods. We assume these spaces satisfy the following conditions:

 ∇×V(K)⊂W(K) (3a) ∇⋅V(K)⊂Q(K) (3b) ∇×W(K)+∇Q(K)⊂V(K) (3c) n∂K×W(K)⊂N(∂K) (3d) γ∂KQ(K)+V(K)⋅n∂K⊂M(∂K). (3e)

All the HDG variants we will study in this paper satisfy (3) and we assume these conditions hold throughout the paper. We now give the HDG scheme under this general setting:
Find such that

 (wh,r)Th−(uh,∇×r)Th−⟨ˆuh,r×n⟩∂Th =0, (4a) (∇×wh,v)Th+⟨τtPN(uh−ˆuh),v⟩∂Th (4b) −(ph,∇⋅v)Th+⟨ˆph,v⋅n⟩∂Th =(f,v)Th, (4c) (∇⋅uh,q)Th+⟨τn(ph−ˆph),q⟩∂Th =0, (4d) −⟨n×wh+τt(uh−ˆuh),η⟩∂Th∖Γ =0, (4e) −⟨ˆuh,η⟩Γ =−⟨g×n,η⟩Γ, (4f) −⟨uh⋅n+τn(ph−ˆph),μ⟩∂Th∖Γ =0, (4g) −⟨ˆph,μ⟩Γ =0, (4h)

for all . In the above equations (4), the two stabilization functions and we assume for all . It is obvious that (4) is a square system. We remark that the unique solvability can be deduced as a consequence of the convergence of the numerical scheme, which we will study in Section 4.

### 2.3 Projections and remainders

The key in our analysis is finding projections satisfying the following Assumption 2.1. Here, under this general setting, we shall just assume the projection exists and proceed the analysis. We remark that these projections are not unique in most cases and our target is to find the projections that can well fit the structures of the approximation spaces and therefore give sharp estimates.

###### Assumption 2.1 (Projection assumption).

For all , there exists a projection

 ΠK:H1(K)3×H1(K)3×H1(K) →W(K)×V(K)×Q(K) (w,u,p) ↦(ΠKw,ΠKu,ΠKp),

such that

 (ΠKw−w,∇×v)K =⟨n×w−PN(n×w),v⟩∂K ∀v∈V(K), (5a) (ΠKu−u,v)K =0 ∀v∈∇×W(K)+∇Q(K), (5b) (ΠKp−p,∇⋅v)K =0 ∀v∈V(K). (5c)

Note that if we have , then (5a) becomes

 (ΠKw−w,∇×v)K=0∀v∈V(K).

In this case, Assumption 2.1 holds obviously as a result of (3a)-(3c), since the projection to satisfy (5). In addition, we have used , , and to represent the first, second, and third component of the projection , respectively. Hence can depend on and as well, and this clarification works similarly for and .

For all the HDG variants we will study in this paper, Assumption 2.1 is satisfied. Namely, we can explicitly construct projections that satisfy (5a)-(5c). We will do this in Section 4.

We next define two operators associated to the projection .

###### Definition 2.1 (Boundary remainders).

For all , we define two operators as follows:

 δΠK±τt:H1(K)3×H1(K)3 →N(∂K) (w,u) ↦n×ΠKw−PN(n×w)±τt(PNΠKu−PNu), (6a) δΠK±τn:H1(K)3×H1(K) →M(∂K) (u,p) ↦ΠKu⋅n−PM(u⋅n)±τn(ΠKp−PMp). (6b)

We call the curl-curl boundary remainder and the grad-div boundary remainder.

By (3d) and (3e), it is easy to see that the above definition is valid. The boundary remainder operators can be regarded as an indicator for how much the projection resembles an HDG projection or a mixed method projection. Consider the grad-div boundary remainder . If we let the second-third component of , namely , to be replaced by the HDG projection with stabilization function , then (holds by definition; see [14]); if the second-third component is replaced by the Raviart-Thomas (RT) and the Brezzi-Douglas-Marini (BDM) projections [6, 34] ([29, 30] by Nédélec for case), we have . On the other hand, if the first-second component of is replaced by edge element associated projections ( projections [29, 30]), then we have .

The following Lemma gives two identities further relating the projection and its associated two boundary remainders.

###### Lemma 2.1 (Weak-commutativity).

For all , denote by and for simplicity. Then

 (∇×(ΠKw−w),v)K±⟨τt(PNΠKu−PNu),v⟩∂K =⟨δΠK±τt,v⟩∂K∀v∈V(K), (7a) (∇⋅(ΠKu−u),q)K±⟨τn(ΠKp−PMp),q⟩∂K =⟨δΠK±τn,q⟩∂K∀q∈Q(K). (7b)
###### Proof.

First note that

 (∇×(ΠKw−w),v)K±⟨τt(PNΠKu−PNu),v⟩∂K =⟨n×(ΠKw−w)±τt(PNΠKu−PNu),v⟩∂K+(ΠKw−w,∇×v)K,

for all . Equation (7a) now follows by using (5a). Equation (7b) can be similarly obtained by using (5b) and (3e). ∎

### 2.4 Estimates

Energy estimates. To proceed with the analysis, we assume Assumption 2.1 is satisfied so that we have a projection satisfying (5a)-(5c) for each . We next define the elementwise projections and associated boundary remainders:

 (Πw,Πu,Πp)=∏K∈ThΠK(w,u,p),δΠτt=∏K∈ThδΠKτt(w,u),δΠτn=∏K∈ThδΠKτn(u,p).

We also define the error terms to simplify notation:

 εwh=Πw−wh, εuh=Πu−uh, εph=Πp−p, ˆεuh=PNu−ˆuh, ˆεph=PMp−ˆph.

Note that

 (εwh,εuh,εph,ˆεuh,ˆεph) ∈Wh×Vh×Qh×Nh×Mh, (δΠτt,δΠτn) ∈∏K∈ThN(∂K)×∏K∈ThM(∂K).

For the following two Propositions (Props. 2.1 and 2.2), we put their proofs in the appendix. Once the HDG variants are specified, we can immediately obtain the error estimates of and by using these two propositions.

###### Proposition 2.1 (Energy identity).

The following energy identity holds

 (εwh,εwh)Th+⟨τt(PNεuh−ˆεuh),PNεuh−ˆεuh⟩∂Th+⟨τn(εph−ˆεph),εph−ˆεph⟩∂Th =(Πw−w,εwh)Th+⟨δΠτt,PNεuh−ˆεuh⟩∂Th+⟨δΠτn,εph−ˆεph⟩∂Th. (8)

From the above identity, we can obtain an estimate for .

Duality estimates. To estimate , we consider the dual equations

 w∗+∇×u∗ =0 in Ω, (9a) −∇×w∗−∇p∗ =θ in Ω, (9b) −∇⋅u∗ =0 in Ω, (9c) n×u∗ =0 on Γ, (9d) p∗ =0 on Γ. (9e)
###### Assumption 2.2 (Elliptic regularity).

The following inequality holds

 ∥w∗∥1,Ω+∥u∗∥2,Ω+∥p∗∥1,Ω≤Creg∥θ∥Ω, (10)

for any , where is a constant depending only on .

We remark that this regularity assumption becomes true if is assumed to be convex additionally. Its proof can be obtained by using [20, Theorem 3.5] and then the identity to transform the original formulation (9) to a Poisson’s equation.

Let be another projection satisfying Assumption 2.1. Note that it is allowed to choose . Define

 (Π∗w∗,Π∗u∗,Π∗p∗)=∏K∈ThΠ∗K(w∗,u∗,p∗), δΠ∗−τt=∏K∈ThδΠ∗K−τt(w∗,u∗),δΠ∗−τn=∏K∈ThδΠ∗K−τn(u∗,p∗).
###### Proposition 2.2 (Duality identity).

The following identity holds

 (Πw−w,Π∗w∗)Th+⟨δΠτt,Π∗u∗−PNu∗⟩∂Th+⟨δΠτn,Π∗p∗−PMp∗⟩∂Th =(Π∗w∗−w∗,Πw−wh)Th+(θ,εuh)Th−⟨δΠ∗−τt,PNεuh−ˆεuh⟩∂Th−⟨δΠ∗−τn,εph−ˆεph⟩∂Th. (11)

Let and proceed, we can obtain an estimate for . We will do this in Section 4 when the approximation spaces are specified.

## 3 Projections

In this section, we give a collection of projections which will become the building blocks for constructing projections satisfying Assumption 2.1. Some of these projections are well known while some are newly devised. For those known, we review their constructions and convergence properties. For those new, we prove their optimal convergence under certain shape regularity conditions of the element. We categorize the projections into two groups: (1) Projections for polyhedral element; (2) Projections for simplex element.

### 3.1 Projections for polyhedral element

In this subsection, we focus on one element , which we assume to be a star-shaped polyhedron (we remark that is also allowed to be a simplex). We define the shape-regularity constant of as any constant satisfying the following conditions (see [5, 16, 23]):

• Chunkiness condition. is star-shaped with respect to a ball with radius and .

• Simplex condition. admits a simplex decomposition such that for any simplex , if is the diameter of and is the inradius, then .

• Local quasi-uniformity. Let and be the areas of the largest and smallest face of respectively, then .

projection. For , the orthogonal projection (or projection)

 Πk:L2(K)3 →Pk(K)3, u ↦Πku,

is defined by solving

 (Πku−u,v)K=0∀v∈Pk(K)3. (12)

We have (see [16])

 h1/2K∥Πku−u∥∂K+∥Πku−u∥K≤ChmK|u|m,K, (13)

where and depends only on and the shape-regularity of .

Curl+ projection. We denote by the homogeneous polynomial space of degree and denote by the surface gradient on face . Define

 N(∂K)=∏F∈EKPk(F)t⊕∇F˜Pk+2(F),

and let be the projection to . For , the curl+ projection

 Πck:H1(K)3 →Pk(K)3, w ↦Πckw,

is defined by

 (Πckw−w,r)K =0∀r∈∇×Pk(K)3⊕(∇×Pk+1(K)3)⊥k, (14a) (Πckw−w,∇×v)K =⟨(n×w)−PN(n×w),v⟩∂K∀v∈(Pk(K)3⊕∇˜Pk+2(K))⊥k+1, (14b)

where means taking orthogonal complement in .

By (14a) and (14b), we obtain

 (Πckw−w,∇×v)K =⟨(n×w)−PN(n×w),v⟩∂K∀v∈Pk+1(K)3. (15)

This can be easily proved by decompose , where and . In addition, note that if , then

 (Πckw−w,P0(K)3)K=0, (16)

which can be derived easily from (14a).

###### Theorem 3.1.

The projection is well defined and

 h1/2K∥Πckw−w∥∂K+∥Πckw−w∥K≤ChmK|w|m,K, (17)

where and depends only on and the shape-regularity of .

###### Proof.

See appendix. ∎

This projection will be the key in our analysis of the two HDG variants using Lehrenfeld-Schöberl type stabilization function (variants and ).

### 3.2 Projections for simplex element

In this subsection, we focus on one simplex element in .

HDG projection. Let to shorten notation. For , the HDG projection (see [14])

 ΠHk,τK:H1(K)3×H1(K) →Pk(K)3×Pk(K), (u,p) ↦(ΠHk,τKu,ΠHk,τKp),

is defined by solving

 (ΠHk,τKu−u,r)K =0∀r∈Pk−1(K)3, (18a) (ΠHk,τKp−p,v)K =0∀v∈Pk−1(K), (18b) ⟨(ΠHk,τKu−u)⋅n+τK(ΠHk,τKp−p),μ⟩∂K =0∀μ∈Rk(∂K), (18c)

where and it satisfies either or .

###### Theorem 3.2 ([14]).

For the projection , we have

 h1/2K∥u−ΠHk,τKu∥∂K+∥u−ΠHk,τKu∥K ≤C(hsK|u|s,K+τsecKhtK|p|t,K), (19a) h1/2K∥p−ΠHk,τKp∥∂K+∥p−ΠHk,τKp∥K ≤C(htK|p|t,K+hsKτmaxK|∇⋅u|s−1,K), (19b)

with , where and are the largest and the second largest values of on the faces of respectively.

BDM-H projection. For , define

 ΠBk,τK:H1(K)3 →Pk(K)3, u ↦ΠBk,τKu,

by solving

 (ΠBk,τKu−u,r)K =0∀r∈Nk−2(K), (20a) ⟨(ΠBk,τKu−u)⋅n+τK(Πk−1p−p),μ⟩∂K =0∀μ∈Rk(∂K), (20b)

where satisfying either or , and is the Nédélec space with .

###### Proposition 3.1.

For the projection , we have

 h1/2K∥ΠBk,τKu−u∥∂K+∥ΠBk,τKu−u∥K≤C(hsK|u|s,K+τmaxKhtK|p|t,K), (21)

where , , is the largest value of , and depends only on and the shape-regularity of .

See appendix. ∎

## 4 Unified error analysis

In this section, we specify those approximations spaces and stabilization functions in the general setting proposed in Section 2 and consider four variants (see Table 2 for an overview). Depending on the choices of the approximations spaces and the types of meshes, we construct different projections. All these projections satisfy Assumption 2.1 and therefore we can easily obtain the error estimates of and by using Proposition 2.1 and Proposition 2.2. We prove that all the variants are optimal in Theorem 4.1.

### 4.1 Four variants and the corresponding projections

We first introduce some notation. Let be the collection of all simplex elements in and let be those non-simplex elements. We denote by with as the collections of the boundaries of the simplex and the non-simplex elements respectively.

Variant : . To the best of our knowledge, this variant has not been considered before. We require .

Let be two fixed constants. For each , we choose the stabilization functions such that and . We choose the projection

 ΠK(w,u,p) :=(Πkw,ΠBk+1,τnu,Πkp), (22a) Π∗K(w∗,u∗,p∗) :=(Πkw∗,ΠBk+1,−τnu∗,Πkp∗). (22b)

For each , we choose the stabilization functions such that and . We choose the projection

 ΠK(w,u,p) :=(Πkw,Πk+1u,Πkp), (23a) Π∗K(w∗,u∗,p∗) :=(Πkw∗,Πk+1u∗,Πkp∗). (23b)

It is easy to verify that the projections