An hp-Hybrid High-Order Method for Variable Diffusion on General Meshes

Joubine Aghili; Daniele A. Di Pietro; Berardo Ruffini

doi:10.1515/cmam-2017-0009

Publicly Available Published by De Gruyter June 17, 2017

An hp-Hybrid High-Order Method for Variable Diffusion on General Meshes

Joubine Aghili , Daniele A. Di Pietro and Berardo Ruffini

From the journal Computational Methods in Applied Mathematics

https://doi.org/10.1515/cmam-2017-0009

Abstract

In this work, we introduce and analyze an hp-hybrid high-order (hp-HHO) method for a variable diffusion problem. The proposed method is valid in arbitrary space dimension and for fairly general polytopal meshes. Variable approximation degrees are also supported. We prove hp-convergence estimates for both the energy- and L2-norms of the error, which are the first of this kind for Hybrid High-Order methods. These results hinge on a novel hp-approximation lemma valid for general polytopal elements in arbitrary space dimension. The estimates are additionally fully robust with respect to the heterogeneity of the diffusion coefficient, and show only a mild dependence on the square root of the local anisotropy, improving previous results for HHO methods. The expected exponential convergence behavior is numerically demonstrated on a variety of meshes for both isotropic and strongly anisotropic diffusion problems.

Keywords: Hybrid High-Order Methods; Discontinuous Skeletal Methods; Polytopal Methods; Variable Diffusion

MSC 2010: 65N30; 65N08; 76R50

1 Introduction

Over the last few years, discretization technologies have appeared that support arbitrary approximation orders on general polytopal meshes. In this work, we focus on a particular instance of such technologies, the so-called Hybrid High-Order (HHO) methods originally introduced in [18, 17]; see also [21] for an introduction. So far, the literature on HHO methods has focused on the h-version of the method with uniform polynomial degree. Our goal is to provide a first example of variable-degree hp-HHO method and carry out a full hp-convergence analysis valid for fairly general meshes and arbitrary space dimension. Let Ω⊂ℝd, d≥1, denote a bounded connected polytopal domain with Lipschitz boundary. We consider the variable diffusion model problem: Find u:Ω→ℝ such that

(1.1){-∇⋅(𝜿⁢∇⁡u)=fin ⁢Ω,u=0on ⁢∂⁡Ω,

where 𝜿 is a uniformly positive, symmetric, tensor-valued field on Ω, while f∈L2⁢(Ω) denotes a volumetric source. For the sake of simplicity, we assume that 𝜿 is piecewise constant on a partition PΩ of Ω into polytopes. The weak formulation of problem (1.1) reads: Find u∈U:=H01⁢(Ω) such that

(1.2)(𝜿⁢∇⁡u,∇⁡v)=(f,v) for all ⁢v∈U,

where we have used the notation (⋅,⋅) for the usual inner products of both L2⁢(Ω) and L2⁢(Ω)d. Here, the scalar-valued field u represents a potential, and the vector-valued field 𝜿⁢∇⁡u the corresponding flux.

For a given polytopal mesh 𝒯h={T} of Ω, the hp-HHO discretization of problem (1.2) is based on two sets of degrees of freedom (DOFs):

skeletal DOFs, consisting in (d-1)-variate polynomials of total degree pF≥0 on each mesh face F, and
elemental DOFs, consisting in d-variate polynomials of degree pT on each mesh element T, where pT denotes the lowest degree of skeletal DOFs on the boundary of T.

Skeletal DOFs are globally coupled, and can be alternatively interpreted as traces of the potential on the mesh faces or as Lagrange multipliers enforcing the continuity of the normal flux at the discrete level; cf. [1, 12] for further insight. Elemental DOFs, on the other hand, are auxiliary DOFs that can be locally eliminated by static condensation, as detailed in [12, Section 2.4] for the case where pF=p for all mesh faces F.

Two key ingredients are devised locally from skeletal and elemental DOFs attached to each mesh element T:

a reconstruction of the potential of degree (pT+1) (i.e., one degree higher than elemental DOFs in T) obtained solving a small Neumann problem, and
a stabilization term penalizing face-based residuals and polynomially consistent up to degree (pT+1).

The local contributions obtained from these two ingredients are then assembled following a standard, finite element-like procedure. The resulting discretization has several appealing features, the most prominent of which are summarized hereafter:

It is valid for fairly general polytopal meshes.
The construction is dimension-independent, which can significantly ease the practical implementation.
It enables the local adaptation of the approximation order, a highly desirable feature when combined with a regularity estimator (whose development will be addressed in a separate work).
It exhibits only a moderate dependence on the diffusion coefficient 𝜿.
It has a moderate computational cost thanks to the possibility of eliminating elemental DOFs locally via static condensation.
Parallel implementations can be simplified by the fact that processes communicate via skeletal unknowns only.

The seminal works on the p- and hp-conforming finite element method on standard meshes date back to the early 80s; cf. [7, 6, 5]. Starting from the late 90s, nonconforming methods on standard meshes supporting arbitrary-order have received a fair amount of attention; a (by far) nonexhaustive list of contributions focusing on scalar diffusive problems similar to the one considered here includes [29, 10, 28, 24, 31]. The possibility of refining both in h and in p on general meshes is, on the other hand, a much more recent research topic. We cite, in particular, hp-composite [3, 25] and polyhedral [9] discontinuous Galerkin methods, and the two-dimensional virtual element method proposed in [8].

The main novelty of the paper are hp-energy- and L2-estimates for HHO methods summarized in Section 3.2. These are the first results of this kind for HHO methods, and among the first for discontinuous skeletal methods in general (a prominent example of discontinuous skeletal methods are the hybridizable discontinuous Galerkin methods of [13]; cf. [12] for a precise study of their relation with HHO methods).

The cornerstone of the analysis is a novel extension of the classical Babuška–Suri hp-approximation results to regular mesh sequences in the sense of [16, Chapter 1] and arbitrary space dimension d≥1; see Lemma 2.3 below. Similar results had been derived in [8] for d=2 and, under assumptions on the mesh different from the ones proposed here, in [9] for d∈{2,3}. A relevant difference with respect to classical approximation results based on the Bramble–Hilbert lemma (see, e.g., [23, Theorem 1.103]) is that we cannot rely here on the notion of reference element, since general polyhedral shapes are allowed. As a result, it turns out that a key point consists in showing that the regularity assumptions adapted to polyhedral meshes imply uniform bounds for the Lipschitz constant of mesh elements.

The resulting energy-norm estimate confirms the characteristic h-superconvergence behavior of HHO methods, whereas we have a more standard scaling as (pT+1)-pT with respect to the polynomial degree pT of elemental DOFs. This scaling is analogous to the best available results for discontinuous Galerkin (dG) methods on rectangular meshes based on polynomials of degree pT, cf. [24] (on more general meshes, the scaling for the symmetric interior penalty dG method is pT-(pT-1/2), half a power less than for the hp-HHO method studied here). Classically, when elliptic regularity holds, the h-convergence order can be increased by 1 for the L2-norm. In our error estimates, the dependence on the diffusion coefficient is carefully tracked, showing full robustness with respect to its heterogeneity and only a moderate dependence with a power of 12 on its local anisotropy when the error in the energy-norm is considered. These results improve the ones in [19] where, considering a different norm for the error, a power of 1 was found. Numerical experiments confirm the expected exponentially convergent behavior for both isotropic and strongly anisotropic diffusion coefficients on a variety of two-dimensional meshes.

The rest of the paper is organized as follows. In Section 2 we introduce the notation and state the basic results required in the analysis including, in particular, Lemma 2.3 (whose proof is detailed in Appendix A). In Section 3 we formulate the hp-HHO method, state our main results, and provide some numerical examples. The proofs of the main results, preceded by the required preparatory material, are collected in Section 4.

2 Setting

In this section we introduce the notation and state the basic results required in the analysis.

2.1 Mesh and Notation

Let ℋ⊂ℝ*+ denote a countable set of meshsizes having 0 as its unique accumulation point. We consider mesh sequences (𝒯h)h∈ℋ where, for all h∈ℋ, 𝒯h={T} is a finite collection of nonempty disjoint open polytopal elements such that Ω¯=⋃T∈𝒯hT¯ and h=maxT∈𝒯h⁡hT (hT stands for the diameter of T). A hyperplanar closed connected subset F of Ω¯ is called a face if it has positive (d-1)-dimensional measure and

either there exist distinct T1,T2∈𝒯h such that F=∂⁡T1∩∂⁡T2 (and F is an interface)
or there exists T∈𝒯h such that F=∂⁡T∩∂⁡Ω (and F is a boundary face).

Interfaces are collected in ℱhi, boundary faces in ℱhb, and we let ℱh:=ℱhi∪ℱhb. It is assumed that the set of faces ℱh partitions the mesh skeleton in the sense that distinct faces have disjoint interiors and that ⋃T∈𝒯h∂⁡T=⋃F∈ℱhF. For all T∈𝒯h, the set ℱT:={F∈ℱh∣F⊂∂⁡T} collects the faces lying on the boundary of T and, for all F∈ℱT, we denote by 𝒏T⁢F the normal vector to F pointing out of T.

The following assumptions on the mesh will be kept throughout the exposition.

Assumption 2.1 (Admissible Mesh Sequence).

We assume that (𝒯h)h∈ℋ is admissible in the sense of [16, Chapter 1], i.e., for all h∈ℋ, 𝒯h admits a matching simplicial submesh 𝔗h and there exists a real number ϱ>0 (the mesh regularity parameter) independent of h such that the following conditions hold:

ϱ⁢hS≤rS for all h∈ℋ and all simplices S∈𝔗h of diameter hS and inradius rS.
ϱ⁢hT≤hS for all h∈ℋ, all T∈𝒯h, and all S∈𝔗h such that S⊂T.
Every mesh element T∈𝒯h is star-shaped with respect to every point of a ball of radius ϱ⁢hT.

Assumption 2.2 (Compliant Mesh Sequence).

We assume that the mesh sequence is compliant with the partition PΩ on which the diffusion tensor 𝜿 is piecewise constant, so that jumps only occur at interfaces and, for all T∈𝒯h, we set

𝜿T:=𝜿∈|Tℙ0(T)d×d.

In what follows, for all T∈𝒯h, κ¯T and κ¯T denote the largest and smallest eigenvalue of 𝜿T, respectively, and λ𝜿,T:=κ¯T/κ¯T the local anisotropy ratio.

2.2 Basic Results

For a subset X of Ω and an index q, Hq⁢(X) denotes the Hilbert space of functions which are in L2⁢(X) together with their weak derivatives of order ≤q, equipped with the usual inner product (⋅,⋅)q,X and associated norm ∥⋅∥q,X. When q=0, we recover the Lebesgue space L2⁢(X), and the subscript 0 is omitted from both the inner product and the norm. The subscript X is also omitted when X=Ω. For a given integer l≥0, we denote by ℙl⁢(X) the space spanned by the restriction to X of d-variate polynomials of degree ≤l. For further use, we also introduce the L2-projector πXl:L1⁢(X)→ℙl⁢(X) such that, for all w∈L1⁢(X),

(2.1)(πXl⁢w-w,v)X=0 for all ⁢v∈ℙl⁢(X).

We recall hereafter a few known results on admissible mesh sequences and refer to [16, Chapter 1] and [14, 15] for a more comprehensive collection. By [16, Lemma 1.41], there exists an integer N∂≥(d+1) (possibly depending on d and ϱ) such that the maximum number of faces of one mesh element is bounded,

(2.2)maxh∈ℋ,T∈𝒯h⁢card⁡(ℱT)≤N∂.

The following multiplicative trace inequality, valid for all h∈ℋ, all T∈𝒯h, and all v∈H1⁢(T), is proved in [16, Lemma 1.49]:

(2.3)∥v∥∂⁡T2≤C⁢(∥v∥T⁢∥∇⁡v∥T+hT-1⁢∥v∥T2),

where C only depends on d and ϱ. We also note the following local Poincaré’s inequality valid for all T∈𝒯h and all v∈H1⁢(T) such that (v,1)T=0:

(2.4)∥v∥T≤CP⁢hT⁢∥∇⁡v∥T,

where CP=π-1 when T is convex, while it can be estimated in terms of ϱ for nonconvex elements (cf., e.g., [34]).

The following functional analysis results lie at the heart of the hp-analysis carried out in Section 4.

Lemma 2.3 (Approximation).

There is a real number C>0 (possibly depending on d and ϱ) such that, for all h∈H, all T∈Th, all integers l≥1, all s≥0, and all v∈Hs+1⁢(T), there exists a polynomial ΠTl⁢v∈Pl⁢(T) satisfying, for all 0≤q≤s+1,

(2.5)∥v-ΠTl⁢v∥q,T≤C⁢hTmin⁡(l,s)-q+1ls+1-q⁢∥v∥s+1,T.

Proof.

See Appendix A. ∎

Lemma 2.4 (Discrete Trace Inequality).

There is a real number C>0 (possibly depending on d and ϱ) such that, for all h∈H, all T∈Th, all integers l≥1, and all v∈Pl⁢(T), it holds

(2.6)∥v∥∂⁡T≤C⁢lhT1/2⁢∥v∥T.

Proof.

When all meshes in the sequence (𝒯h)h∈ℋ are simplicial and conforming, the proof of (2.6) can be found in [30, Theorem 4.76] for d=2; for d≥2 the proof is analogous. The extension to admissible mesh sequences in the sense of Assumption 2.1 can be done following the reasoning in [16, Lemma 1.46]. ∎

3 Discretization

In this section, we formulate the hp-HHO method, state our main results, and provide some numerical examples.

3.1 The hp-HHO Method

We present in this section an extension of the classical HHO method of [18] accounting for variable polynomial degrees. Let a vector p¯h=(pF)F∈ℱh∈ℕℱh of skeletal polynomial degrees be given. For all T∈𝒯h, we denote by p¯T=(pF)F∈ℱT the restriction of p¯h to ℱT, and we introduce the following local space of DOFs:

(3.1)U¯Tp¯T:=ℙpT(T)×(×F∈ℱTℙpF(F)),pT:=minF∈ℱTpF.

The standard HHO notation v¯T=(vT,(vF)F∈ℱT) is used for a generic element of U¯Tp¯T. We define the local potential reconstruction operator

rTpT+1:U¯Tp¯T→ℙpT+1⁢(T)

such that, for all v¯T∈U¯Tp¯T and all w∈ℙpT+1⁢(T),

(3.2)(𝜿T⁢∇⁡rTpT+1⁢v¯T,∇⁡w)T=-(vT,∇⋅(𝜿T⁢∇⁡w))T+∑F∈ℱT(vF,𝜿T⁢∇⁡w⋅𝒏T⁢F)F

and

(rTpT+1⁢v¯T-vT,1)T=0.

Note that computing rTpT+1⁢v¯T according to (3.2) requires to invert the 𝜿T-weighted stiffness matrix of ℙk+1⁢(T), which can be efficiently accomplished by a Cholesky solver.

We define on U¯Tp¯T×U¯Tp¯T the local bilinear form aT such that

(3.3)aT⁢(u¯T,v¯T):=(𝜿T⁢∇⁡rTpT+1⁢u¯T,∇⁡rTpT+1⁢v¯T)T+sT⁢(u¯T,v¯T),

(3.4)sT⁢(u¯T,v¯T):=∑F∈ℱTκFhT⁢(δT⁢Fp¯T⁢u¯T,δT⁢Fp¯T⁢v¯T)F,

where, for all F∈ℱT, we let κF:=𝜿T⁢𝒏T⁢F⋅𝒏T⁢F and the face-based residual operator

δT⁢Fp¯T:U¯Tp¯T→ℙpF⁢(F)

is such that, for all v¯T∈U¯Tp¯T,

(3.5)δT⁢Fp¯T⁢v¯T:=πFpF⁢(vF-rTpT+1⁢v¯T+πTpT⁢rTpT+1⁢v¯T-vT).

The first contribution in aT is in charge of consistency, whereas the second ensures stability by a least-square penalty of the face-based residuals δT⁢Fp¯T. This subtle form for δT⁢Fp¯T ensures that the residual vanishes when its argument is the interpolate of a function in ℙpT+1⁢(T), and is required for high-order h-convergence (a detailed motivation is provided in [18, Remark 6]).

The global space of DOFs and its subspace with strongly enforced boundary conditions are defined, respectively, as

U¯hp¯h:=(×T∈𝒯hℙpT(T))×(×F∈ℱhℙpF(F)),U¯h,0p¯h:={v¯h∈U¯hp¯h∣vF≡0 for all F∈ℱhb}.

Note that interface DOFs in U¯hp¯h are single-valued, meaning that they match from one element to the adjacent one. We use the notation

v¯h=((vT)T∈𝒯h,(vF)F∈ℱh)

for a generic DOF vector in U¯hp¯h and, for all T∈𝒯h, we denote by v¯T∈U¯Tp¯T its restriction to T. For further use, we also introduce the global interpolator

I¯hp¯h:H1⁢(Ω)→U¯hp¯h

such that, for all v∈H1⁢(Ω),

(3.6)I¯hp¯h⁢v:=((πTpT⁢v)T∈𝒯h,(πFpF⁢v)F∈ℱh),

and denote by I¯Tp¯T its restriction to T∈𝒯h.

The hp-HHO discretization of problem (1.2) consists in seeking u¯h∈U¯h,0p¯h such that

(3.7)ah⁢(u¯h,v¯h)=lh⁢(v¯h) for all ⁢v¯h∈U¯h,0p¯h,

where the global bilinear form ah on U¯hp¯h×U¯hp¯h and the linear form lh on U¯hp¯h are assembled element-wise setting

ah⁢(u¯h,v¯h):=∑T∈𝒯haT⁢(u¯T,v¯T),lh⁢(v¯h):=∑T∈𝒯h(f,vT)T.

Remark 3.1 (Static Condensation).

Using a standard static condensation procedure, it is possible to eliminate element-based DOFs locally and solve (3.7) by inverting a system in the skeletal unknowns only. For the sake of conciseness, we do not repeat the details here and refer instead to [12, Section 2.4]. Accounting for the strong enforcement of boundary conditions, the size of the system after static condensation is

(3.8)Ndof=∑F∈ℱhi(pF+d-1pF).

Remark 3.2 (Finite Element Interpretation).

A finite element interpretation of the scheme (3.7) is possible following the extension proposed in [12, Remark 3] of the ideas originally developed in [4] in the context of nonconforming virtual element methods. For all F∈ℱhi, we denote by [⋅]F the usual jump operator (the sign is irrelevant), which we extend to boundary faces F∈ℱhb setting [φ]F:=φ. Let

𝔘h,0p¯h:={𝔳h∈L2(Ω)∣𝔳h∈|T𝔘Tp¯T for all T∈𝒯h and πFpF([𝔳T]F)=0 for all F∈ℱh},

where, for all T∈𝒯h, we have introduced the local space

𝔘Tp¯T:={𝔳T∈H1(T)∣∇⋅(𝜿T∇𝔳T)∈ℙpT(T) and 𝜿T∇𝔳T⋅|F𝒏T⁢F∈ℙpF(F) for all F∈ℱT}.

It can be proved that, for all T∈𝒯h,

I¯Tp¯T:𝔘Tp¯T→U¯Tp¯T

is an isomorphism. Thus, the triplet (T,𝔘Tp¯T,I¯Tp¯T) defines a finite element in the sense of Ciarlet [11]. Additionally, problem (3.7) can be reformulated as the nonconforming finite element method: Find 𝔲h∈𝔘h,0p¯h such that

𝔞h⁢(𝔲h,𝔳h)=lh⁢(𝔳h) for all ⁢𝔳h∈𝔘h,0p¯h,

where

𝔞h⁢(𝔲h,𝔳h):=ah⁢(I¯hp¯h⁢𝔲h,I¯hp¯h⁢𝔳h),lh⁢(𝔳h):=lh⁢(I¯hp¯h⁢𝔳h),

and it can be proved that 𝔲h is the unique element of 𝔘h,0p¯h such that u¯h=I¯hp¯h⁢𝔲h with u¯h unique solution to (3.7).

3.2 Main Results

We next state our main results. The proofs are postponed to Section 4. For all T∈𝒯h, we denote by ∥⋅∥a,T and |⋅|s,T the seminorms defined on U¯Tp¯T by the bilinear forms aT and sT, respectively, and by ∥⋅∥a,h the seminorm defined by the bilinear form ah on U¯hp¯h. We also introduce the penalty seminorm |⋅|s,h such that, for all v¯h∈U¯hp¯h,

|v¯h|s,h2:=∑T∈𝒯h|v¯T|s,T2.

Note that ∥⋅∥a,h is a norm on the subspace U¯h,0p¯h with strongly enforced boundary conditions (the arguments are analogous to those of [17, Proposition 5]). We will also need the global reconstruction operator

rhp¯h:U¯hp¯h→L2⁢(Ω)

such that, for all v¯h∈U¯hp¯h,

(rhp¯hv¯h)=|TrTpT+1v¯T for all T∈𝒯h.

Finally, for the sake of conciseness, throughout the rest of the paper we denote by a≲b the inequality a≤C⁢b with real number C>0 independent of h, p¯h, 𝜿 and, for local inequalities on a mesh element T∈𝒯h, also of T.

Our first estimate concerns the error measured in energy-like norms.

Theorem 3.3 (Energy Error Estimate).

Let u∈U and u¯h∈U¯h,0p¯h denote the unique solutions of problems (1.2) and (3.7), respectively, and set

(3.9)u¯^h:=I¯hp¯h⁢u.

Assuming that u∈|THpT+2(T) for all T∈Th, it holds

(3.10)∥u¯h-u¯^h∥a,h≲(∑T∈𝒯hκ¯T⁢λ𝜿,T⁢hT2⁢(pT+1)(pT+1)2⁢pT⁢∥u∥pT+2,T2)1/2.

Consequently, denoting by ∇h the broken gradient on Th, we have that

(3.11)∥𝜿1/2⁢∇h⁡(u-rhp¯h⁢u¯h)∥2+|u¯h|s,h2≲∑T∈𝒯hκ¯T⁢λ𝜿,T⁢hT2⁢(pT+1)(pT+1)2⁢pT⁢∥u∥pT+2,T2.

Proof.

See Section 4.3. ∎

In (3.10) and (3.11), we observe the characteristic improved h-convergence of HHO methods (cf. [18]), whereas, in terms of p-convergence, we have a more standard scaling as (pT+1)-pT (i.e., half a power more than discontinuous Galerkin methods based on polynomials of degree pT, cf., e.g., [28]). In (3.11), we observe that the left-hand side has the same convergence rate (both in h and in p) as the interpolation error

∥𝜿1/2⁢∇h⁡(u-rhp¯h⁢u¯^h)∥2+|u¯^h|s,h2,

as can be verified combining (4.2) and (4.4) below. Note that, in this case, the p-convergence is limited by the second term, which measures the discontinuity of the potential reconstruction at interfaces. An inspection of formulas (3.10) and (3.11) also shows that the method is fully robust with respect to the heterogeneity of the diffusion coefficient, while only a moderate dependence (with a power of 12) is observed with respect to its local anisotropy ratio.

For the sake of completeness, we also provide an estimate of the L2-error between the piecewise polynomial fields uh and u^h such that

uh:=|TuT and u^h:=|Tu^T=πTpTu for all T∈𝒯h.

To this end, we need elliptic regularity in the following form: For all g∈L2⁢(Ω), the unique element z∈U such that

(3.12)(𝜿⁢∇⁡z,∇⁡v)=(g,v) for all ⁢v∈U,

satisfies the a priori estimate

(3.13)∥z∥2≲κ¯-1⁢∥g∥L2⁢(Ω),κ¯:=minT∈𝒯h⁡κ¯T.

The following result is proved in Section 4.4.

Figure 1

Convergence with p-refinement (regular test case).

$(a) ∥u¯h-u¯^h∥a,h\lVert\underline{u}_{h}-\widehat{\underline{u}}_{h}\rVert_{\mathrm{a},h} vs. NdofN_{\mathrm{dof}}$

(a)

∥u¯h-u¯^h∥a,h vs. Ndof

$(b) ∥uh-u^h∥\lVert u_{h}-\widehat{u}_{h}\rVert vs. NdofN_{\mathrm{dof}}$

(b)

∥uh-u^h∥ vs. Ndof

$(c) ∥u¯h-u¯^h∥a,h\lVert\underline{u}_{h}-\widehat{\underline{u}}_{h}\rVert_{\mathrm{a},h} vs. p$

(c)

∥u¯h-u¯^h∥a,h vs. p

$(d) ∥uh-u^h∥\lVert u_{h}-\widehat{u}_{h}\rVert vs. p$

(d)

∥uh-u^h∥ vs. p

Figure 2

Convergence with p-refinement (Le Potier’s test case).

$(a) ∥u¯h-u¯^h∥a,h\lVert\underline{u}_{h}-\widehat{\underline{u}}_{h}\rVert_{\mathrm{a},h} vs. NdofN_{\mathrm{dof}}$

(a)

∥u¯h-u¯^h∥a,h vs. Ndof

$(b) ∥uh-u^h∥\lVert u_{h}-\widehat{u}_{h}\rVert vs. NdofN_{\mathrm{dof}}$

(b)

∥uh-u^h∥ vs. Ndof

$(c) ∥u¯h-u¯^h∥a,h\lVert\underline{u}_{h}-\widehat{\underline{u}}_{h}\rVert_{\mathrm{a},h} vs. p$

(c)

∥u¯h-u¯^h∥a,h vs. p

$(d) ∥uh-u^h∥\lVert u_{h}-\widehat{u}_{h}\rVert vs. p$

(d)

∥uh-u^h∥ vs. p

Theorem 3.4 (L2-Error Estimate).

Under the assumptions of Theorem 3.3, and further assuming elliptic regularity (3.13) and that f∈HpT+ΔT⁢(T) for all T∈Th with ΔT=1 if pT=0 while ΔT=0 otherwise, it holds

κ¯⁢∥uh-u^h∥≲κ¯1/2⁢λ𝜿⁢h⁢(∑T∈𝒯hλ𝜿,T⁢κ¯T⁢hT2⁢(pT+1)(pT+1)2⁢pT⁢∥u∥pT+2,T2)1/2+(∑T∈𝒯hhT2⁢(pT+2)(pT+ΔT)2⁢(pT+2)⁢∥f∥pT+ΔT2)1/2

with λ𝛋:=maxT∈Th⁡λ𝛋,T and κ¯:=maxT∈Th⁡κ¯T.

3.3 Numerical Examples

We close this section with some numerical examples. The h-convergence properties of the method (3.7) have been numerically investigated in [18, Section 4]. To illustrate its p-convergence properties, we solve on the unit square domain Ω=(0,1)2 the homogeneous Dirichlet problem with exact solution u=sin⁡(π⁢x1)⁢sin⁡(π⁢x2) and right-hand side f chosen accordingly. We consider two values for the diffusion coefficients:

𝜿1=𝑰2,𝜿2=((x2-x¯2)2+ϵ⁢(x1-x¯1)2-(1-ϵ)⁢(x1-x¯1)⁢(x2-x¯2)-(1-ϵ)⁢(x1-x¯1)⁢(x2-x¯2)(x1-x¯1)2+ϵ⁢(x2-x¯2)2),

where 𝑰2 denotes the identity matrix of dimension 2, 𝒙¯:=-(0.1,0.1), and ϵ=1⋅10-2. The choice 𝜿=𝜿1 (“regular” test case) yields a homogeneous isotropic problem, while the choice 𝜿=𝜿2 (“Le Potier’s” test case [27]) corresponds to a highly anisotropic problem where the principal axes of the diffusion tensor vary at each point of the domain. Figures 1 and 2 depict the energy- and L2-errors as a function of the number of skeletal DOFs Ndof (cf. (3.8)) when pF=p for all F∈ℱh and p∈{0,…,9} for the proposed choices for 𝜿 on the meshes of Figure 3. In all the cases, the expected exponentially convergent behavior is observed. Interestingly, the best performance in terms of error vs. Ndof is obtained for the Cartesian and Voronoi meshes. A comparison of the results for the two values of the diffusion coefficients allows to appreciate the robustness of the method with respect to anisotropy.

Figure 3

Meshes considered in the p-convergence test of Section 3.3. The triangular, Cartesian, refined, and staggered meshes originate from the FVCA5 benchmark [26]; the hexagonal mesh was originally introduced in [20]; the Voronoi mesh was obtained using the PolyMesher algorithm of [33].

(a)

Triangular

(b)

Cartesian

(c)

Refined

(d)

Staggered

(e)

Hexagonal

(f)

Voronoi

4 Convergence Analysis

In this section we prove the results stated in Section 3.2.

4.1 Consistency of the Potential Reconstruction

Preliminary to the convergence analysis is the study of the approximation properties of the potential reconstruction rTpT+1 defined by (3.2) when its argument is the interpolate of a regular function. Let a mesh element T∈𝒯h be fixed. For any integer l≥1, we define the oblique elliptic projector ϖ𝜿,Tl:H1⁢(T)→ℙl⁢(T) such that, for all v∈H1⁢(T), (ϖ𝜿,Tl⁢v-v,1)T=0 and it holds

(4.1)(𝜿T⁢∇⁡(ϖ𝜿,Tl⁢v-v),∇⁡w)T=0 for all ⁢w∈ℙl⁢(T).

Proposition 4.1 (Characterization of rTpT+1∘I¯Tp¯T).

It holds, for all T∈Th,

rTpT+1∘I¯Tp¯T=ϖ𝜿,TpT+1.

Proof.

For a generic v∈H1⁢(T), letting v¯T=I¯Tp¯T⁢v in (3.2), we infer, for all w∈ℙpT+1⁢(T),

(𝜿T⁢∇⁡(rTpT+1∘I¯Tp¯T)⁡v,∇⁡w)T=-(πTpT⁢v,∇⋅(𝜿T⁢∇⁡w))T+∑F∈ℱT(πFpF⁢v,𝜿T⁢∇⁡w⋅𝒏T⁢F)F

=-(v,∇⋅(𝜿T⁢∇⁡w))T+∑F∈ℱT(v,𝜿T⁢∇⁡w⋅𝒏T⁢F)F

=(𝜿T⁢∇⁡v,∇⁡w)T,

where we have used the fact that ∇⋅(𝜿T⁢∇⁡w)∈ℙpT-1⁢(T)⊂ℙpT⁢(T) and (𝜿T∇w)⋅|F𝒏T⁢F∈ℙpT(F)⊂ℙpF(F) (cf. definition (3.1) of pT) to cancel the projectors in the second line, and an integration by parts to conclude. ∎

We next study the approximation properties of ϖ𝜿,Tl, from which those of (rTpT+1∘I¯Tp¯T) follow in the light of Proposition 4.1.

Lemma 4.2 (Approximation Properties of ϖκ,Tl).

For all integers l≥1, all mesh elements T∈Th, all 0≤s≤l, and all v∈Hs+1⁢(T), it holds

∥𝜿T1/2⁢∇⁡(v-ϖ𝜿,Tl⁢v)∥T+hT1/2l⁢∥𝜿T1/2⁢∇⁡(v-ϖ𝜿,Tl⁢v)∥∂⁡T+κ¯T1/2hT⁢∥v-ϖ𝜿,Tl⁢v∥T+κ¯T1/2hT1/2⁢∥v-ϖ𝜿,Tl⁢v∥∂⁡T

(4.2)≲κ¯T1/2⁢hTmin⁡(l,s)ls⁢∥v∥s+1,T.

Proof.

By definition (4.1) of ϖ𝜿,Tl, it holds

(4.3)∥𝜿T1/2∇(v-ϖ𝜿,Tlv)∥T2=minw∈ℙl⁢(T)∥𝜿T1/2∇(v-w)∥T2≤κ¯T∥∇(v-ΠTlv)∥T2.

Hence, using (2.5) with q=1, it is readily inferred that

∥𝜿T1/2⁢∇⁡(v-ϖ𝜿,Tl⁢v)∥T≲κ¯T1/2⁢hTmin⁡(l,s)ls⁢∥v∥s+1,T.

To prove the second bound in (4.2), use the triangle inequality to infer

∥𝜿T1/2⁢∇⁡(v-ϖ𝜿,Tl⁢v)∥∂⁡T≤∥𝜿T1/2⁢∇⁡(v-ΠTl⁢v)∥∂⁡T+∥𝜿T1/2⁢∇⁡(ΠTl⁢v-ϖ𝜿,Tl⁢v)∥∂⁡T=:𝔗1+𝔗2.

For the first term, the multiplicative trace inequality (2.3) combined with (2.5) (with q=1,2) gives

𝔗1≲κ¯T1/2⁢hTmin⁡(l,s)-1/2ls-1/2⁢∥v∥s+1,T.

For the second term, we have

𝔗2≲lhT1/2⁢∥𝜿T1/2⁢∇⁡(ΠTl⁢v-ϖ𝜿,Tl⁢v)∥T

≤lhT1/2⁢(∥𝜿T1/2⁢∇⁡(ΠTl⁢v-v)∥T+∥𝜿T1/2⁢∇⁡(v-ϖ𝜿,Tl⁢v)∥T)

≲κ¯T1/2⁢lhT1/2⁢∥∇⁡(ΠTl⁢v-v)∥T

≲κ¯T1/2⁢hTmin⁡(l,s)-1/2ls-1⁢∥v∥s+1,T,

where we have used the discrete trace inequality (2.6) in the first line, the triangle inequality in the second line, the estimate (4.3) in the third, and the approximation result (2.5) with q=1 to conclude. To obtain the third bound in (4.2), after recalling that (v-ϖ𝜿,Tl⁢v,1)T=0, we apply the local Poincaré inequality (2.4) to infer

∥v-ϖ𝜿,Tl⁢v∥T≲hTκ¯T1/2⁢∥𝜿T1/2⁢∇⁡(v-ϖ𝜿,Tl⁢v)∥T≲κ¯T1/2κ¯T1/2⁢hTmin⁡(l,s)+1ls⁢∥v∥s+1,T,

where the conclusion follows from the first bound in (4.2). Finally, to obtain the last bound, we use the multiplicative trace inequality (2.3) to infer

∥v-ϖ𝜿,Tl⁢v∥∂⁡T2≲κ¯T-1/2⁢∥v-ϖ𝜿,Tl⁢v∥T⁢∥𝜿T1/2⁢∇⁡(v-ϖ𝜿,Tl⁢v)∥T+hT-1⁢∥v-ϖ𝜿,Tl⁢v∥T2,

and use the first and third bound in (4.2) to estimate the various terms. ∎

4.2 Consistency of the Stabilization Term

The consistency properties of the stabilization bilinear form sT defined by (3.4) are summarized in the following lemma.

Lemma 4.3 (Consistency of the Stabilization Term).

For all T∈Th, all 0≤q≤pT, and all v∈Hq+2⁢(T), it holds

(4.4)|I¯Tp¯T⁢v|s,T≲κ¯T1/2⁢λ𝜿,T1/2⁢hTmin⁡(pT,q)+1(pT+1)q⁢∥v∥q+2,T.

Proof.

Let T∈𝒯h and v∈Hq+2⁢(T) and set, for the sake of brevity,

vˇT:=(rTpT+1∘I¯Tp¯T)⁢v=ϖ𝜿,TpT+1⁢v

(cf. Proposition 4.1). For all F∈ℱT, recalling the definitions of the face residual δT⁢Fp¯T (cf. (3.5)) and of the local interpolator I¯Tp¯T (cf. (3.6)), together with the fact that pT≤pF by definition (3.1), we get

δT⁢Fp¯T⁢I¯Tp¯T⁢v=πFpF⁢(vˇT-v)-πTpT⁢(vˇT-v).

Using the triangle inequality and the L2⁢(F)-stability of πFpF, we infer,

(4.5)∥δT⁢Fp¯TI¯Tp¯Tv∥F≤∥vˇT-v∥F+∥πTpT(vˇT-v)∥F=:𝔗1+𝔗2.

For the first term, the approximation properties (4.2) of ϖ𝜿,TpT+1 (with l=pT+1 and s=q+1) readily yield

(4.6)𝔗1≲λ𝜿,T1/2⁢hTmin⁡(pT,q)+3/2(pT+1)q⁢∥v∥q+2,T.

For the second term, on the other hand, the discrete trace inequality (2.6) followed by the L2⁢(T)-stability of πTpT and (4.2) (with l=pT+1 and s=q+1) gives

(4.7)𝔗2≲(pT+1)hT1/2⁢∥πTpT⁢(vˇT-v)∥T≤(pT+1)hT1/2⁢∥vˇT-v∥T≲λ𝜿,T1/2⁢hTmin⁡(pT,q)+3/2(pT+1)q⁢∥v∥q+2,T.

The bound (4.4) follows using (4.6)–(4.7) in the right-hand side of (4.5), squaring the resulting inequality, multiplying it by κF/hT, summing over F∈ℱT, and using the bound (2.2) on card⁡(ℱT). ∎

4.3 Energy Error Estimate

Proof of Theorem 3.3.

We start by noting the following abstract error estimate:

(4.8)∥u¯h-u¯^h∥a,h≤supv¯h∈U¯h,0p¯h,∥v¯h∥a,h=1⁡Eh⁢(v¯h),

with nonconformity error

(4.9)Eh⁢(v¯h):=lh⁢(v¯h)-ah⁢(u¯^h,v¯h).

To prove (4.8), it suffices to observe that

∥u¯h-u¯^h∥a,h2=ah⁢(u¯h-u¯^h,u¯h-u¯^h)

=ah⁢(u¯h,u¯h-u¯^h)-ah⁢(u¯^h,u¯h-u¯^h)

=lh⁢(u¯h-u¯^h)-ah⁢(u¯^h,u¯h-u¯^h),

where we have used the definition of the ∥⋅∥a,h-norm in the first line, the linearity of ah in its first argument in the second line, and the discrete problem (3.7) in the third. The conclusion follows dividing both sides by ∥u¯h-u¯^h∥a,h, using linearity, and passing to the supremum.

We next bound the nonconformity error Eh⁢(v¯h) for a generic vector of DOFs v¯h∈U¯h,0p¯h. A preliminary step consists in finding a more appropriate rewriting for Eh⁢(v¯h). Observing that f=-∇⋅(𝜿⁢∇⁡u) a.e. in Ω, integrating by parts element-by-element, and using the continuity of the normal component of 𝜿⁢∇⁡u across interfaces together with the strongly enforced boundary conditions in U¯h,0p¯h to insert vF into the second term in parentheses, we infer

(4.10)lh⁢(v¯h)=∑T∈𝒯h((𝜿T⁢∇⁡u,∇⁡vT)T+∑F∈ℱT(𝜿T⁢∇⁡u⋅𝒏T⁢F,vF-vT)F).

Setting, for the sake of conciseness (cf. Proposition 4.1),

(4.11)uˇT:=rTpT+1⁢u¯^T=ϖ𝜿,TpT+1⁢u,

and using the definition (3.2) of rTpT+1⁢v¯T with w=uˇT, we have

(4.12)ah⁢(u¯^h,v¯h)=∑T∈𝒯h((𝜿T⁢∇⁡uˇT,∇⁡vT)T+∑F∈ℱT(𝜿T⁢∇⁡uˇT⋅𝒏T⁢F,vF-vT)F+sT⁢(u¯^T,v¯T)).

Subtracting (4.12) from (4.10), and observing that the first terms inside the summations cancel out owing to (4.1), we have

(4.13)Eh⁢(v¯h)=∑T∈𝒯h(∑F∈ℱT(𝜿T⁢∇⁡(uˇT-u)⋅𝒏T⁢F,vF-vT)F+sT⁢(u¯^T,v¯T)).

Denote by 𝔗1⁢(T) and 𝔗2⁢(T) the two summands in parentheses. Using the Cauchy–Schwarz inequality followed by the approximation properties (4.2) of uˇT (with l=s=pT+1) and (4.16) below, we have for the first term

|𝔗1⁢(T)|≤hT1/2⁢∥𝜿T1/2⁢∇⁡(uˇT-u)∥∂⁡T×(∑F∈ℱTκFhT⁢∥vF-vT∥F2)1/2

(4.14)≲κ¯T1/2⁢λ𝜿,T1/2⁢hTpT+1(pT+1)pT⁢∥u∥pT+2,T⁢∥v¯T∥a,T.

For the second term, the Cauchy–Schwarz inequality followed by (4.4) (with q=pT) readily yields

(4.15)|𝔗2⁢(T)|≲κ¯T1/2⁢λ𝜿,T1/2⁢hTpT+1(pT+1)pT⁢∥u∥pT+2,T⁢|v¯T|s,T≲κ¯T1/2⁢λ𝜿,T1/2⁢hTpT+1(pT+1)pT⁢∥u∥pT+2,T⁢∥v¯T∥a,T.

Using (4.14)–(4.15) to estimate the right-hand side of (4.13), applying the Cauchy–Schwarz inequality, and passing to the supremum yields (3.10). To prove (3.11), it suffices to observe that, inserting u¯^h and using the triangle inequality, we obtain

∥𝜿1/2⁢∇h⁡(u-rhp¯h⁢u¯h)∥2+|u¯h|s,h2≲∥𝜿1/2⁢∇h⁡(u-rhp¯h⁢u¯^h)∥2+|u¯^h|s,h2+∥u¯h-u¯^h∥a,h2.

and (3.11) follows using the estimates (4.2), (4.4), and (3.10) to bound the terms in the right-hand side. ∎

Remark 4.4 (Nonconformity Error).

The name nonconformity error for the quantity (4.9) is reminiscent of the Finite Element literature. As a matter of fact Eh corresponds to the contribution that appears in the estimate provided by the Second Strang Lemma when nonconforming methods are considered; see, e.g., [23, Lemma 2.25] and, in particular, the third term in the right-hand side of [23, (2.20)].

Proposition 4.5 (Estimate of Boundary Difference Seminorm).

It holds, for all v¯T∈U¯Tp¯T,

(4.16)∑F∈ℱTκFhT⁢∥vF-vT∥F2≲λ𝜿,T⁢∥v¯T∥a,T2.

Proof.

Let T∈𝒯h, v¯T∈U¯Tp¯T, and set, for the sake of brevity, vˇT:=rTpT+1⁢v¯T. We have, for all F∈ℱT,

∥vF-vT∥F=∥πFpF⁢(vF-vT)∥F

=∥πFpF⁢(vF-vˇT+πTpT⁢vˇT-vT+vˇT-πTpT⁢vˇT)∥F

(4.17)≤∥δT⁢Fp¯T⁢v¯T∥F+∥vˇT-πTpT⁢vˇT∥F,

where we have used the fact that pT≤pF (cf. (3.1)) to infer that vT∈|FℙpF(F) and thus insert πFpF in the first line, added and subtracted (vˇT-πTpT⁢vˇT) in the second line, used the triangle inequality together with the definition (3.5) of the face-based residual δT⁢Fp¯T and the L2⁢(F)-stability of πFpF in the third. To conclude, we observe that, if pT=0, then the discrete trace inequality (2.6) followed by Poincaré’s inequality yield

∥vˇT-πT0⁢vˇT∥F≲hT1/2⁢κ¯T-1/2⁢∥𝜿⁢∇⁡vˇT∥T

while, if pT≥1, then

∥vˇT-πTpT⁢vˇT∥F=∥vˇT-πT0⁢vˇT-πTpT⁢(vˇT-πT0⁢vˇT)∥F

≲pT+1hT1/2⁢∥vˇT-πT0⁢vˇT-πTpT⁢(vˇT-πT0⁢vˇT)∥T

≲pT+1hT1/2⁢hTpT⁢∥vˇT-πT0⁢vˇT∥1,T

≲hT1/2⁢κ¯T-1/2⁢∥𝜿T1/2⁢∇⁡vˇT∥T,

where we have inserted ±πT0⁢vˇT in the first line, used the discrete trace inequality (2.6) in the second line, the L2⁢(T)-optimality of πTpT together with the approximation properties (2.5) (with l=pT and q=s=0) in the third line, and we have concluded observing that (pT+1)/pT≤2 and using the local Poincaré inequality (2.4) to infer

∥vˇT-πT0⁢vˇT∥1,T≲hT1/2⁢∥∇⁡vˇT∥T.

Plugging the above bounds for ∥vˇT-πTpT⁢vˇT∥F into (4.17), squaring the resulting inequality, multiplying it by κF/hT, summing over F∈ℱT, and recalling the bound (2.2) on card⁡(ℱT), we deduce (4.16). ∎

4.4 L2-Error Estimate

Proof of Theorem 3.4.

We let z∈U solve (3.12) with g=u^h-u and set z¯^h:=I¯hp¯h⁢z and, for all T∈𝒯h (cf. Proposition 4.1),

(4.18)zˇT:=rTpT+1⁢z¯^T=ϖ𝜿,TpT+1⁢z.

For the sake of brevity, we also let e¯h:=u¯^h-u¯h∈U¯h,0p¯h (recall the definition (3.9) of u¯^h), so that u^T-uT=eT for all T∈𝒯h. We start by observing that

(4.19)∥eh∥2=-(∇⋅(𝜿⁢z),eh)=∑T∈𝒯h((𝜿T⁢∇⁡z,∇⁡eT)T+∑F∈ℱT(𝜿T⁢∇⁡z⋅𝒏T⁢F,eF-eT)F),

where we have used the fact that -∇⋅(𝜿⁢z)=eh a.e. in Ω followed by element-by-element partial integration together with the continuity of the normal component of 𝜿T⁢∇⁡z across interfaces and the strongly enforced boundary conditions in U¯h,0p¯h to insert eF into the last term.

In view of adding and subtracting ah⁢(e¯h,z¯^h) to the right-hand side of (4.19), we next provide two useful reformulations of this quantity. First, we have

ah⁢(e¯h,z¯^h)=ah⁢(u¯^h,z¯^h)-ah⁢(u¯h,z¯^h)+(f,z)-(𝜿⁢∇⁡u,∇⁡z)

=∑T∈𝒯h((𝜿T⁢∇⁡uˇT,∇⁡zˇT)T-(𝜿T⁢∇⁡u,∇⁡z)T+sT⁢(u¯^T,z¯^T)+(f,z-πTpT⁢z)T)

(4.20)=∑T∈𝒯h((𝜿T⁢∇⁡(uˇT-u),∇⁡(zˇT-z))T+sT⁢(u¯^T,z¯^T)+(f-πTpT⁢f,z-πT1-ΔT⁢z)T),

where we have added the quantity (f,z)-(𝜿⁢∇⁡u,∇⁡z)=0 (cf. (1.2)) in the first line, we have passed to the second line using the definition (3.3) of aT (with u¯T=u¯^T and v¯T=z¯^T) together with the discrete problem (3.7) to infer

ah⁢(u¯h,z¯^h)=(f,z^h)=∑T∈𝒯h(f,πTpT⁢z)T,

and we have concluded using the definitions (4.1) of ϖ𝜿,TpT+1 (together with (4.11) and (4.18)) and (2.1) of πTpT and πT1-ΔT. Second, using the definition (3.2) of rTpT+1 (with v¯T=e¯T), we obtain

(4.21)ah⁢(e¯h,z¯^h)=∑T∈𝒯h((𝜿T⁢∇⁡z,∇⁡eT)T+∑F∈ℱT(𝜿T⁢∇⁡zˇT,eF-eT)F+sT⁢(z¯^T,e¯T)),

where we have additionally used the fact that zˇT=ϖ𝜿,TpT+1⁢z (cf. (4.18)) together with the definition (4.1) of ϖ𝜿,TpT+1 to replace zˇT by z in the first term in parentheses.

Thus, adding (4.20) and subtracting (4.21) from (4.19), we obtain after rearranging

(4.22)∥eh∥2=∑T∈𝒯h(𝔗1⁢(T)+𝔗2⁢(T)+𝔗3⁢(T)),

with

𝔗1⁢(T):=∑F∈ℱT(𝜿T⁢∇⁡(z-zˇT)⋅𝒏T⁢F,eF-eT)F+sT⁢(z¯^T,e¯T),

𝔗2⁢(T):=(𝜿T⁢∇⁡(uˇT-u),∇⁡(zˇT-z))T+sT⁢(u¯^T,z¯^T),

𝔗3⁢(T):=(f-πTpT⁢f,z-πT1-ΔT⁢z)T.

Using the Cauchy–Schwarz inequality, the approximation properties (4.2) of zˇT (with l=pT+1 and s=1) together with the consistency properties (4.4) of sT (with q=0) for the first factor, and the bound (4.16) for the second factor, we get

|𝔗1⁢(T)|≤(hT⁢∥𝜿T1/2⁢∇⁡(z-zˇT)∥∂⁡T2+|z¯^T|s,T2)1/2×(∑F∈ℱTκFhT⁢∥eF-eT∥F2+|e¯T|s,T2)1/2

(4.23)≲κ¯T1/2⁢λ𝜿,T⁢hT⁢∥e¯T∥a,T⁢∥z∥2,T.

For the second term, the Cauchy–Schwarz inequality followed by the approximation properties (4.2) of uˇT (with q=pT) and zˇT (with q=1), and the consistency properties (4.4) of sT (with q=pT and q=0 for the first and second factor, respectively) yield

|𝔗2⁢(T)|≤(∥𝜿1/2⁢∇⁡(uˇT-u)∥T2+|u¯^T|s,T2)1/2×(∥𝜿1/2⁢∇⁡(zˇT-z)∥T2+|z¯^T|s,T2)1/2

(4.24)≲κ¯T⁢λ𝜿,T⁢hTpT+2(pT+1)pT⁢∥u∥pT+2,T⁢∥z∥2,T.

Finally, for the third term we have, when pT=0,

(4.25)|𝔗3⁢(T)|≤∥f-πT0⁢f∥T⁢∥z-πT0⁢z∥T≲hT2⁢∥f∥1,T⁢∥z∥1,T≤hT2⁢∥f∥1,T⁢∥z∥2,T,

while, when pT≥1,

|𝔗3⁢(T)|≤∥f-πTpT⁢f∥T⁢∥z-πT1⁢z∥T

≤∥f-ΠTpT⁢f∥T⁢∥z-ΠT1⁢z∥T

(4.26)≲hTpT+2pTpT+2⁢∥f∥pT⁢∥z∥2,T≤hTpT+2pTpT+2⁢∥f∥pT,T⁢∥z∥2,T,

where we have used the optimality of πTpT in the L2⁢(T)-norm to pass to the second line and the approximation properties (2.5) of ΠTpT to conclude. Using (4.23)–(4.26) to bound the right-hand side of (4.22), and recalling the energy error estimate (3.10) and elliptic regularity (3.13) concludes the proof. ∎

Funding source: Agence Nationale de la Recherche

Award Identifier / Grant number: project HHOMM (ANR-15-CE40-0005)

Funding statement: We acknowledge the support of Agence Nationale de la Recherche through project HHOMM (ANR-15-CE40-0005).

A Proof of Lemma 2.3

Let K^⊂ℝd be a L-Lipschitz set (that is, such that its boundary can be locally parametrized by means of L-Lipschitz functions) with diam⁢(K^)=1, and fix r0>1 and a d-cube R⁢(r0) containing K^. In [6, Lemma 3.1], the following is shown: Suppose that K^ is a cube, then given a function v∈Hs+1⁢(K^), there exists a polynomial ΠK^l⁢v∈ℙl⁢(K^) that satisfies

(A.1)∥v-ΠK^l⁢v∥q,K^≲1ls+1-q⁢∥v∥s+1,K^

for q≤s+1. Actually, by analyzing the proof of that result, we can reinforce that statement not only for cubes, but for sets K^ which are regular enough, in the sense that there has to exist an extension operator E:Hs+1⁢(K^)→Hs+1⁢(R⁢(2⁢r0)) such that

(A.2)∥E⁢(v)∥s+1,R⁢(2⁢r0)≤C⁢∥v∥s+1,R⁢(K^),E⁢(v)=0⁢ on ⁢R⁢(2⁢r0)∖R⁢(32⁢r0).

So the problem is to show the existence of such an extension (in any dimension d≥1). The minimal assumptions such that this is doable are studied in [32, Theorems 5 and 5′, p. 181]. Namely, by means of a careful inspection of [32, Theorems 5 and 5′] and their proofs, and in particular formulas (25), (30) and the end of the proof of Theorem 5 at p. 192, we get that the constant C in (A.2) depends on the Lipschitz constant L and on the (minimal) number of L-Lipschitz coverings of K^, that is, the number of open sets which cover ∂⁡K^ and in each of whom ∂⁡K^ can be parametrized by means of an L-Lipschitz function.^[1] Thus, to get the hp-estimate (2.5), we have to show two facts:

first, that (2.5) follows from (A.1) whenever we replace K^ with an element T of the mesh;
second, that (A.1) actually holds with a uniform bounding constant for the elements of the mesh.

(i) Proof of (2.5) for regular elements. Assume for the moment that any T∈𝒯h satisfies the regularity assumptions needed to have the extension (A.2), and thus (A.1), with T^:=T/hT in place of K^. We additionally assume, without loss of generality, that the barycenter of T (and thus of T^) is . Then, by homogeneity, we get that, for every f∈Hr⁢(T), letting f^⁢(𝒙):=f⁢(𝒙/λ),

(A.3)∥f∥r,T≲λd/2-r⁢∥f^∥r,Tλ.

Thus, setting r=s+1, λ=hT, and f=v-𝔮 with 𝔮 a generic polynomial of degree l, we get by (A.1) (applied to v-𝔮 and ΠTl⁢(v-𝔮) in place of v and ΠTl⁢v, respectively),

∥v-ΠTl⁢v∥q,T=∥(v-𝔮)-ΠTl⁢(v-𝔮)∥q,T≲hTd/2-qls+1-q⁢∥v^-𝔮^∥s+1,T^.

Now, [22, Theorem 3.2] states that

|f-ΠT^l⁢f|i,T^≲|f|l,T^

for any i≤l+1. Thus by the definition of the Hq-norm we conclude that

∥v-ΠTl⁢v∥q,T≲hTd/2-qls+1-q⁢(∑i=min⁡(l,s)+1s+1|v^|i,T^2)1/2≲hTmin⁡(l,s)-q+1ls+1-q⁢∥v∥s+1,T.

Notice that we exploited again (A.3) to scale back on T, and the fact that hT<1.

(ii) Proof of regularity under Assumption 2.1. To conclude the proof, we are left to show that Assumption 2.1 entails a uniform bound only in terms of ϱ for the Lipschitz constant of every element T∈𝒯h. To this aim, consider 𝒙∈∂⁡T. Then, 𝒙∈S for some (convex) element of the submesh S∈𝔗h contained in T. Since the sets S form a submesh of 𝒯h, a uniform bound on the Lipschitz regularity of ∂⁡S would imply a bound on the Lipschitz regularity of ∂⁡T. Thus, we focus on the regularity of S. Since S is convex, we can cover ∂⁡S by means of 2⁢(d+1) open sets Ui, such that ∂⁡S∩Ui admits a local convex (and thus Lipschitz) parametrization ϕi, i.e., there exists an orthogonal coordinate system such that ∂⁡S∩Ui is the graph of a Lipschitz function ϕi:Ii⊂ℝd-1→ℝ. This bound on the number of open sets Ui is crucial to get [32, Theorem 5] to work (clearly, thanks to (2.2), the bound on the number of Lipschitz coverings of T is bounded by a constant 2d⁢N∂=c⁢(d,ϱ)). We claim that each ϕi is 1/ϱ-Lipschitz.

Suppose that 𝒙∈Ui=:U and set ϕ:=ϕi. Up to a rotation and a rescaling, we can suppose that 𝒙= and ϕ⁢(𝒙)=ϕ⁢()=0. Let now rs be the inradius of S and hS be its diameter. By Assumption 2.1, we know that hsrS≤1ϱ. Let Brs be a ball contained in S of radius rS. Up to a further rotation of center 𝒙= of the coordinate system, we can suppose that BrS is centered on the xd-axis. In place of ϕ:I→ℝ, it is useful to consider its Lipschitz extension ϕ~:ℝd-1→ℝ defined by, denoting by |⋅| the usual Euclidian norm,

ϕ~⁢(𝒙):=inf⁡{ϕ⁢(𝒚)+Lip⁢(ϕ)⁢|𝒚-𝒙|∣𝒚∈ℝd-1},

We know that ϕ~ is Lipschitz on ℝd-1 and that Lip⁢(ϕ~)=Lip⁢(ϕ) (see for instance [2, Proposition 2.12]). Moreover, it is clear that ϕ~ is convex on ℝd-1. The fact that BrS⊂S and it is centered on the xd-axis (without loss of generality, we can suppose that its center is 𝝃=(,ξd) with ξd>0) translates into the fact that BrS is contained in the epigraph of ϕ~ and its center has distance from ∈ℝd at most hS. Let now 𝒑∈∂⁡ϕ~⁢(0), where ∂⁡ϕ~ is the subdifferential of ϕ~. Then, for every 𝒚∈ℝd-1 we have

ϕ~⁢(𝒚)≥𝒑⋅𝒚.

By choosing 𝒚=λ⁢𝒑, with λ≠0, we get the inequality

(A.4)ϕ~⁢(λ⁢𝒑)λ⁢|𝒑|≥|𝒑|.

Since the epigraph of ϕ~ contains BrS, which is centered at a height less than hS on the xd-axis, and by the convexity of ϕ~, we have that the truncated cone

C={(𝒙′,xn)∈ℝd-1×ℝ∣hS≥xn≥hSrS⁢|𝒙′|}

is contained in the epigraph of ϕ~ (see Figure 4). Then from (A.4) we get

|𝒑|2≤ϕ~⁢(𝒑)≤hSrS⁢|𝒑|

and so, by Assumption 2.1,

|𝒑|≤hSrS≤1ϱ.

Let now 𝒚∈ℝd-1. Then

|ϕ~⁢(𝒚)-ϕ~⁢(𝒙)|≤|𝒑⋅(𝒚-𝒙)|≤|𝒑|⁢|𝒚-𝒙|≤ϱ-1⁢|𝒚-𝒙|.

Since 𝒙 is arbitrary, this shows that Lip⁢(ϕ)=Lip⁢(ϕ~)=ϱ-1, and so that ∂⁡S is ϱ-1-Lipschitz.

Figure 4

Illustration for point (ii) in the proof of Lemma 2.3.

References

[1] J. Aghili, S. Boyaval and D. A. Di Pietro, Hybridization of mixed high-order methods on general meshes and application to the Stokes equations, Comput. Methods Appl. Math. 15 (2015), no. 2, 111–134. 10.1515/cmam-2015-0004Search in Google Scholar

[2] L. Ambrosio, N. Fusco and D. Pallara, Functions of Bounded Variation and Free Discontinuity Problems, Oxford Math. Monogr., Oxford University Press, Oxford, 2000. Search in Google Scholar

[3] P. F. Antonietti, S. Giani and P. Houston, hp-version composite discontinuous Galerkin methods for elliptic problems on complicated domains, SIAM J. Sci. Comput. 35 (2013), no. 3, A1417–A1439. 10.1137/120877246Search in Google Scholar

[4] B. Ayuso de Dios, K. Lipnikov and G. Manzini, The nonconforming virtual element method, ESAIM Math. Model Numer. Anal. (M2AN) 50 (2016), no. 3, 879–904. 10.1051/m2an/2015090Search in Google Scholar

[5] I. Babuška and M. Suri, The h-p version of the finite element method with quasi-uniform meshes, RAIRO Modél. Math. Anal. Numér. 21 (1987), no. 2, 199–238. 10.1051/m2an/1987210201991Search in Google Scholar

[6] I. Babuška and M. Suri, The optimal convergence rate of the p-version of the finite element method, SIAM J. Numer. Anal. 24 (1987), no. 4, 750–776. 10.1137/0724049Search in Google Scholar

[7] I. Babuška, B. A. Szabo and I. N. Katz, The p-version of the finite element method, SIAM J. Numer. Anal. 18 (1981), no. 3, 515–545. 10.1137/0718033Search in Google Scholar

[8] L. Beirão da Veiga, A. Chernov, L. Mascotto and A. Russo, Basic principes of hp virtual elements on quasiuniform meshes, Math. Models Methods Appl. Sci. 26 (2016), no. 8, 1567–1598. 10.1142/S021820251650038XSearch in Google Scholar

[9] A. Cangiani, E. H. Georgoulis and P. Houston, hp-version discontinuous Galerkin methods on polygonal and polyhedral meshes, Math. Models Methods Appl. Sci. 24 (2014), no. 10, 2009–2041. 10.1007/978-3-319-67673-9Search in Google Scholar

[10] P. Castillo, B. Cockburn, D. Scötzau and C. Schwab, Optimal a priori error estimates for the hp-version of the local discontinuous Galerkin method for convection-diffusion problems, Math. Comp. 71 (2001), no. 238, 455–478. 10.1090/S0025-5718-01-01317-5Search in Google Scholar

[11] P. G. Ciarlet, The Finite Element Method for Elliptic Problems, Classics Appl. Math. 40, Society for Industrial and Applied Mathematics, Philadelphia, 2002. 10.1137/1.9780898719208Search in Google Scholar

[12] B. Cockburn, D. A. Di Pietro and A. Ern, Bridging the hybrid high-order and hybridizable discontinuous Galerkin methods, ESAIM Math. Model. Numer. Anal. (M2AN) 50 (2016), no. 3, 635–650. 10.1051/m2an/2015051Search in Google Scholar

[13] B. Cockburn, J. Gopalakrishnan and R. Lazarov, Unified hybridization of discontinuous Galerkin, mixed, and continuous Galerkin methods for second order elliptic problems, SIAM J. Numer. Anal. 47 (2009), no. 2, 1319–1365. 10.1137/070706616Search in Google Scholar

[14] D. A. Di Pietro and J. Droniou, A hybrid high-order method for Leray–Lions elliptic equations on general meshes, Math. Comp. 86 (2016), no. 307, 2159–2191. 10.1090/mcom/3180Search in Google Scholar

[15] D. A. Di Pietro and J. Droniou, Ws,p-approximation properties of elliptic projectors on polynomial spaces, with application to the error analysis of a hybrid high-order discretisation of Leray–Lions problems, Math. Models Methods Appl. Sci. 27 (2017), no. 5, 879–908. 10.1142/S0218202517500191Search in Google Scholar

[16] D. A. Di Pietro and A. Ern, Mathematical Aspects of Discontinuous Galerkin Methods, Math. Appl. (Berlin) 69, Springer, Berlin, 2012. 10.1007/978-3-642-22980-0Search in Google Scholar

[17] D. A. Di Pietro and A. Ern, A hybrid high-order locking-free method for linear elasticity on general meshes, Comput. Methods Appl. Mech. Engrg. 283 (2015), 1–21. 10.1016/j.cma.2014.09.009Search in Google Scholar

[18] D. A. Di Pietro, A. Ern and S. Lemaire, An arbitrary-order and compact-stencil discretization of diffusion on general meshes based on local reconstruction operators, Comput. Methods Appl. Math. 14 (2014), no. 4, 461–472. 10.1515/cmam-2014-0018Search in Google Scholar

[19] D. A. Di Pietro, A. Ern and S. Lemaire, A review of hybrid high-order methods: Formulations, computational aspects, comparison with other methods, Building Bridges: Connections and Challenges in Modern Approaches to Numerical Partial Differential Equations (Durham 2014), Lecture Notes in Comput. Sci. 114, Springer, Cham (2016), 205–236. 10.1007/978-3-319-41640-3_7Search in Google Scholar

[20] D. A. Di Pietro and S. Lemaire, An extension of the Crouzeix–Raviart space to general meshes with application to quasi-incompressible linear elasticity and Stokes flow, Math. Comp. 84 (2015), no. 291, 1–31. 10.1090/S0025-5718-2014-02861-5Search in Google Scholar

[21] D. A. Di Pietro and R. Tittarelli, An introduction to hybrid high-order methods, Lectures from the Fall 2016 Thematic Quarter at Institut Henri Poincaré, Springer, to appear. 10.1007/978-3-319-94676-4_4Search in Google Scholar

[22] T. Dupont and R. Scott, Polynomial approximation of functions in Sobolev spaces, Math. Comp. 34 (1980), no. 150, 441–463. 10.1090/S0025-5718-1980-0559195-7Search in Google Scholar

[23] A. Ern and J.-L. Guermond, Theory and Practice Of Finite Elements, Appl. Math. Sci. 159, Springer, New York, 2004. 10.1007/978-1-4757-4355-5Search in Google Scholar

[24] E. H. Georgoulis and E. Süli, Optimal error estimates for the hp-version interior penalty discontinuous Galerkin finite element method, IMA J. Numer. Anal. 25 (2005), 205–220. 10.1093/imanum/drh014Search in Google Scholar

[25] S. Giani and P. Houston, hp-adaptive composite discontinuous Galerkin methods for elliptic problems on complicated domains, Numer. Methods Partial Differential Equations 30 (2014), no. 4, 1342–1367. 10.1002/num.21872Search in Google Scholar

[26] R. Herbin and F. Hubert, Benchmark on discretization schemes for anisotropic diffusion problems on general grids, Finite Volumes for Complex Applications V, John Wiley & Sons, New York (2008), 659–692. Search in Google Scholar

[27] C. Le Potier, A finite volume method for the approximation of highly anisotropic diffusion operators on unstructured meshes, Finite Volumes for Complex Applications IV, ISTE, London (2005), 401–412. Search in Google Scholar

[28] I. Perugia and D. Schötzau, A hp-analysis of the local discontinuous Galerkin method for diffusion problems, J. Sci. Comput. 17 (2002), no. 1–4, 561–571. 10.1023/A:1015118613130Search in Google Scholar

[29] B. Rivière, M. F. Wheeler and V. Girault, Improved energy estimates for interior penalty, constrained and discontinuous Galerkin methods for elliptic problems. Part I, Comput. Geosci. 3 (1999), 337–360. 10.1023/A:1011591328604Search in Google Scholar

[30] C. Schwab, p- and hp-FEM. Theory and Application to Solid and Fluid Mechanics, Oxford University Press, Oxford, 1998. Search in Google Scholar

[31] B. Stamm and T. P. Wihler, hp-optimal discontinuous Galerkin methods for linear elliptic problems, Math. Comp. 79 (2010), no. 272, 2117–2133. 10.1090/S0025-5718-10-02335-5Search in Google Scholar

[32] E. M. Stein, Singular Integrals and Differentiability Properties of Functions, Princeton Math. Ser. 30, Princeton University Press, Princeton, 1970. 10.1515/9781400883882Search in Google Scholar

[33] C. Talischi, G. H. Paulino, A. Pereira and I. F. M. Menezes, PolyMesher: A general-purpose mesh generator for polygonal elements written in Matlab, Struct. Multidiscip. Optim. 45 (2012), no. 3, 309–328. 10.1007/s00158-011-0706-zSearch in Google Scholar

[34] M. Vohralík, A posteriori error estimates for lowest-order mixed finite element discretizations of convection-diffusion-reaction equations, SIAM J. Numer. Anal. 45 (2007), no. 4, 1570–1599. 10.1137/060653184Search in Google Scholar

Received: 2017-4-4

Revised: 2017-4-6

Accepted: 2017-4-7

Published Online: 2017-6-17

Published in Print: 2017-7-1

An hp-Hybrid High-Order Method for Variable Diffusion on General Meshes

Abstract

1 Introduction

2 Setting

2.1 Mesh and Notation

Assumption 2.1 (Admissible Mesh Sequence).

Assumption 2.2 (Compliant Mesh Sequence).

2.2 Basic Results

Lemma 2.3 (Approximation).

Proof.

Lemma 2.4 (Discrete Trace Inequality).

Proof.

3 Discretization

3.1 The hp-HHO Method

Remark 3.1 (Static Condensation).

Remark 3.2 (Finite Element Interpretation).

3.2 Main Results

Theorem 3.3 (Energy Error Estimate).

Proof.

Theorem 3.4 (L2-Error Estimate).

3.3 Numerical Examples

4 Convergence Analysis

4.1 Consistency of the Potential Reconstruction

Proposition 4.1 (Characterization of rTpT+1∘I¯Tp¯T).

Proof.

Lemma 4.2 (Approximation Properties of ϖκ,Tl).

Proof.

4.2 Consistency of the Stabilization Term

Lemma 4.3 (Consistency of the Stabilization Term).

Proof.

4.3 Energy Error Estimate

Proof of Theorem 3.3.

Remark 4.4 (Nonconformity Error).

Proposition 4.5 (Estimate of Boundary Difference Seminorm).

Proof.

4.4 L2-Error Estimate

Proof of Theorem 3.4.

A Proof of Lemma 2.3

References

Journal and Issue

Articles in the same Issue