Publicly Available Published by De Gruyter May 24, 2022

DPG Methods for a Fourth-Order div Problem

Thomas Führer , Pablo Herrera and Norbert Heuer

From the journal Computational Methods in Applied Mathematics

https://doi.org/10.1515/cmam-2021-0246

Abstract

We study a fourth-order div problem and its approximation by the discontinuous Petrov–Galerkin method with optimal test functions. We present two variants, based on first and second-order systems. In both cases, we prove well-posedness of the formulation and quasi-optimal convergence of the approximation. Our analysis includes the fully-discrete schemes with approximated test functions, for general dimension and polynomial degree in the first-order case, and for two dimensions and lowest-order approximation in the second-order case. Numerical results illustrate the performance for quasi-uniform and adaptively refined meshes.

Keywords: Fourth-Order div Equation; Fourth-Order Elliptic PDE; Discontinuous Petrov–Galerkin Method; Ultraweak Formulation; Trace Operators; Optimal Test Functions; Fortin Operator

MSC 2010: 35J35; 65N30; 35J67

1 Introduction

There has been a sustained interest in the numerical analysis of fourth-order problems. Important examples are some models for the mechanics of thin structures [13, 24, 25], and more recently the quad-curl problem related with Maxwell transmission eigenvalue problems [2, 27] and magneto-hydrodynamics with hyper-resistivity [29, 5]. In this paper, we study a fourth-order problem with the ( ∇ ⁡ div ) 2 -operator which appears in strain-gradient elasticity; see, e.g., [26]. In two dimensions, this operator is nothing but the quad-curl operator whose approximation has been studied, e.g., in [1, 23, 31], to give a few references. In [12], a mixed method for the very ( ∇ ⁡ div ) 2 -operator is considered.

The motivation of this work is to continue developing DPG techniques for fourth-order problems. The DPG framework, short for discontinuous Petrov–Galerkin method with optimal test functions, has been established by Demkowicz and Gopalakrishnan [8] to aim at automatic discrete stability of (Petrov) Galerkin schemes for well-posed variational formulations. This is particularly relevant for fluid flow and singularly perturbed problems; cf., e.g., [6, 10, 15, 22, 30]. The DPG approach has important advantages over traditional and mixed Galerkin settings. First, any well-posed variational formulation is suitable. Second, this includes ultraweak formulations where conforming approximations are easier to realize since only trace variables are affected. Third, being a minimum residual method, DPG has built-in a posteriori error estimation and adaptivity; cf. [3, 7, 9].

The development of stable formulations and characterization of traces is at the heart of our studies of DPG formulations for thin structures, [14, 17, 18, 19]. Here, we continue the DPG development for higher-order problems, considering a fourth-order div problem. We present two variational formulations of ultraweak type based on systems of first and second order, and analyze their approximations by the DPG method. As we will see, the analysis of the first-order case is relatively straightforward, though including that of a saddle point system with non-zero diagonal blocks. On the other hand, the second-order formulation reduces the number of unknowns (though static condensation of field variables is an option) but complicates the analysis of traces typically appearing in ultraweak formulations. This is particularly the case when analyzing the fully discrete setting that includes the approximation of optimal test functions; cf. [21]. In the end, the approximation of traces is identical in both settings.

Let us finish the introduction with presenting our model problem and giving a description of what follows. Let Ω ⊂ R d ( d ≥ 2 ) be a bounded, simply connected Lipschitz domain with boundary Γ = ∂ ⁡ Ω . For a given vector function 𝒇, we consider the problem of finding the vector function 𝒖 that satisfies

(1.1) ( ∇ ⁡ div ) 2 ⁢ u + u = f in ⁢ Ω , div ⁡ u = u ⋅ n = 0 on ⁢ Γ .

Here, 𝒏 indicates the exterior unit normal vector on Γ, and below will be used generically along the boundary of any Lipschitz domain. For ease of presentation, we have selected homogeneous boundary conditions. This is not essential, and in fact, one of our numerical examples requires non-homogeneous traces on a part of the boundary.

In the next section, we recall and study trace operators and spaces that stem from the ( ∇ ⁡ div ) 2 -operator. For the first-order system, these are classical operators, recalled in § 2.1, whereas in § 2.2, we introduce operators needed for the second-order system. The first-order formulation, along with the general DPG method and a specific discrete setting, is introduced and analyzed in Section 3. Well-posedness of the formulation and quasi-optimal convergence of the DPG method with optimal test functions and with approximated test functions are stated in Theorems 7, 8 and Corollary 9. Here, the fully discrete case considers a general space dimension and arbitrary polynomial degrees (though uniform for simplicity); see § 3.1. Proofs are collected in § 3.2.2. Section 4 considers the second-order system, using the same setup as in Section 3. Well-posedness of the formulation and quasi-optimal convergence of the schemes are stated by Theorems 12, 13 and Corollary 14. The main difference in this case is that the fully discrete setting is only analyzed for two space dimensions and lowest-order approximation; see § 4.1. Proofs for the setup of the second-order system are collected in § 4.2.2. In Section 5, we present numerical results for both DPG schemes, considering an example with smooth solution and an example with corner singularity. In the latter case, we also consider adaptive variants of both schemes. In all cases, we observe the expected convergence orders.

Throughout the paper, a ≲ b means that a ≤ c ⁢ b with a generic constant c > 0 that is independent of involved functions and diameters of elements, except otherwise noted. Similarly, we use the notation a ≳ b , and a ≃ b means that a ≲ b and b ≲ a .

2 Trace Spaces and Norms

We consider a mesh T = { T } consisting of general non-intersecting Lipschitz elements so that Ω ¯ = ⋃ { T ¯ ; T ∈ T } , and denote by S = { ∂ ⁡ T ; T ∈ T } its skeleton.

Throughout this section, let ω ⊂ Ω denote a Lipschitz domain, for simplicity simply connected, with boundary ∂ ⁡ ω . We use the standard L 2 spaces L 2 ⁢ ( ω ) , L 2 ⁢ ( ω ) = ( L 2 ⁢ ( ω ) ) d of scalar and vector fields, respectively, with generic norm ∥ ⋅ ∥ ω and L 2 ⁢ ( ω ) -bilinear form ( ⋅ , ⋅ ) ω . We drop the index 𝜔 when ω = Ω . We also use the standard Sobolev spaces H 1 ⁢ ( ω ) and H 0 1 ⁢ ( ω ) with norm ∥ ⋅ ∥ 1 , ω . Furthermore,

H ⁢ ( div , ω ) := { v ∈ L 2 ⁢ ( ω ) ; div ⁡ v ∈ L 2 ⁢ ( ω ) } , H 0 ⁢ ( div , ω ) := { v ∈ H ⁢ ( div , ω ) ; ( v ⋅ n ) | ∂ ⁡ ω = 0 } , H ⁢ ( ∇ ⁡ div , ω ) := { v ∈ L 2 ⁢ ( ω ) ; ∇ ⁡ div ⁡ v ∈ L 2 ⁢ ( ω ) } , H 0 ⁢ ( ∇ ⁡ div , ω ) := { v ∈ H ⁢ ( ∇ ⁡ div , ω ) ; ( v ⋅ n ) | ∂ ⁡ ω = ( div ⁡ v ) | ∂ ⁡ ω = 0 }

with respective norms

∥ v ∥ div , ω 2 := ∥ v ∥ ω 2 + ∥ div ⁡ v ∥ ω 2 , ∥ v ∥ ∇ ⁡ div , ω 2 := ∥ v ∥ ω 2 + ∥ ∇ ⁡ div ⁡ v ∥ ω 2 ,

and dropping 𝜔 when ω = Ω . We also need the corresponding product spaces, denoted in the same way but replacing 𝜔 with 𝒯, e.g., H ⁢ ( ∇ ⁡ div , T ) = Π T ∈ T ⁢ H ⁢ ( ∇ ⁡ div , T ) with canonical product norm ∥ ⋅ ∥ ∇ ⁡ div , T . The L 2 ⁢ ( T ) -dualities are denoted by ( ⋅ , ⋅ ) T .

2.1 Canonical Trace Operators

We define operators

tr grad : { H 1 ⁢ ( T ) → H ⁢ ( div , T ) ′ , v ↦ ⟨ tr grad ⁡ ( v ) , τ ⟩ S ∇ := ( v , div ⁡ τ ) T + ( ∇ ⁡ v , τ ) T ,

tr div : { H ⁢ ( div , T ) → H 1 ⁢ ( T ) ′ , τ ↦ ⟨ tr div ⁡ ( τ ) , v ⟩ S div := ⟨ tr grad ⁡ ( v ) , τ ⟩ S ∇ .

Analogously, we define the canonical trace operators for ω ⊂ Ω (instead of 𝒯) with support on ∂ ⁡ ω ,

tr ω grad : H 1 ⁢ ( ω ) → H ⁢ ( div , ω ) ′ and tr ω div : H ⁢ ( div , ω ) → H 1 ⁢ ( ω ) ′ .

Of course, given v ∈ H 1 ⁢ ( ω ) and τ ∈ H ⁢ ( div , ω ) , tr ω grad ⁡ ( v ) = v | ∂ ⁡ ω and tr ω div ⁡ ( τ ) = ( τ ⋅ n ) | ∂ ⁡ ω are the classical traces living in the canonical Sobolev spaces H 1 / 2 ⁢ ( ∂ ⁡ ω ) = tr ω grad ⁡ ( H 1 ⁢ ( ω ) ) and H - 1 / 2 ⁢ ( ∂ ⁡ ω ) = tr ω div ⁡ ( H ⁢ ( div , ω ) ) , respectively.

The statements of the following proposition are proved in [4, Theorem 2.3, Remark 2.5].

Proposition 1

The following statements hold.

If v ∈ H 1 ⁢ ( T ) , then
v ∈ H 0 1 ⁢ ( Ω ) ⇔ ⟨ tr div ⁡ ( τ ) , v ⟩ S div = 0 for all ⁢ τ ∈ H ⁢ ( div , Ω ) , v ∈ H 1 ⁢ ( Ω ) ⇔ ⟨ tr div ⁡ ( τ ) , v ⟩ S div = 0 for all ⁢ τ ∈ H 0 ⁢ ( div , Ω ) .
If τ ∈ H ⁢ ( div , T ) , then
τ ∈ H 0 ⁢ ( div , Ω ) ⇔ ⟨ tr grad ⁡ ( v ) , τ ⟩ S ∇ = 0 for all ⁢ v ∈ H 1 ⁢ ( Ω ) , τ ∈ H ⁢ ( div , Ω ) ⇔ ⟨ tr grad ⁡ ( v ) , τ ⟩ S ∇ = 0 for all ⁢ v ∈ H 0 1 ⁢ ( Ω ) .

The canonical trace spaces are

H 1 / 2 ⁢ ( ∂ ⁡ ω ) := tr ω grad ⁡ ( H 1 ⁢ ( ω ) ) , H 1 / 2 ⁢ ( S ) := tr grad ⁡ ( H 1 ⁢ ( Ω ) ) , H 0 1 / 2 ⁢ ( S ) := tr grad ⁡ ( H 0 1 ⁢ ( Ω ) ) , H - 1 / 2 ⁢ ( ∂ ⁡ ω ) := tr ω div ⁡ ( H ⁢ ( div , ω ) ) , H - 1 / 2 ⁢ ( S ) := tr div ⁡ ( H ⁢ ( div , Ω ) ) , H 0 - 1 / 2 ⁢ ( S ) := tr div ⁡ ( H 0 ⁢ ( div , Ω ) )

with respective operator norms

∥ v ^ ∥ ( div , ω ) ′ := sup ∥ τ ∥ div , ω = 1 ⁡ ⟨ v ^ , τ ⟩ ∂ ⁡ ω ∇ ( v ^ ∈ H 1 / 2 ⁢ ( ∂ ⁡ ω ) ) , ∥ v ^ ∥ ( div , T ) ′ := sup ∥ τ ∥ div , T = 1 ⁡ ⟨ v ^ , τ ⟩ S ∇ ( v ^ ∈ H 1 / 2 ⁢ ( S ) ) , ∥ τ ^ ∥ ( 1 , ω ) ′ := sup ∥ v ∥ 1 , ω = 1 ⁡ ⟨ τ ^ , v ⟩ ∂ ⁡ ω div ( τ ^ ∈ H - 1 / 2 ⁢ ( ∂ ⁡ ω ) ) , ∥ τ ^ ∥ ( 1 , T ) ′ := sup ∥ v ∥ 1 , T = 1 ⁡ ⟨ τ ^ , v ⟩ S div ( τ ^ ∈ H - 1 / 2 ⁢ ( S ) ) .

Here, the dualities between the trace spaces and their corresponding test spaces are taken to be consistent with the trace operator. For instance, given v ^ ∈ H 1 / 2 ⁢ ( S ) and τ ∈ H ⁢ ( div , T ) ,

⟨ v ^ , τ ⟩ S ∇ := ⟨ tr grad ⁡ ( v ) , τ ⟩ S ∇ with ⁢ v ∈ H 1 ⁢ ( Ω ) ⁢ such that ⁢ tr grad ⁡ ( v ) = v ^ .

Alternatively, the canonical trace norms are

∥ v ^ ∥ 1 / 2 , ∂ ⁡ ω := inf { ∥ v ∥ 1 , ω ; v ∈ H 1 ( ω ) , tr ω grad ( v ) = v ^ } ( v ^ ∈ H 1 / 2 ⁢ ( ∂ ⁡ ω ) ) , ∥ v ^ ∥ 1 / 2 , S := inf { ∥ v ∥ 1 ; v ∈ H 1 ( Ω ) , tr grad ( v ) = v ^ } ( v ^ ∈ H 1 / 2 ⁢ ( S ) ) , ∥ τ ^ ∥ - 1 / 2 , ∂ ⁡ ω := inf { ∥ τ ∥ div , ω ; τ ∈ H ( div , ω ) , tr ω div ( τ ) = τ ^ } ( τ ^ ∈ H - 1 / 2 ⁢ ( ∂ ⁡ ω ) ) , ∥ τ ^ ∥ - 1 / 2 , S := inf { ∥ τ ∥ div ; τ ∈ H ( div , Ω ) , tr div ( τ ) = τ ^ } ( τ ^ ∈ H - 1 / 2 ⁢ ( S ) ) .

The equivalence of the corresponding trace norms will be relevant. In fact, they are pairwise identical as stated by the next proposition.

Proposition 2

Proposition 2 ([4, Lemma 2.1])

The identities

∥ v ^ ∥ 1 / 2 , ∂ ⁡ ω = ∥ v ^ ∥ ( div , ω ) ′ ( v ^ ∈ H 1 / 2 ⁢ ( ∂ ⁡ ω ) ) , ∥ v ^ ∥ 1 / 2 , S = ∥ v ^ ∥ ( div , T ) ′ ( v ^ ∈ H 1 / 2 ⁢ ( S ) ) , ∥ τ ^ ∥ - 1 / 2 , ∂ ⁡ ω = ∥ τ ^ ∥ ( 1 , ω ) ′ ( τ ^ ∈ H - 1 / 2 ⁢ ( ∂ ⁡ ω ) ) , ∥ τ ^ ∥ - 1 / 2 , S = ∥ τ ^ ∥ ( 1 , T ) ′ ( τ ^ ∈ H - 1 / 2 ⁢ ( S ) )

hold true.

2.2 A grad-div Trace Operator

Proposition 3

Let ω ⊂ Ω be a Lipschitz domain. The norms ∥ ⋅ ∥ ∇ ⁡ div , ω and ∥ div ⋅ ∥ ω + ∥ ⋅ ∥ ∇ ⁡ div , ω are equivalent in H ⁢ ( ∇ ⁡ div , ω ) , with the non-trivial constant depending on 𝜔.

Proof

Noting that ∇ ⁡ div ⁡ v ∈ L 2 ⁢ ( ω ) implies that div ⁡ v ∈ L 2 ⁢ ( Ω ) , it follows that H ⁢ ( ∇ ⁡ div , ω ) is Banach with respect to both norms, ∥ ⋅ ∥ ∇ ⁡ div , ω and ∥ ⋅ ∥ ∇ ⁡ div , ω + ∥ div ⋅ ∥ ω . The trivial bound ∥ ⋅ ∥ ∇ ⁡ div , ω ≤ ∥ ⋅ ∥ ∇ ⁡ div , ω + ∥ div ⋅ ∥ ω then implies the equivalence of both norms. ∎

Analogously to the canonical traces, we define traces induced by the operator ∇ ⁡ div ,

tr ∇ ⁡ div : { H ⁢ ( ∇ ⁡ div , T ) → H ⁢ ( ∇ ⁡ div , T ) ′ , v ↦ ⟨ tr ∇ ⁡ div ⁡ ( v ) , τ ⟩ S ∇ ⁡ div := ( v , ∇ ⁡ div ⁡ τ ) T - ( ∇ ⁡ div ⁡ v , τ ) T ,

and analogously for ω ⊂ Ω (instead of 𝒯), tr ω ∇ ⁡ div : H ⁢ ( ∇ ⁡ div , ω ) → H ⁢ ( ∇ ⁡ div , ω ) ′ with support on ∂ ⁡ ω , and duality ⟨ ⋅ , ⋅ ⟩ ∂ ⁡ ω ∇ ⁡ div . By Proposition 3, v ∈ H ⁢ ( ∇ ⁡ div , ω ) implies that v ∈ H ⁢ ( div , ω ) and div ⁡ v ∈ H 1 ⁢ ( ω ) . Therefore, tr ω ∇ ⁡ div ⁡ ( v ) has the two well-posed classical trace components v ⋅ n | ∂ ⁡ ω and ( div ⁡ v ) | ∂ ⁡ ω , and using the trace operators tr grad , tr div introduced in § 2.1, we can write

(2.1) ⟨ tr ω ∇ ⁡ div ⁡ ( v ) , τ ⟩ ∂ ⁡ ω ∇ ⁡ div = ⟨ tr ω div ⁡ ( v ) , div ⁡ τ ⟩ ∂ ⁡ ω div - ⟨ tr ω grad ⁡ ( div ⁡ v ) , τ ⟩ ∂ ⁡ ω ∇ for all ⁢ v , τ ∈ H ⁢ ( ∇ ⁡ div , ω ) .

Though, the separate canonical components are not uniformly bounded with respect to the domain 𝜔 under consideration. In the following, tr ω ∇ ⁡ div ⁡ ( v ) = ( g 1 , g 2 ) for v ∈ H ⁢ ( ∇ ⁡ div , ω ) means that g 1 = tr ω div ⁡ ( v ) = ( v ⋅ n ) | ∂ ⁡ ω and g 2 = tr ω grad ⁡ ( div ⁡ v ) = ( div ⁡ v ) | ∂ ⁡ ω .

We proceed with considering the trace operator tr ∇ ⁡ div and its induced trace spaces. We define local trace spaces

H ∇ ⁡ div ⁢ ( ∂ ⁡ ω ) := tr ω ∇ ⁡ div ⁡ ( H ⁢ ( ∇ ⁡ div , ω ) )

and the “skeleton” variants

H ∇ ⁡ div ⁢ ( S ) := tr ∇ ⁡ div ⁡ ( H ⁢ ( ∇ ⁡ div , Ω ) ) , H 0 ∇ ⁡ div ⁢ ( S ) := tr ∇ ⁡ div ⁡ ( H 0 ⁢ ( ∇ ⁡ div , Ω ) )

provided with operator norms

∥ v ^ ∥ ( ∇ ⁡ div , ω ) ′ := sup ∥ τ ∥ ∇ ⁡ div , ω = 1 ⁡ ⟨ v ^ , τ ⟩ ∂ ⁡ ω ∇ ⁡ div ( v ^ ∈ H ∇ ⁡ div ⁢ ( ∂ ⁡ ω ) ) , ∥ v ^ ∥ ( ∇ ⁡ div , T ) ′ := sup ∥ τ ∥ ∇ ⁡ div , T = 1 ⁡ ⟨ v ^ , τ ⟩ S ∇ ⁡ div ( v ^ ∈ H ∇ ⁡ div ⁢ ( S ) ) .

As before, the dualities are defined to be consistent with the definition of the trace operator. The canonical trace norms are

∥ v ^ ∥ - 1 / 2 , 1 / 2 , ∂ ⁡ ω := inf { ∥ v ∥ ∇ ⁡ div , ω ; v ∈ H ( ∇ div , ω ) , tr ω ∇ ⁡ div ( v ) = v ^ } ( v ^ ∈ H ∇ ⁡ div ⁢ ( ∂ ⁡ ω ) ) , ∥ v ^ ∥ - 1 / 2 , 1 / 2 , S := inf { ∥ v ∥ ∇ ⁡ div ; v ∈ H ( ∇ div , Ω ) , tr ∇ ⁡ div ( v ) = v ^ } ( v ^ ∈ H ∇ ⁡ div ⁢ ( S ) ) ,

and in Proposition 5, we will see that ∥ ⋅ ∥ - 1 / 2 , 1 / 2 , S = ∥ ⋅ ∥ ( ∇ ⁡ div , T ) ′ .

The following proposition relates the combined ∇ ⁡ div traces with the canonical ones. We formulate it for an arbitrary Lipschitz domain but only will use it for Ω.

Proposition 4

Let ω ⊂ Ω be a Lipschitz domain. The spaces

( H ∇ ⁡ div ⁢ ( ∂ ⁡ ω ) , ∥ ⋅ ∥ - 1 / 2 , 1 / 2 , ∂ ⁡ ω ) , ( H - 1 / 2 ⁢ ( ∂ ⁡ ω ) , ∥ ⋅ ∥ - 1 / 2 , ∂ ⁡ ω ) × ( H 1 / 2 ⁢ ( ∂ ⁡ ω ) , ∥ ⋅ ∥ 1 / 2 , ∂ ⁡ ω )

are identical with equivalent norms and equivalence constants that may depend on 𝜔. More specifically, the inclusion H - 1 / 2 ⁢ ( ∂ ⁡ ω ) × H 1 / 2 ⁢ ( ∂ ⁡ ω ) ⊂ H ∇ ⁡ div ⁢ ( ∂ ⁡ ω ) is due to the existence of a bounded extension operator E ω : H - 1 / 2 ⁢ ( ∂ ⁡ ω ) × H 1 / 2 ⁢ ( ∂ ⁡ ω ) → H ⁢ ( ∇ ⁡ div , ω ) with tr ω ∇ ⁡ div ∘ E ω = id .

Proof

We already have seen trace identification (2.1) which shows the algebraic inclusion

H ∇ ⁡ div ⁢ ( ∂ ⁡ ω ) ⊂ H - 1 / 2 ⁢ ( ∂ ⁡ ω ) × H 1 / 2 ⁢ ( ∂ ⁡ ω ) .

Now, to show the other inclusion and its continuity, let ( g 1 , g 2 ) ∈ H - 1 / 2 ⁢ ( ∂ ⁡ ω ) × H 1 / 2 ⁢ ( ∂ ⁡ ω ) be given. We extend g 1 := ( g 1 , 0 ) ↦ v 1 ∈ H ⁢ ( ∇ ⁡ div , ω ) and g 2 := ( 0 , g 2 ) ↦ v 2 ∈ H ⁢ ( ∇ ⁡ div , ω ) such that

tr ω ∇ ⁡ div ⁡ ( v j ) = ( tr ω div ⁡ ( v j ) , tr ω grad ⁡ ( div ⁡ v j ) ) = g j ,

and with boundedness ∥ v 1 ∥ ∇ ⁡ div , ω ≲ ∥ g 1 ∥ - 1 / 2 , ∂ ⁡ ω and ∥ v 2 ∥ ∇ ⁡ div , ω ≲ ∥ g 2 ∥ 1 / 2 , ∂ ⁡ ω ( j = 1 , 2 ). The extension

( g 1 , g 2 ) ↦ E ω ⁢ ( g 1 , g 2 ) := v := v 1 + v 2

then satisfies v ∈ H ⁢ ( ∇ ⁡ div , ω ) with tr ω ∇ ⁡ div ⁡ ( v ) = ( g 1 , g 2 ) and the bound

∥ tr ω ∇ ⁡ div ⁡ ( v ) ∥ - 1 / 2 , 1 / 2 , ∂ ⁡ ω ≤ ∥ v ∥ ∇ ⁡ div , ω ≲ ∥ g 1 ∥ - 1 / 2 , ∂ ⁡ ω + ∥ g 2 ∥ 1 / 2 , ∂ ⁡ ω .

The equivalence of both norms then follows since H ∇ ⁡ div ⁢ ( ∂ ⁡ ω ) = H - 1 / 2 ⁢ ( ∂ ⁡ ω ) × H 1 / 2 ⁢ ( ∂ ⁡ ω ) are Banach spaces.

Extension of g 1 . We define ϕ ∈ H 0 1 ⁢ ( ω ) and ψ ∈ H 1 ⁢ ( ω ) as the solutions to

- Δ ⁢ ϕ = α in ⁢ ω , - Δ ⁢ ψ = ϕ in ⁢ ω , ∂ n ⁡ ψ = - g 1 on ⁢ ∂ ⁡ ω ,

where α ∈ R is chosen so that ( ϕ , 1 ) ω - ⟨ g 1 , 1 ⟩ ∂ ⁡ ω div = 0 . Then we select v 1 := - ∇ ⁡ ψ and observe that

curl ⁡ v 1 = ( ∇ ⁡ div ) 2 ⁢ v 1 = 0 and tr ω ∇ ⁡ div ⁡ ( v 1 ) = ( v 1 ⋅ n , div ⁡ v 1 ) | ∂ ⁡ ω = ( g 1 , 0 ) .

Furthermore,

∥ v 1 ∥ ω 2 = ∥ ∇ ⁡ ψ ∥ ω 2 ≲ ∥ ϕ ∥ ω 2 + ∥ g 1 ∥ - 1 / 2 , ∂ ⁡ ω 2 , ∥ ∇ ⁡ div ⁡ v 1 ∥ ω 2 + ∥ ϕ ∥ ω 2 = ∥ ∇ ⁡ ϕ ∥ ω 2 + ∥ ϕ ∥ ω 2 ≲ ∥ ∇ ⁡ ϕ ∥ ω 2 .

It is not difficult to see that ∥ ∇ ⁡ ϕ ∥ ω ≃ | α | ≲ ∥ g 1 ∥ - 1 / 2 , ∂ ⁡ ω . This finishes the extension of g 1 .

Extension of g 2 . We define ϕ ∈ H 1 ⁢ ( ω ) and ψ ∈ H 1 ⁢ ( ω ) as the solutions to

- Δ ⁢ ϕ = 0 in ⁢ ω , ϕ = g 2 on ⁢ ∂ ⁡ ω , - Δ ⁢ ψ = ϕ in ⁢ ω , ∂ n ⁡ ψ = β on ⁢ ∂ ⁡ ω ,

where β ∈ R is chosen so that ( ϕ , 1 ) ω + ⟨ β , 1 ⟩ ∂ ⁡ ω div = 0 . Let us denote the extension g 1 = ( g 1 , 0 ) ↦ v 1 from the previous step as E ω 1 ⁢ ( g 1 ) . Then we define v 2 := E ω 1 ⁢ ( β ) - ∇ ⁡ ψ and observe that curl ⁡ v 2 = ( ∇ ⁡ div ) 2 ⁢ v 2 = 0 and tr ω ∇ ⁡ div ⁡ ( v 2 ) = ( 0 , g 2 ) . Furthermore,

∥ v 2 ∥ ω 2 ≤ 2 ⁢ ( ∥ E ω 1 ⁢ ( β ) ∥ ω 2 + ∥ ∇ ⁡ ψ ∥ ω 2 ) ≲ | β | 2 + ∥ ϕ ∥ ω 2 ≲ | β | 2 + ∥ g 2 ∥ 1 / 2 , ∂ ⁡ ω 2 , ∥ ∇ ⁡ div ⁡ v 2 ∥ ω 2 ≤ 2 ⁢ ( ∥ ∇ ⁡ div ⁡ E ω 1 ⁢ ( β ) ∥ ω 2 + ∥ ∇ ⁡ ϕ ∥ ω 2 ) ≲ | β | 2 + ∥ g 2 ∥ 1 / 2 , ∂ ⁡ ω 2 .

The bound | β | ≲ ∥ ϕ ∥ ω ≲ ∥ g 2 ∥ 1 / 2 , ∂ ⁡ ω finally proves the required boundedness of the extension v 2 . ∎

The following propositions will be needed to analyze our second-order ultraweak formulation of problem (1.1).

Proposition 5

The identity

∥ v ^ ∥ - 1 / 2 , 1 / 2 , S = ∥ v ^ ∥ ( ∇ ⁡ div , T ) ′ ( v ^ ∈ H ∇ ⁡ div ⁢ ( S ) )

holds true.

Proof

The proof of this result is by now standard; cf. [4, proofs of Lemma 2.2, Theorem 2.3] and [17, proofs of Lemma 3.2, Proposition 3.5]. We only recall the essential steps, following the setting from [17]. The bound ∥ ⋅ ∥ ( ∇ ⁡ div , T ) ′ ≤ ∥ ⋅ ∥ - 1 / 2 , 1 / 2 , S is immediate by definition of the trace operator. For the other inequality, let v ^ ∈ H ∇ ⁡ div ⁢ ( S ) be given. We define τ ∈ H ⁢ ( ∇ ⁡ div , T ) as the solution to

( τ , δ ⁢ τ ) + ( ∇ ⁡ div ⁡ τ , ∇ ⁡ div ⁡ δ ⁢ τ ) T = ⟨ v ^ , δ ⁢ τ ⟩ S ∇ ⁡ div for all ⁢ δ ⁢ τ ∈ H ⁢ ( ∇ ⁡ div , T )

and τ ↦ v ∈ H ⁢ ( ∇ ⁡ div , Ω ) as the solution to

( v , δ ⁢ v ) + ( ∇ ⁡ div ⁡ v , ∇ ⁡ div ⁡ δ ⁢ v ) = ⟨ tr ∇ ⁡ div ⁡ ( δ ⁢ v ) , τ ⟩ S ∇ ⁡ div for all ⁢ δ ⁢ v ∈ H ⁢ ( ∇ ⁡ div , Ω ) .

It follows that v = ∇ ⁡ div ⁡ τ (piecewise on 𝒯), leading to tr ∇ ⁡ div ⁡ ( v ) = v ^ and ∥ τ ∥ ∇ ⁡ div , T 2 = ∥ v ∥ ∇ ⁡ div 2 = ⟨ v ^ , τ ⟩ S ∇ ⁡ div , so that

∥ v ^ ∥ ( ∇ ⁡ div , T ) ′ ≥ ⟨ v ^ , τ ⟩ S ∇ ⁡ div ∥ τ ∥ ∇ ⁡ div , T = ∥ v ∥ ∇ ⁡ div ≥ ∥ v ^ ∥ - 1 / 2 , 1 / 2 , S . ∎

Proposition 6

If v ∈ H ⁢ ( ∇ ⁡ div , T ) , then

v ∈ H 0 ⁢ ( ∇ ⁡ div , Ω ) ⇔ ⟨ tr ∇ ⁡ div ⁡ ( τ ) , v ⟩ S ∇ ⁡ div = 0 for all ⁢ τ ∈ H ⁢ ( ∇ ⁡ div , Ω ) , v ∈ H ⁢ ( ∇ ⁡ div , Ω ) ⇔ ⟨ tr ∇ ⁡ div ⁡ ( τ ) , v ⟩ S ∇ ⁡ div = 0 for all ⁢ τ ∈ H 0 ⁢ ( ∇ ⁡ div , Ω ) .

Proof

The statements can be proved by standard techniques. Let us recall them briefly. In both cases, the direction “⇒” holds by integration by parts, taking into account the inherent regularities of H ⁢ ( ∇ ⁡ div , Ω ) mentioned previously. The direction “⇐” without boundary condition holds by the distributional definition of the derivatives,

∇ ⁡ div ⁡ v ⁢ ( ϕ ) := ( v , ∇ ⁡ div ⁡ ϕ ) = ( ∇ ⁡ div ⁡ v , ϕ ) T + ⟨ tr ∇ ⁡ div ⁡ ( v ) , ϕ ⟩ S ∇ ⁡ div = ( ∇ ⁡ div ⁡ v , ϕ ) T for all ⁢ ϕ ∈ C 0 ∞ ⁢ ( Ω ) d ,

giving ∇ ⁡ div ⁡ v ∈ L 2 ⁢ ( Ω ) . Now, to show the boundary condition in the first statement, let v ∈ H ⁢ ( ∇ ⁡ div , Ω ) be given. We use Proposition 4 to extend an arbitrary datum

g = ( g 1 , g 2 ) ∈ H - 1 / 2 ⁢ ( Γ ) × H 1 / 2 ⁢ ( Γ )

to E Ω ⁢ ( g ) ∈ H ⁢ ( ∇ ⁡ div , Ω ) . By assumption, the fact that E Ω is a right inverse of tr Ω ∇ ⁡ div and relation (2.1), we find that

0 = ⟨ tr ∇ ⁡ div ⁡ ( E Ω ⁢ ( g ) ) , v ⟩ S ∇ ⁡ div = ⟨ tr Ω ∇ ⁡ div ⁡ ( E Ω ⁢ ( g ) ) , v ⟩ Γ ∇ ⁡ div = ⟨ g 1 , div ⁡ v ⟩ Γ div - ⟨ g 2 , v ⟩ Γ ∇ = ⟨ g 1 , div ⁡ v ⟩ Γ - ⟨ g 2 , v ⋅ n ⟩ Γ ,

where ⟨ ⋅ , ⋅ ⟩ Γ denotes the duality between H - 1 / 2 ⁢ ( Γ ) and H 1 / 2 ⁢ ( Γ ) (in any order) with L 2 ⁢ ( Γ ) as pivot space. We conclude that div ⁡ v = v ⋅ n = 0 on Γ. ∎

3 First-Order Variational Formulation and DPG Method

Denoting u 1 := u , we write (1.1) as the first-order system (3.1)

(3.1a) ∇ ⁡ u 4 + u 1 = f ,

(3.1b) div ⁡ u 3 - u 4 = 0 ,

(3.1c) ∇ ⁡ u 2 - u 3 = 0 ,

(3.1d) div ⁡ u 1 - u 2 = 0

in Ω subject to u 1 ⋅ n = u 2 = 0 on Γ. Testing (3.1a)–(3.1d) with piecewise smooth functions v 1 , v 2 , v 3 , v 4 , respectively, and employing the trace operators tr grad , tr div , we obtain

( u 1 , v 1 - ∇ ⁡ v 4 ) T - ( u 2 , v 4 + div ⁡ v 3 ) T - ( u 3 , v 3 + ∇ ⁡ v 2 ) T - ( u 4 , v 2 + div ⁡ v 1 ) T + ⟨ tr div ⁡ ( u 1 ) , v 4 ⟩ S div + ⟨ tr grad ⁡ ( u 2 ) , v 3 ⟩ S ∇ + ⟨ tr div ⁡ ( u 3 ) , v 2 ⟩ S div + ⟨ tr grad ⁡ ( u 4 ) , v 1 ⟩ S ∇ = ( f , v 1 ) .

We replace the traces by independent trace variables

u ^ 1 = tr div ⁡ ( u 1 ) , u ^ 2 = tr grad ⁡ ( u 2 ) , u ^ 3 = tr div ⁡ ( u 3 ) , u ^ 4 = tr grad ⁡ ( u 4 )

and define the spaces

U := U 0 × U ^ with ⁢ U 0 := L 2 ⁢ ( Ω ) × L 2 ⁢ ( Ω ) × L 2 ⁢ ( Ω ) × L 2 ⁢ ( Ω ) , U ^ := H 0 - 1 / 2 ⁢ ( S ) × H 0 1 / 2 ⁢ ( S ) × H - 1 / 2 ⁢ ( S ) × H 1 / 2 ⁢ ( S )

and

V ⁢ ( T ) := H ⁢ ( div , T ) × H 1 ⁢ ( T ) × H ⁢ ( div , T ) × H 1 ⁢ ( T )

with (squared) norms

∥ ( u 1 , u 2 , u 3 , u 4 , u ^ 1 , u ^ 2 , u ^ 3 , u ^ 4 ) ∥ U 2 := ∥ ( u 1 , u 2 , u 3 , u 4 ) ∥ 2 + ∥ ( u ^ 1 , u ^ 2 , u ^ 3 , u ^ 4 ) ∥ U ^ 2 := ∥ u 1 ∥ 2 + ∥ u 2 ∥ 2 + ∥ u 3 ∥ 2 + ∥ u 4 ∥ 2 + ∥ u ^ 1 ∥ - 1 / 2 , S 2 + ∥ u ^ 2 ∥ 1 / 2 , S 2 + ∥ u ^ 3 ∥ - 1 / 2 , S 2 + ∥ u ^ 4 ∥ 1 / 2 , S 2

and

∥ ( v 1 , v 2 , v 3 , v 4 ) ∥ V ⁢ ( T ) 2 := ∥ v 1 ∥ div , T 2 + ∥ v 2 ∥ 1 , T 2 + ∥ v 3 ∥ div , T 2 + ∥ v 4 ∥ 1 , T 2 ,

respectively. Then our first variational formulation of (1.1) reads: find u = ( u 1 , u 2 , u 3 , u 4 , u ^ 1 , u ^ 2 , u ^ 3 , u ^ 4 ) ∈ U such that

(3.2) b ⁢ ( u , v ) = L ⁢ ( v ) for all ⁢ v = ( v 1 , v 2 , v 3 , v 4 ) ∈ V ⁢ ( T ) ,

in operator form u ∈ U : B ⁢ u = L in V ⁢ ( T ) ′ . Here,

b ⁢ ( u , v ) := ( u 1 , v 1 - ∇ ⁡ v 4 ) T - ( u 2 , v 4 + div ⁡ v 3 ) T - ( u 3 , v 3 + ∇ ⁡ v 2 ) T - ( u 4 , v 2 + div ⁡ v 1 ) T + ⟨ u ^ 1 , v 4 ⟩ S div + ⟨ u ^ 2 , v 3 ⟩ S ∇ + ⟨ u ^ 3 , v 2 ⟩ S div + ⟨ u ^ 4 , v 1 ⟩ S ∇

and L ⁢ ( v ) := ( f , v 1 ) .

Theorem 7

For given f ∈ L 2 ⁢ ( Ω ) , there exists a unique solution u = ( u 1 , u 2 , u 3 , u 4 , u ^ 1 , u ^ 2 , u ^ 3 , u ^ 4 ) ∈ U to (3.2). It is uniformly bounded, ∥ u ∥ U ≲ ∥ f ∥ with a hidden constant that is independent of 𝒇 and 𝒯. Furthermore, u := u 1 solves (1.1).

For a proof of this theorem, we refer to § 3.2.2.

In order to define our approximation scheme, we consider a (family of) discrete subspace(s) U h ⊂ U (being set on a sequence of meshes 𝒯) and introduce the trial-to-test operator T : U → V ⁢ ( T ) by

⟨ ⟨ T ⁢ ( u ) , v ⟩ ⟩ V ⁢ ( T ) = b ⁢ ( u , v ) for all ⁢ v ∈ V ⁢ ( T ) .

Here, ⟨ ⟨ ⋅ , ⋅ ⟩ ⟩ V ⁢ ( T ) denotes the inner product of V ⁢ ( T ) that induces the norm defined previously. Then our first DPG scheme for problem (1.1) is: find u h ∈ U h such that

(3.3) b ⁢ ( u h , T ⁢ δ ⁢ u ) = L ⁢ ( T ⁢ δ ⁢ u ) for all ⁢ δ ⁢ u ∈ U h .

The following theorem states its quasi-optimal approximation property.

Theorem 8

Let f ∈ L 2 ⁢ ( Ω ) be given. For any finite-dimensional subspace U h ⊂ U , there exists a unique solution u h ∈ U h to (3.3). It satisfies the quasi-optimal error estimate ∥ u - u h ∥ U ≲ ∥ u - w ∥ U for all w ∈ U h with a hidden constant that is independent of 𝒯.

A proof of this theorem is given in § 3.2.2.

3.1 Fully Discrete Method

In practice, the action of the trial-to-test operator 𝔗 has to be approximated. To this end, we replace the test space by a sufficiently rich finite-dimensional subspace V h ⁢ ( T ) ⊂ V ⁢ ( T ) and use the approximated operator T h : U h → V h ⁢ ( T ) defined by ⟨ ⟨ T h ⁢ ( u ) , v ⟩ ⟩ V ⁢ ( T ) = b ⁢ ( u , v ) for all v ∈ V h ⁢ ( T ) . The fully discrete DPG method then is

(3.4) u h ∈ U h : b ⁢ ( u h , T h ⁢ δ ⁢ u ) = L ⁢ ( T h ⁢ δ ⁢ u ) for all ⁢ δ ⁢ u ∈ U h .

Discrete stability and quasi-optimal convergence of this scheme is guaranteed by the existence of a (sequence of) Fortin operator(s) Π F : V ⁢ ( T ) → V h ⁢ ( T ) , i.e., Π F is uniformly bounded with respect to the underlying (sequence of) space(s) V h ⁢ ( T ) and satisfies b ⁢ ( u h , v - Π F ⁢ v ) = 0 for all u h ∈ U h , v ∈ V ⁢ ( T ) ; see [21] for details.

In order to be specific, let us fix the discretization. We consider a polygonal domain Ω ⊂ R d ( d ≥ 2 ) and meshes 𝒯 of shape-regular simplicial elements,

sup T ∈ T ⁡ diam ⁢ ( T ) d | T | ⁢ is uniformly bounded .

In the following, F T denotes the set of faces of an element T ∈ T . For a non-negative integer 𝑝, let us define the discrete spaces

P p ⁢ ( T ) := { v : T → R ; v ⁢ is a polynomial of degree ⁢ p } ( T ∈ T ) , P p ⁢ ( T ) := { v ∈ L 2 ⁢ ( Ω ) ; v | T ∈ P p ⁢ ( T ) ⁢ for all ⁢ T ∈ T } , P p ⁢ ( F ) := { v : F → R ; v ⁢ is a polynomial of degree ⁢ p } ( F ∈ F T , T ∈ T ) , P p ⁢ ( F T ) := { v ∈ L 2 ⁢ ( ∂ ⁡ T ) ; v | F ∈ P p ⁢ ( F ) ⁢ for all ⁢ F ∈ F T } ( T ∈ T ) , U div , p ⁢ ( T ) := { v ∈ H ( div , Ω ) ; ( v ⋅ n ) | ∂ ⁡ T ∈ P p ( F T ) for all T ∈ T } , P p ⁢ ( S ) := tr div ⁡ ( U div , p ⁢ ( T ) ) , P 0 p ( S ) := tr div ( U div , p ( T ) ∩ H 0 ( div , Ω ) ) , U grad , p ⁢ ( T ) := { v ∈ H 1 ⁢ ( Ω ) ; v | ∂ ⁡ T ∈ P p ⁢ ( F T ) ⁢ for all ⁢ T ∈ T } , P p , c ⁢ ( S ) := tr grad ⁡ ( U grad , p ⁢ ( T ) ) , P 0 p , c ( S ) := tr grad ( U grad , p ( T ) ∩ H 0 1 ( Ω ) ) .

Then we consider the approximation space

(3.5) U h := U h p := P p ⁢ ( T ) d × P p ⁢ ( T ) × P p ⁢ ( T ) d × P p ⁢ ( T ) × P 0 p ⁢ ( S ) × P 0 p + 1 , c ⁢ ( S ) × P p ⁢ ( S ) × P p + 1 , c ⁢ ( S )

and discrete test space

(3.6) V h ⁢ ( T ) := V h p ⁢ ( T ) := P p + 2 ⁢ ( T ) d × P p + d + 1 ⁢ ( T ) × P p + 2 ⁢ ( T ) d × P p + d + 1 ⁢ ( T ) .

We employ the operators Π p + 2 div and Π p + d + 1 grad from [21, Section 3.4], and using the results reported there, it is clear that

Π F ⁢ ( v ) := ( Π p + 2 div ⁢ v 1 , Π p + d + 1 grad ⁢ v 2 , Π p + 2 div ⁢ v 3 , Π p + d + 1 grad ⁢ v 4 )

for v = ( v 1 , v 2 , v 3 , v 4 ) ∈ V ⁢ ( T ) is a Fortin operator as required for our case. For details, we refer to [21, Lemmas 3.2, 3.3]. Note the increased polynomial degree by d + 1 instead of 𝑑 as in [21] for the approximation of H 1 ⁢ ( T ) . This is necessary due to the presence of the terms ( u 2 , v 4 ) and ( u 4 , v 2 ) in the bilinear form, contrary to the Poisson case in [21] without reaction term.

Combining the established existence of a Fortin operator with Theorem 8, we conclude the quasi-optimal convergence of the fully-discrete DPG scheme; cf. [21, Theorem 2.1].

Corollary 9

Let f ∈ L 2 ⁢ ( Ω ) and p ∈ N 0 be given. Selecting the approximation and test spaces U h , V h ⁢ ( T ) as in (3.5), (3.6), respectively, system (3.4) has a unique solution u h ∈ U h . It satisfies the quasi-optimal error estimate

∥ u - u h ∥ U ≲ ∥ u - w ∥ U for all ⁢ w ∈ U h

with a hidden constant that is independent of 𝒯.

3.2 Proofs for the First-Order Formulation

In this section, we prove Theorems 7, 8. This requires some preliminary steps given next.

3.2.1 Stability of the Adjoint Problem

The adjoint problem to (3.1) with continuous test functions

( v 1 , v 2 , v 3 , v 4 ) ∈ V 0 ⁢ ( Ω ) := H 0 ⁢ ( div , Ω ) × H 0 1 ⁢ ( Ω ) × H ⁢ ( div , Ω ) × H 1 ⁢ ( Ω )

reads as follows: given g 1 , g 3 ∈ L 2 ⁢ ( Ω ) and g 2 , g 4 ∈ L 2 ⁢ ( Ω ) , find ( v 1 , v 2 , v 3 , v 4 ) ∈ V 0 ⁢ ( Ω ) such that

(3.7) v 1 - ∇ ⁡ v 4 = g 1 , v 4 + div ⁡ v 3 = g 2 , v 3 + ∇ ⁡ v 2 = g 3 , v 2 + div ⁡ v 1 = g 4 .

It is a fully non-homogeneous version of (3.1) with some sign changes. We show that this problem is well posed.

Lemma 10

There is a unique solution v = ( v 1 , v 2 , v 3 , v 4 ) ∈ V 0 ⁢ ( Ω ) to (3.7), and the bound

∥ v ∥ V ⁢ ( T ) 2 ≲ ∥ g 1 ∥ 2 + ∥ g 2 ∥ 2 + ∥ g 3 ∥ 2 + ∥ g 4 ∥ 2

holds with a constant that is independent of the given data.

Proof

Eliminating v 1 and v 3 , (3.7) becomes

v 4 + div ⁢ ( g 3 - ∇ ⁡ v 2 ) = g 2 , v 2 + div ⁢ ( g 1 + ∇ ⁡ v 4 ) = g 4 ,

and in variational form: find ( v 2 , v 4 ) ∈ H 0 1 ⁢ ( Ω ) × H 1 ⁢ ( Ω ) such that

(3.8) a ⁢ ( v 2 , v 4 ; δ ⁢ v 2 , δ ⁢ v 4 ) = ( g 2 , δ ⁢ v 2 ) + ( g 3 , ∇ ⁡ δ ⁢ v 2 ) + ( g 4 , δ ⁢ v 4 ) + ( g 1 , ∇ ⁡ δ ⁢ v 4 )

for any ( δ ⁢ v 2 , δ ⁢ v 4 ) ∈ H 0 1 ⁢ ( Ω ) × H 1 ⁢ ( Ω ) , where

a ⁢ ( v 2 , v 4 ; δ ⁢ v 2 , δ ⁢ v 4 ) := ( ∇ ⁡ v 2 , ∇ ⁡ δ ⁢ v 2 ) + ( v 4 , δ ⁢ v 2 ) + ( v 2 , δ ⁢ v 4 ) - ( ∇ ⁡ v 4 , ∇ ⁡ δ ⁢ v 4 ) .

Note that the discarded duality between the traces on Γ of ( g 1 + ∇ ⁡ v 4 ) ⋅ n and δ ⁢ v 4 implies the required homogeneous boundary condition for v 1 := g 1 + ∇ ⁡ v 4 . Problem (3.8) is a mixed system of saddle point form. It is somewhat non-standard due to the presence of the negative semi-definite form - ( ∇ ⁡ v 4 , ∇ ⁡ δ ⁢ v 4 ) . We show its well-posedness. Then v := ( v 1 , v 2 , v 3 , v 4 ) with v 1 := g 1 + ∇ ⁡ v 4 and v 3 := g 3 - ∇ ⁡ v 2 satisfies v ∈ V 0 ⁢ ( Ω ) , solves (3.7) and is bounded as claimed.

Noting the boundedness of the appearing linear and bilinear forms, to show well-posedness of (3.8), it is enough to verify the inf-sup condition

(3.9) sup 0 ≠ ( δ ⁢ v 2 , δ ⁢ v 4 ) ∈ H 0 1 ⁢ ( Ω ) × H 1 ⁢ ( Ω ) ⁡ a ⁢ ( v 2 , v 4 ; δ ⁢ v 2 , δ ⁢ v 4 ) ∥ ∇ ⁡ δ ⁢ v 2 ∥ + ∥ δ ⁢ v 4 ∥ 1 ≳ ∥ ∇ ⁡ v 2 ∥ + ∥ v 4 ∥ 1 for all ⁢ ( v 2 , v 4 ) ∈ H 0 1 ⁢ ( Ω ) × H 1 ⁢ ( Ω )

and injectivity

a ⁢ ( δ ⁢ v 2 , δ ⁢ v 4 ; v 2 , v 4 ) = 0 ⁢ for all ⁢ ( δ ⁢ v 2 , δ ⁢ v 4 ) ∈ H 0 1 ⁢ ( Ω ) × H 1 ⁢ ( Ω ) ⟹ v 2 = v 4 = 0

for ( v 2 , v 4 ) ∈ H 0 1 ⁢ ( Ω ) × H 1 ⁢ ( Ω ) . The injectivity is immediate since a ⁢ ( v 2 , - v 4 ; v 2 , v 4 ) = ∥ ∇ ⁡ v 2 ∥ 2 + ∥ ∇ ⁡ v 4 ∥ 2 = 0 implies that v 2 = 0 , v 4 = c ∈ R , and a ⁢ ( δ ⁢ v 2 , 0 ; 0 , c ) = ( δ ⁢ v 2 , c ) = 0 for all δ ⁢ v 2 ∈ H 0 1 ⁢ ( Ω ) means that c = 0 .

We are left with proving (3.9). To this end, we define

∥ ( v 2 , v 4 ) ∥ a 2 := ∥ ∇ ⁡ v 2 ∥ 2 + ∥ ∇ ⁡ v 4 ∥ 2 + | a ⁢ ( v 2 , v 4 ; w , 0 ) | 2 , ( v 2 , v 4 ) ∈ H 0 1 ⁢ ( Ω ) × H 1 ⁢ ( Ω ) ,

where w ∈ H 0 1 ⁢ ( Ω ) is a fixed function solving - Δ ⁢ w = 1 in Ω. One verifies that ∥ ⋅ ∥ a is a norm in H 0 1 ⁢ ( Ω ) × H 1 ⁢ ( Ω ) and that the equivalence

∥ ( v 2 , v 4 ) ∥ a + ∥ v 4 ∥ ≃ ∥ ∇ ⁡ v 2 ∥ + ∥ ∇ ⁡ v 4 ∥ + ∥ v 4 ∥ for all ⁢ ( v 2 , v 4 ) ∈ H 0 1 ⁢ ( Ω ) × H 1 ⁢ ( Ω )

holds true. By Rellich’s compactness of the embedding H 1 ⁢ ( Ω ) → L 2 ⁢ ( Ω ) and the Peetre–Tartar lemma (see, e.g., [20, Theorem I.2.1]), we also have the equivalence

∥ ( v 2 , v 4 ) ∥ a ≃ ∥ ∇ ⁡ v 2 ∥ + ∥ ∇ ⁡ v 4 ∥ + ∥ v 4 ∥ for all ⁢ ( v 2 , v 4 ) ∈ H 0 1 ⁢ ( Ω ) × H 1 ⁢ ( Ω ) .

The inf-sup property follows by selecting δ ⁢ v 2 := v 2 + a ⁢ ( v 2 , v 4 ; w , 0 ) ⁢ w , δ ⁢ v 4 := - v 4 , calculating

a ⁢ ( v 2 , v 4 ; δ ⁢ v 2 , δ ⁢ v 4 ) = ∥ ∇ ⁡ v 2 ∥ 2 + ∥ ∇ ⁡ v 4 ∥ 2 + | a ⁢ ( v 2 , v 4 ; w , 0 ) | 2 = ∥ ( v 2 , v 4 ) ∥ a 2

and bounding

∥ ∇ ⁡ δ ⁢ v 2 ∥ 2 + ∥ ∇ ⁡ δ ⁢ v 4 ∥ 2 + ∥ δ ⁢ v 4 ∥ 2 ≲ ∥ ∇ ⁡ v 2 ∥ 2 + ∥ v 4 ∥ 1 2 + | a ⁢ ( v 2 , v 4 ; w , 0 ) | 2 ⁢ ∥ ∇ ⁡ w ∥ 2 ≲ ∥ ∇ ⁡ v 2 ∥ 2 + ∥ v 4 ∥ 1 2 .

The hidden constants only depend on Ω. ∎

Lemma 11

The adjoint operator B * : V ⁢ ( T ) → U ′ is injective.

Proof

The proof is standard. We have to show that, if v ∈ V ⁢ ( T ) satisfies b ⁢ ( δ ⁢ u , v ) = 0 for all δ ⁢ u ∈ U , then v = 0 . Selecting functions δ ⁢ u with only one of the trace components non-zero and all the field variables zero, Proposition 1 shows that v ∈ V 0 ⁢ ( Ω ) . Furthermore, 𝔳 solves (3.7) with homogeneous data so that v = 0 by Lemma 10. ∎

3.2.2 Proofs of Theorems 7, 8

The well-posedness of (3.2) follows by standard arguments.

Boundedness of the functional and bilinear form. By definition of the norms in 𝔘 and V ⁢ ( T ) , and using Proposition 2, the uniform boundedness of b ⁢ ( ⋅ , ⋅ ) and 𝐿 is immediate.
Injectivity. This is Lemma 11.
Inf-sup condition. By [4, Theorem 3.3], the inf-sup condition
(3.10) sup 0 ≠ v ∈ V ⁢ ( T ) ⁡ b ⁢ ( u ; v ) ∥ v ∥ V ⁢ ( T ) ≳ ∥ u ∥ U for all ⁢ u = ( u 0 , u ^ ) ∈ U 0 × U ^ = U
follows from the inf-sup conditions
(3.11) sup 0 ≠ v ∈ V ⁢ ( T ) ⁡ b ⁢ ( 0 , u ^ ; v ) ∥ v ∥ V ⁢ ( T ) ≳ ∥ u ^ ∥ U ^ for all ⁢ u ^ ∈ U ^ ,

(3.12) sup 0 ≠ v ∈ V 0 ⁢ ( Ω ) ⁡ b ⁢ ( u 0 , 0 ; v ) ∥ v ∥ V ⁢ ( T ) ≳ ∥ u 0 ∥ for all ⁢ u 0 ∈ U 0 .
We note that Proposition 1 implies that V 0 ⁢ ( Ω ) is the correct test space in (3.12). It corresponds to Y 0 in [4]; see [4, relation (17)]. Relation (3.11) is true (with constant 1) by Proposition 2, and inf-sup condition (3.12) holds due to Lemma 10: selecting data ( g 1 , g 2 , g 3 , g 4 ) := ( u 1 , - u 2 , - u 3 , - u 4 ) in (3.7) with solution v * ∈ V 0 ⁢ ( Ω ) , we bound
sup 0 ≠ v ∈ V 0 ⁢ ( Ω ) ⁡ b ⁢ ( u 0 , 0 ; v ) ∥ v ∥ V ⁢ ( T ) ≥ ( u 1 , g 1 ) - ( u 2 , g 2 ) - ( u 3 , g 3 ) - ( u 4 , g 4 ) ∥ v * ∥ V ⁢ ( T ) ≳ ∥ u 0 ∥ .

We conclude that there is a unique and stable solution u = ( u 0 , u ^ ) ∈ U 0 × U ^ of (3.2). It remains to note that u 0 solves (3.1), satisfies the homogeneous boundary conditions, and

u ^ = ( tr div ⁡ ( u 1 ) , tr grad ⁡ ( u 2 ) , tr div ⁡ ( u 3 ) , tr grad ⁡ ( u 4 ) ) .

Furthermore, u := u 1 solves (1.1). This finishes the proof of Theorem 7.

Theorem 8 follows in the standard way. Note that the DPG scheme minimizes the residual ∥ B ⁢ ( u - u h ) ∥ V ⁢ ( T ) ′ . Then its quasi-optimal convergence in the norm ∥ ⋅ ∥ U is due to the equivalence of this norm and the residual norm ∥ B ⋅ ∥ V ⁢ ( T ) ′ . Specifically, the boundedness ∥ B ⋅ ∥ V ⁢ ( T ) ′ ≲ ∥ ⋅ ∥ U is due to the boundedness of b ⁢ ( ⋅ , ⋅ ) , and the inverse estimate is proved by inf-sup property (3.10).

4 Second-Order Variational Formulation and DPG Method

In this section, we use the same notation as in Section 3 for spaces, norms, bilinear form, functional 𝐿, operator 𝐵 and trial-to-test operator 𝔗.

Denoting w := - ∇ ⁡ div ⁡ u , we write (1.1) as the second-order system (4.1)

(4.1a) - ∇ ⁡ div ⁡ w + u = f ,

(4.1b) ∇ ⁡ div ⁡ u + w = 0

in Ω subject to u ⋅ n = div ⁡ u = 0 on Γ. Testing (4.1a) and (4.1b) with piecewise smooth functions 𝒗 and - τ , respectively, and employing the trace operator tr ∇ ⁡ div , we obtain

( u , v - ∇ ⁡ div ⁡ τ ) T - ( w , τ + ∇ ⁡ div ⁡ v ) T + ⟨ tr ∇ ⁡ div ⁡ ( u ) , τ ⟩ S ∇ ⁡ div + ⟨ tr ∇ ⁡ div ⁡ ( w ) , v ⟩ S ∇ ⁡ div = ( f , v ) .

We introduce independent trace variables u ^ = tr ∇ ⁡ div ⁡ ( u ) , w ^ = tr ∇ ⁡ div ⁡ ( w ) and define the spaces

U := U 0 × U ^ := ( L 2 ⁢ ( Ω ) × L 2 ⁢ ( Ω ) ) × ( H 0 ∇ ⁡ div ⁢ ( S ) × H ∇ ⁡ div ⁢ ( S ) ) , V ⁢ ( T ) := H ⁢ ( ∇ ⁡ div , T ) × H ⁢ ( ∇ ⁡ div , T )

with (squared) norms

∥ ( u , w , u ^ , w ^ ) ∥ U 2 := ∥ ( u , w ) ∥ 2 + ∥ ( u ^ , w ^ ) ∥ U ^ 2 := ∥ u ∥ 2 + ∥ w ∥ 2 + ∥ u ^ ∥ - 1 / 2 , 1 / 2 , S 2 + ∥ w ^ ∥ - 1 / 2 , 1 / 2 , S 2 , ∥ ( v , τ ) ∥ V ⁢ ( T ) 2 := ∥ v ∥ ∇ ⁡ div , T 2 + ∥ τ ∥ ∇ ⁡ div , T 2 ,

respectively. Then our second variational formulation of (1.1) reads: find u = ( u , w , u ^ , w ^ ) ∈ U such that

(4.2) b ⁢ ( u , v ) = L ⁢ ( v ) for all ⁢ v = ( v , τ ) ∈ V ⁢ ( T ) ,

in operator form u ∈ U : B ⁢ u = L in V ⁢ ( T ) ′ . Here,

b ⁢ ( u , v ) := ( u , v - ∇ ⁡ div ⁡ τ ) T - ( w , τ + ∇ ⁡ div ⁡ v ) T + ⟨ u ^ , τ ⟩ S ∇ ⁡ div + ⟨ w ^ , v ⟩ S ∇ ⁡ div

and L ⁢ ( v ) := ( f , v ) .

Theorem 12

For given f ∈ L 2 ⁢ ( Ω ) , there exists a unique solution u = ( u , w , u ^ , w ^ ) ∈ U to (4.2). It is uniformly bounded, ∥ u ∥ U ≲ ∥ f ∥ with a hidden constant that is independent of 𝒇 and 𝒯. Furthermore, 𝒖 solves (1.1).

For a proof of this theorem, we refer to § 4.2.2.

To define the discretization scheme, we consider as before a (family of) discrete subspace(s) U h ⊂ U and use the trial-to-test operator T : U → V ⁢ ( T ) defined by ⟨ ⟨ T ⁢ ( u ) , v ⟩ ⟩ V ⁢ ( T ) = b ⁢ ( u , v ) for all v ∈ V ⁢ ( T ) , where ⟨ ⟨ ⋅ , ⋅ ⟩ ⟩ V ⁢ ( T ) denotes the inner product of the new space V ⁢ ( T ) ,

⟨ ⟨ ( v 1 , v 2 ) , ( δ ⁢ v 1 , δ ⁢ v 2 ) ⟩ ⟩ V ⁢ ( T ) := ∑ j = 1 2 ( v j , δ ⁢ v j ) + ( ∇ ⁡ div ⁡ v j , ∇ ⁡ div ⁡ δ ⁢ v j ) T

for ( v 1 , v 2 ) , ( δ ⁢ v 1 , δ ⁢ v 2 ) ∈ V h ⁢ ( T ) . The second DPG scheme for problem (1.1) is: find u h ∈ U h such that

(4.3) b ⁢ ( u h , T ⁢ δ ⁢ u ) = L ⁢ ( T ⁢ δ ⁢ u ) for all ⁢ δ ⁢ u ∈ U h .

Again, it delivers a quasi-best approximation.

Theorem 13

Let f ∈ L 2 ⁢ ( Ω ) be given. For any finite-dimensional subspace U h ⊂ U , there exists a unique solution u h ∈ U h to (4.3). It satisfies the quasi-optimal error estimate ∥ u - u h ∥ U ≲ ∥ u - w ∥ U for all w ∈ U h with a hidden constant that is independent of 𝒯.

A proof of this theorem is given in § 4.2.2.

4.1 Fully Discrete Method

In the case of the second-order formulation, we restrict our fully discrete analysis to two space dimensions ( d = 2 ) and lowest-order approximations. Using the notation for (piecewise) polynomial spaces from § 3.1, our approximation space is

(4.4) U h := P 0 ⁢ ( T ) 2 × P 0 ⁢ ( T ) 2 × P 0 0 ⁢ ( S ) × P 0 1 , c ⁢ ( S ) × P 0 ⁢ ( S ) × P 1 , c ⁢ ( S ) ,

and as enriched test space, we choose

(4.5) V h ⁢ ( T ) := P 3 ⁢ ( T ) 2 × P 3 ⁢ ( T ) 2 .

Here, 𝒯 denotes a (family of) shape-regular triangular mesh(es). The notation for our approximation scheme is as before

(4.6) u h ∈ U h : b ⁢ ( u h , T h ⁢ δ ⁢ u ) = L ⁢ ( T h ⁢ δ ⁢ u ) for all ⁢ δ ⁢ u ∈ U h ,

where the approximated trial-to-test operator T h : U h → V h ⁢ ( T ) is defined accordingly.

The existence of a Fortin operator in this case is due to Lemma 18 (given in § 4.2.3 below), showing the existence of a local variant Π T ∇ ⁡ div : H ⁢ ( ∇ ⁡ div , T ) → P 3 ⁢ ( T ) 2 . We simply define Π ∇ ⁡ div : H ⁢ ( ∇ ⁡ div , T ) → P 3 ⁢ ( T ) 2 as ( Π ∇ ⁡ div ⁢ v ) | T := Π T ∇ ⁡ div ⁢ v | T ( T ∈ T ), and Π F ⁢ ( v ) := ( Π ∇ ⁡ div ⁢ v , Π ∇ ⁡ div ⁢ τ ) for v = ( v , τ ) ∈ V ⁢ ( T ) . The statements of Lemma 18 imply that this operator has the properties of a Fortin operator. Then, analogously to Corollary 9, we conclude the quasi-optimal convergence of the fully discrete DPG scheme.

Corollary 14

Let f ∈ L 2 ⁢ ( Ω ) be given. Selecting the approximation and test spaces U h , V h ⁢ ( T ) as in (4.4), (4.5), respectively, system (4.6) has a unique solution u h ∈ U h . It satisfies the quasi-optimal error estimate

∥ u - u h ∥ U ≲ ∥ u - w ∥ U for all ⁢ w ∈ U h

with a hidden constant that is independent of 𝒯.

4.2 Proofs for the Second-Order Formulation

The steps in this section are analogous to the ones in § 3.2.2.

4.2.1 Stability of the Adjoint Problem

The adjoint problem to (4.1) with continuous test functions ( v , τ ) ∈ V 0 ⁢ ( Ω ) := H 0 ⁢ ( ∇ ⁡ div , Ω ) × H ⁢ ( ∇ ⁡ div , Ω ) reads as follows: given g 1 , g 2 ∈ L 2 ⁢ ( Ω ) , find ( v , τ ) ∈ V 0 ⁢ ( Ω ) such that

(4.7) v - ∇ ⁡ div ⁡ τ = g 1 , τ + ∇ ⁡ div ⁡ v = g 2 .

It is a well-posed, fully non-homogeneous version of (4.1).

Lemma 15

There is a unique solution v = ( v , τ ) ∈ V 0 ⁢ ( Ω ) to (4.7), and the bound ∥ v ∥ V ⁢ ( T ) 2 ≲ ∥ g 1 ∥ 2 + ∥ g 2 ∥ 2 holds with a constant that is independent of the given data.

Proof

Eliminating 𝒗, (4.7) becomes τ + ∇ ⁡ div ⁢ ( g 1 + ∇ ⁡ div ⁡ τ ) = g 2 , in variational form: find τ ∈ H ⁢ ( ∇ ⁡ div , Ω ) such that

(4.8) ( ∇ ⁡ div ⁡ τ , ∇ ⁡ div ⁡ δ ⁢ τ ) + ( τ , δ ⁢ τ ) = ( g 2 , δ ⁢ τ ) - ( g 1 , ∇ ⁡ div ⁡ δ ⁢ τ ) for all ⁢ δ ⁢ τ ∈ H ⁢ ( ∇ ⁡ div , Ω ) .

As for the first-order system (see the proof of Lemma 10), the discarded duality between the trace

tr ∇ ⁡ div ⁡ ( g 1 + ∇ ⁡ div ⁡ τ )

and δ ⁢ τ implies the homogeneous boundary condition for v := g 1 + ∇ ⁡ div ⁡ τ . By the Lax–Milgram lemma, problem (4.8) is well posed. Then ( v , τ ) with v := g 1 + ∇ ⁡ div ⁡ τ solves (4.7) and satisfies the claimed bound. ∎

Lemma 16

The adjoint operator B * : V ⁢ ( T ) → U ′ is injective.

Proof

The proof is identical to the one of Lemma 11, using Proposition 6 and Lemma 15 instead of Proposition 1 and Lemma 10, respectively. ∎

4.2.2 Proofs of Theorems 12, 13

To show the well-posedness of (4.2), we repeat the arguments from § 3.2.2. The boundedness of the functional and bilinear form is identical to the previous case. The injectivity of B * holds by Lemma 16. In order to prove the inf-sup condition, it is enough to show

(4.9) sup 0 ≠ v ∈ V ⁢ ( T ) ⁡ b ⁢ ( 0 , u ^ ; v ) ∥ v ∥ V ⁢ ( T ) ≳ ∥ u ^ ∥ U ^ for all ⁢ u ^ ∈ U ^ ,

(4.10) sup 0 ≠ v ∈ V 0 ⁢ ( Ω ) ⁡ b ⁢ ( u 0 , 0 ; v ) ∥ v ∥ V ⁢ ( T ) ≳ ∥ u 0 ∥ for all ⁢ u 0 ∈ U 0 .

The test space in (4.10) is the appropriate one due to Proposition 6.

Relation (4.9) is true (with constant 1) by Proposition 5, and inf-sup condition (4.10) can be seen by using Lemma 15 with data ( g 1 , g 2 ) := ( u , - w ) in (4.7). We conclude that there is a unique and stable solution u = ( u 0 , u ^ ) ∈ U 0 × U ^ of (4.2). It remains to note that ( u , w ) solves (4.1), 𝒖 satisfies the homogeneous boundary conditions, and u ^ = tr ∇ ⁡ div ⁡ ( u ) , w ^ = tr ∇ ⁡ div ⁡ ( w ) . Furthermore, 𝒖 solves (1.1). This finishes the proof of Theorem 12.

The proof of Theorem 13 is identical to the one of Theorem 8.

4.2.3 Fortin Operator

Let T ^ denote the reference triangle with vertices ( 0 , 0 ) , ( 1 , 0 ) , ( 0 , 1 ) . We start with a preliminary Fortin operator Π ~ T ^ ∇ ⁡ div on the reference element. The corresponding operator on an element T ∈ T is not uniformly bounded with respect to diam ⁢ ( T ) so that it has to be adapted. This will be studied in Lemma 18.

Lemma 17

There is a linear bounded operator Π ~ T ^ ∇ ⁡ div : H ⁢ ( ∇ ⁡ div , T ^ ) → P 3 ⁢ ( T ^ ) 2 that satisfies, for any w ∈ P 0 ⁢ ( T ^ ) 2 and w ^ ∈ P 0 ⁢ ( ∂ ⁡ T ^ ) × P 1 , c ⁢ ( ∂ ⁡ T ^ ) ,

( w , ( 1 - Π ~ T ^ ∇ ⁡ div ) ⁢ v ) T ^ = ⟨ w ^ , ( 1 - Π ~ T ^ ∇ ⁡ div ) ⁢ v ⟩ ∂ ⁡ T ^ ∇ ⁡ div = 0 for all ⁢ v ∈ H ⁢ ( ∇ ⁡ div , T ^ ) .

Proof

For the construction of this operator, we employ the strategy of Nagaraj, Petrides and Demkowicz [28]; see also [16]. For given v ∈ H ⁢ ( ∇ ⁡ div , T ^ ) , we define Π ~ T ^ ∇ ⁡ div ⁢ v := v * , where

( v * , w , v ^ ) ∈ P 3 ⁢ ( T ^ ) 2 × P 0 ⁢ ( T ^ ) 2 × ( P 0 ⁢ ( ∂ ⁡ T ^ ) × P 1 , c ⁢ ( ∂ ⁡ T ^ ) )

is the solution to the mixed problem

⟨ ⟨ v * , δ ⁢ v ⟩ ⟩ ∇ ⁡ div , T ^ + ( w , δ ⁢ v ) T ^ + ⟨ v ^ , δ ⁢ v ⟩ ∂ ⁡ T ^ ∇ ⁡ div = 0 for all ⁢ δ ⁢ v ∈ P 3 ⁢ ( T ^ ) 2 , ( δ ⁢ u , v * ) T ^ = ( δ ⁢ u , v ) T ^ for all ⁢ δ ⁢ u ∈ P 0 ⁢ ( T ^ ) 2 , ⟨ δ ⁢ u ^ , v * ⟩ ∂ ⁡ T ^ ∇ ⁡ div = ⟨ δ ⁢ u ^ , v ⟩ ∂ ⁡ T ^ ∇ ⁡ div for all ⁢ δ ⁢ u ^ ∈ P 0 ⁢ ( ∂ ⁡ T ^ ) × P 1 , c ⁢ ( ∂ ⁡ T ^ ) .

Using a canonical basis, one can show that the matrix of this problem is invertible, thus giving the result. ∎

Lemma 18

For T ∈ T , there is a linear operator Π T ∇ ⁡ div : H ⁢ ( ∇ ⁡ div , T ) → P 3 ⁢ ( T ) 2 that satisfies, for any w ∈ P 0 ⁢ ( T ) 2 and w ^ ∈ P 0 ⁢ ( ∂ ⁡ T ) × P 1 , c ⁢ ( ∂ ⁡ T ) ,

( w , ( 1 - Π T ∇ ⁡ div ) ⁢ v ) T = ( w , ∇ ⁡ div ⁢ ( 1 - Π T ∇ ⁡ div ) ⁢ v ) T = ⟨ w ^ , ( 1 - Π T ∇ ⁡ div ) ⁢ v ⟩ ∂ ⁡ T ∇ ⁡ div = 0 for all ⁢ v ∈ H ⁢ ( ∇ ⁡ div , T ) ,

and which is uniformly bounded with respect to diam ⁢ ( T ) ,

∥ Π T ∇ ⁡ div ⁢ v ∥ ∇ ⁡ div , T ≲ ∥ v ∥ ∇ ⁡ div , T for all ⁢ v ∈ H ⁢ ( ∇ ⁡ div , T ) , T ∈ T .

Proof

In this proof, Π p denotes the L 2 -projection onto polynomials of degree 𝑝 (on the domain under consideration). Let T ∈ T and v ∈ H ⁢ ( ∇ ⁡ div , T ) be given. Furthermore, let T ^ ∋ x ^ ↦ B T ⁢ x ^ + a T ∈ T be an affine mapping with a T ∈ R 2 , B T ∈ R 2 × 2 , and J T := det ⁡ ( B T ) .

We will define Π T ∇ ⁡ div ⁢ v by splitting 𝒗 into three components and treat them differently. To this end, we need three operators: ℐ in step (i) below, Π 3 div in step (ii), already used in § 3.1, and Π ~ T ∇ ⁡ div in step (iii), based on Π ~ T ^ ∇ ⁡ div .

(i) Following Ern et al. [11, Theorem 3.2], there is an interpolation operator

I : H ⁢ ( div , T ) → R ⁢ T 0 ⁢ ( T ) := { z : T → R 2 ; there exist ⁢ a , b , c ∈ R : z = ( a , b ) T + c ⁢ x ⁢ for all ⁢ x ∈ T }

with div ∘ I = Π 0 ∘ div on H ⁢ ( div , T ) and boundedness

(4.11) ∥ I ⁢ z ∥ T ≲ ∥ z ∥ T + diam ⁢ ( T ) ⁢ ∥ ( 1 - Π 0 ) ⁢ div ⁡ z ∥ T for all ⁢ z ∈ H ⁢ ( div , T )

uniformly with respect to diam ⁢ ( T ) . Note that ∥ ∇ ⁡ div ⁡ T ⁢ z ∥ T = 0 for all z ∈ H ⁢ ( ∇ ⁡ div , T ) .

Then we write v = I ⁢ v + w with w := ( 1 - I ) ⁢ v and use the Helmholtz decomposition w = ∇ ⁡ ϕ + curl ⁡ ψ , where ϕ ∈ H 0 1 ⁢ ( T ) solves Δ ⁢ ϕ = div ⁡ w , and ψ ∈ H ⁢ ( curl , T ) (with obvious definition of the space). Integration by parts, the commuting property of ℐ, the Bramble–Hilbert lemma and Poincare’s inequality prove that

∥ ∇ ⁡ ϕ ∥ T 2 = - ( div ⁡ w , ϕ ) T = - ( ( 1 - Π 0 ) ⁢ div ⁡ v , ϕ ) T ≲ ∥ ( 1 - Π 0 ) ⁢ div ⁡ v ∥ T ⁢ ∥ ϕ ∥ T ≲ diam ⁢ ( T ) 2 ⁢ ∥ ∇ ⁡ div ⁡ v ∥ T ⁢ ∥ ∇ ⁡ ϕ ∥ T ,

(4.12) i.e., ∥ ∇ ⁡ ϕ ∥ T ≲ diam ⁢ ( T ) 2 ⁢ ∥ ∇ ⁡ div ⁡ v ∥ T .

(ii) We continue to use the Fortin operator Π p div : H ⁢ ( div , T ) → P p ⁢ ( T ) 2 from [21, Lemma 3.3] already employed in § 3.1, in this case for p = 3 . It is uniformly bounded in the H ⁢ ( div , T ) -norm and has the properties

( δ ⁢ u , Π 3 div ⁢ z ) T = ( δ ⁢ u , z ) T , ⟨ δ ⁢ u ^ , Π 3 div ⁢ z ⟩ ∂ ⁡ T ∇ = ⟨ δ ⁢ u ^ , z ⟩ ∂ ⁡ T ∇ for all ⁢ z ∈ H ⁢ ( div , T ) ,

for any δ ⁢ u ∈ P 1 ⁢ ( T ) 2 and δ ⁢ u ^ ∈ P 2 , c ⁢ ( ∂ ⁡ T ) . Furthermore, div ∘ Π 3 div = Π 2 ∘ div ; see [21, proof of Lemma 3.3]. On the one hand, this property implies that

(4.13) ⟨ δ ⁢ u ^ , Π 3 div ⁢ curl ⁡ ψ ⟩ ∂ ⁡ T ∇ ⁡ div = ⟨ δ ⁢ u ^ , curl ⁡ ψ ⟩ ∂ ⁡ T ∇ ⁡ div for all ⁢ δ ⁢ u ^ ∈ P 0 ⁢ ( ∂ ⁡ T ) × P 1 , c ⁢ ( ∂ ⁡ T ) .

On the other hand, also using (4.11) and the Bramble–Hilbert lemma, we conclude that

∥ Π 3 div ⁢ curl ⁡ ψ ∥ T ≲ ∥ curl ⁡ ψ ∥ T + ∥ div ⁢ curl ⁡ ψ ∥ T ≤ ∥ w ∥ T ≲ ∥ v ∥ T + diam ⁢ ( T ) 2 ⁢ ∥ ∇ ⁡ div ⁡ v ∥ T , ∥ ∇ ⁡ div ⁡ Π 3 div ⁢ curl ⁡ ψ ∥ T = ∥ ∇ ⁡ Π 2 ⁢ div ⁢ curl ⁡ ψ ∥ T = 0 .

(iii) We transform the operator Π ~ T ^ ∇ ⁡ div from Lemma 17 through the Piola map T ^ → T . Denoting by z ^ ∈ H ⁢ ( ∇ ⁡ div , T ^ ) the preimage of z ∈ H ⁢ ( ∇ ⁡ div , T ) under the Piola transform, we define

Π ~ T ∇ ⁡ div ⁢ z := | J T | - 1 ⁢ B T ⁢ Π ~ T ^ ∇ ⁡ div ⁢ z ^ ∘ F T - 1 .

By scaling and Lemma 17, we bound

∥ Π ~ T ∇ ⁡ div ⁢ z ∥ T ≃ ∥ Π ~ T ^ ∇ ⁡ div ⁢ z ^ ∥ T ^ ≲ ∥ z ∥ T + diam ⁢ ( T ) 2 ⁢ ∥ ∇ ⁡ div ⁡ z ∥ T , ∥ ∇ ⁡ div ⁡ Π ~ T ∇ ⁡ div ⁢ z ∥ T ≃ diam ⁢ ( T ) - 2 ⁢ ∥ ∇ ⁡ div ⁡ Π ~ T ^ ∇ ⁡ div ⁢ z ^ ∥ T ^ ≲ diam ⁢ ( T ) - 2 ⁢ ∥ z ∥ T + ∥ ∇ ⁡ div ⁡ z ∥ T .

These estimates together with (4.12), and noting that ∇ ⁡ div ⁡ w = ∇ ⁡ div ⁢ ( 1 - I ) ⁢ v = ∇ ⁡ div ⁡ v , yield

∥ Π ~ T ∇ ⁡ div ⁢ ∇ ⁡ ϕ ∥ T ≲ ∥ ∇ ⁡ ϕ ∥ T + diam ⁢ ( T ) 2 ⁢ ∥ ∇ ⁡ div ⁡ w ∥ T ≲ diam ⁢ ( T ) 2 ⁢ ∥ ∇ ⁡ div ⁡ v ∥ T , ∥ ∇ ⁡ div ⁡ Π ~ T ∇ ⁡ div ⁢ ∇ ⁡ ϕ ∥ T ≲ diam ⁢ ( T ) - 2 ⁢ ∥ ∇ ⁡ ϕ ∥ T + ∥ ∇ ⁡ div ⁡ w ∥ T ≲ ∥ ∇ ⁡ div ⁡ v ∥ T .

Definition of Π T ∇ ⁡ div ⁢ v . Finally, we define Π T ∇ ⁡ div ⁢ v := I ⁢ v + Π 3 div ⁢ curl ⁡ ψ + Π ~ T ∇ ⁡ div ⁢ ∇ ⁡ ϕ . By steps (i)–(iii), it is clear that this operator is uniformly bounded. It remains to verify the orthogonality relations. By the orthogonality properties of Π 3 div and Π ~ T ∇ ⁡ div , we immediately verify that

( δ ⁢ u , Π T ∇ ⁡ div ⁢ v ) T = ( δ ⁢ u , I ⁢ v + Π 3 div ⁢ curl ⁡ ψ + Π ~ T ∇ ⁡ div ⁢ ∇ ⁡ ϕ ) T = ( δ ⁢ u , I ⁢ v + curl ⁡ ψ + ∇ ⁡ ϕ ) T = ( δ ⁢ u , v ) T

for any δ ⁢ u ∈ P 0 ⁢ ( T ) 2 and, recalling Lemma 17 and (4.13),

⟨ δ ⁢ u ^ , Π T ∇ ⁡ div ⁢ v ⟩ ∂ ⁡ T ∇ ⁡ div = ⟨ δ ⁢ u ^ , I ⁢ v + Π 3 div ⁢ curl ⁡ ψ + Π ~ T ∇ ⁡ div ⁢ ∇ ⁡ ϕ ⟩ ∂ ⁡ T ∇ ⁡ div = ⟨ δ ⁢ u ^ , I ⁢ v + curl ⁡ ψ + ∇ ⁡ ϕ ⟩ ∂ ⁡ T ∇ ⁡ div = ⟨ δ ⁢ u ^ , v ⟩ ∂ ⁡ T ∇ ⁡ div

for any δ ⁢ u ^ ∈ P 0 ⁢ ( ∂ ⁡ T ) × P 1 , c ⁢ ( ∂ ⁡ T ) . The remaining orthogonality relation is checked by using the one we have just seen and trace operator tr T ∇ ⁡ div ,

( δ ⁢ u , ∇ ⁡ div ⁡ Π T ∇ ⁡ div ⁢ v ) T = ( ∇ ⁡ div ⁡ δ ⁢ u , Π T ∇ ⁡ div ⁢ v ) T + ⟨ tr T ∇ ⁡ div ⁡ δ ⁢ u , Π T ∇ ⁡ div ⁢ v ⟩ ∂ ⁡ T ∇ ⁡ div = ⟨ tr T ∇ ⁡ div ⁡ δ ⁢ u , v ⟩ ∂ ⁡ T ∇ ⁡ div = ( δ ⁢ u , ∇ ⁡ div ⁡ v ) T for all ⁢ δ ⁢ u ∈ P 0 ⁢ ( T ) 2

since tr T ∇ ⁡ div ⁡ ( P 0 ⁢ ( T ) 2 ) ⊂ P 0 ⁢ ( ∂ ⁡ T ) × P 1 , c ⁢ ( ∂ ⁡ T ) . ∎

5 Numerical Results

We consider both DPG schemes for two examples in two dimensions, one with smooth solution and one with a corner singularity. For both schemes, we use lowest-order spaces, i.e., discrete spaces U h and V h ⁢ ( T ) as in § 3.1 with p = 0 for the scheme based on the first-order system, and discrete spaces as in § 4.1 for the method based on the second-order system.

Example with Smooth Solution

In this case, we select Ω = ( 0 , 1 ) 2 , the manufactured solution

u ⁢ ( x ) := ( x 1 2 ⁢ ( x 1 - 1 ) 2 ⁢ x 2 2 ⁢ ( x 2 - 1 ) 2 , sin 2 ⁡ ( π ⁢ x 1 ) ⁢ sin 2 ⁡ ( π ⁢ x 2 ) ) T for ⁢ x = ( x 1 , x 2 ) ∈ Ω ,

and calculate 𝒇 accordingly.

Figure 1

Smooth example, quasi-uniform meshes. Left: first-order system, right: second-order system

Considering a sequence of quasi-uniform triangular meshes, we expect both schemes to deliver approximation errors of the order O ⁢ ( h ) = O ⁢ ( dim ⁡ ( U h ) - 1 / 2 ) . This is confirmed in Figure 1, with error curves for the first scheme on the left, and the second on the right. Here and below, 𝜂 stands for the corresponding residual ∥ B ⁢ u h - L ∥ V h ⁢ ( T ) ′ , and note that u 1 = u , w = - u 3 = - ∇ ⁡ div ⁡ u .

Example with Singular Solution

We select the L-shaped domain

Ω = { ( x 1 , x 2 ) ; | x 1 | + | x 2 | < 2 4 } ∖ { ( x 1 , x 2 ) ; | x 1 + 2 4 | + | x 2 | ≤ 2 4 }

with initial mesh as on the left of Figure 2, solution u := curl ⁡ v with v ⁢ ( r , φ ) := r 2 / 3 ⁢ cos ⁡ ( 2 ⁢ φ / 3 ) ( ( r , φ ) denoting polar coordinates centered at the origin so that v = 0 on the boundary edges meeting at the origin), and corresponding right-hand side function 𝒇. This solution satisfies the boundary conditions div ⁡ u = 0 on Γ and u ⋅ n = 0 along the edges meeting at the incoming corner. The non-homogeneous trace of u ⋅ n along the other edges is implemented by approximation (piecewise constant L 2 -projection in H - 1 / 2 ⁢ ( Γ ) and piecewise linear interpolation in H 1 / 2 ⁢ ( Γ ) , with subsequent extensions to elements of the discrete trace spaces).

Figure 2

Initial mesh (left) and an adaptively refined mesh with 1674 elements

Figure 3

Singular example, quasi-uniform meshes. Left: first-order system, right: second-order system

Figure 4

Singular example, adaptively refined meshes. Left: first-order system, right: second-order system

When considering the first-order system (3.1), we are dealing with solution u 1 = u and zero components u 2 , u 3 , u 4 . In the case of the second-order system (4.1), we have the solution 𝒖 and zero component 𝒘. In any case, the first (and only) non-zero component (along with its canonical trace) has a reduced regularity u ∈ H s ⁢ ( Ω ) with s < 2 / 3 , a standard Sobolev space. Therefore, for a sequence of quasi-uniform meshes, we expect a reduced convergence order O ⁢ ( h 2 / 3 ) = O ⁢ ( dim ⁡ ( U h ) - 1 / 3 ) , whereas for a scheme with appropriate mesh refinement at the incoming corner, it should be of optimal order O ⁢ ( dim ⁡ ( U h ) - 1 / 2 ) . The DPG method belongs to the family of minimum residual methods so that a posteriori error estimation is inherent; cf. [3] for details. We use the local residuals measured in V h ⁢ ( T ) ′ , T ∈ T , to steer mesh refinements, based on newest-vertex-bisection and Dörfler (or bulk) marking with parameter 3 / 4 . For both schemes, our numerical results confirm the reduced convergence order when using quasi-uniform meshes (see Figure 3 with first scheme on the left and second on the right) and illustrate optimal convergence of the adaptive variants (Figure 4). An adaptively refined mesh is shown in Figure 2.

Funding source: Fondo Nacional de Desarrollo Científico y Tecnológico

Award Identifier / Grant number: 1190009

Award Identifier / Grant number: 1210391

Funding statement: Supported by ANID-Chile through FONDECYT projects 1190009, 1210391.

References

[1] S. C. Brenner, J. Sun and L.-y. Sung, Hodge decomposition methods for a quad-curl problem on planar domains, J. Sci. Comput. 73 (2017), no. 2–3, 495–513. 10.1007/s10915-017-0449-0Search in Google Scholar

[2] F. Cakoni, D. Colton, P. Monk and J. Sun, The inverse electromagnetic scattering problem for anisotropic media, Inverse Problems 26 (2010), no. 7, Article ID 074004. 10.1088/0266-5611/26/7/074004Search in Google Scholar

[3] C. Carstensen, L. Demkowicz and J. Gopalakrishnan, A posteriori error control for DPG methods, SIAM J. Numer. Anal. 52 (2014), no. 3, 1335–1353. 10.1137/130924913Search in Google Scholar

[4] C. Carstensen, L. Demkowicz and J. Gopalakrishnan, Breaking spaces and forms for the DPG method and applications including Maxwell equations, Comput. Math. Appl. 72 (2016), no. 3, 494–522. 10.1016/j.camwa.2016.05.004Search in Google Scholar

[5] L. Chacón, A. N. Simakov and A. Zocco, Steady-state properties of driven magnetic reconnection in 2D electron magnetohydrodynamics, Phys. Rev. Lett. 99 (2007), Article ID 235001. 10.1103/PhysRevLett.99.235001Search in Google Scholar PubMed

[6] J. Chan, N. Heuer, T. Bui-Thanh and L. Demkowicz, A robust DPG method for convection-dominated diffusion problems II: Adjoint boundary conditions and mesh-dependent test norms, Comput. Math. Appl. 67 (2014), no. 4, 771–795. 10.1016/j.camwa.2013.06.010Search in Google Scholar

[7] L. Demkowicz, T. Führer, N. Heuer and X. Tian, The double adaptivity paradigm: (How to circumvent the discrete inf-sup conditions of Babuška and Brezzi), Comput. Math. Appl. 95 (2021), 41–66. 10.1016/j.camwa.2020.10.002Search in Google Scholar

[8] L. Demkowicz and J. Gopalakrishnan, A class of discontinuous Petrov–Galerkin methods. II. Optimal test functions, Numer. Methods Partial Differential Equations 27 (2011), no. 1, 70–105. 10.1002/num.20640Search in Google Scholar

[9] L. Demkowicz, J. Gopalakrishnan and A. H. Niemi, A class of discontinuous Petrov-Galerkin methods. Part III: Adaptivity, Appl. Numer. Math. 62 (2012), no. 4, 396–427. 10.1016/j.apnum.2011.09.002Search in Google Scholar

[10] L. Demkowicz and N. Heuer, Robust DPG method for convection-dominated diffusion problems, SIAM J. Numer. Anal. 51 (2013), no. 5, 2514–2537. 10.1137/120862065Search in Google Scholar

[11] A. Ern, T. Gudi, I. Smears and M. Vohralík, Equivalence of local- and global-best approximations, a simple stable local commuting projector, and optimal h ⁢ p approximation estimates in H ⁢ ( div ) , IMA J. Numer. Anal. (2021), 10.1093/imanum/draa103. 10.1093/imanum/draa103Search in Google Scholar

[12] R. Fan, Y. Liu and S. Zhang, Mixed schemes for fourth-order DIV equations, Comput. Methods Appl. Math. 19 (2019), no. 2, 341–357. 10.1515/cmam-2018-0003Search in Google Scholar

[13] W. Flügge, Stresses in Shells, Springer, Berlin, 1960. 10.1007/978-3-662-29731-5Search in Google Scholar

[14] T. Führer, A. Haberl and N. Heuer, Trace operators of the bi-Laplacian and applications, IMA J. Numer. Anal. 41 (2021), no. 2, 1031–1055. 10.1093/imanum/draa012Search in Google Scholar

[15] T. Führer and N. Heuer, Robust coupling of DPG and BEM for a singularly perturbed transmission problem, Comput. Math. Appl. 74 (2017), no. 8, 1940–1954. 10.1016/j.camwa.2016.09.016Search in Google Scholar

[16] T. Führer and N. Heuer, Fully discrete DPG methods for the Kirchhoff–Love plate bending model, Comput. Methods Appl. Mech. Engrg. 343 (2019), 550–571. 10.1016/j.cma.2018.08.041Search in Google Scholar

[17] T. Führer, N. Heuer and A. H. Niemi, An ultraweak formulation of the Kirchhoff–Love plate bending model and DPG approximation, Math. Comp. 88 (2019), no. 318, 1587–1619. 10.1090/mcom/3381Search in Google Scholar

[18] T. Führer, N. Heuer and A. H. Niemi, A DPG method for shallow shells, preprint (2021), https://arxiv.org/abs/2107.07624. 10.1007/s00211-022-01308-wSearch in Google Scholar

[19] T. Führer, N. Heuer and F.-J. Sayas, An ultraweak formulation of the Reissner-Mindlin plate bending model and DPG approximation, Numer. Math. 145 (2020), no. 2, 313–344. 10.1007/s00211-020-01116-0Search in Google Scholar

[20] V. Girault and P.-A. Raviart, Finite Element Methods for Navier–Stokes Equations, Springer Ser. Comput. Math. 5, Springer, Berlin, 1986. 10.1007/978-3-642-61623-5Search in Google Scholar

[21] J. Gopalakrishnan and W. Qiu, An analysis of the practical DPG method, Math. Comp. 83 (2014), no. 286, 537–552. 10.1090/S0025-5718-2013-02721-4Search in Google Scholar

[22] N. Heuer and M. Karkulik, A robust DPG method for singularly perturbed reaction-diffusion problems, SIAM J. Numer. Anal. 55 (2017), no. 3, 1218–1242. 10.1137/15M1041304Search in Google Scholar

[23] Q. Hong, J. Hu, S. Shu and J. Xu, A discontinuous Galerkin method for the fourth-order curl problem, J. Comput. Math. 30 (2012), no. 6, 565–578. 10.4208/jcm.1206-m3572Search in Google Scholar

[24] G. Kirchhoff, Über das Gleichgewicht und die Bewegung einer elastischen Scheibe, J. Reine Angew. Math. 40 (1850), 51–88. Search in Google Scholar

[25] A. E. H. Love, The small free vibrations and deformation of a thin elastic shell, Phil. Trans. R. Soc. Lond. A 179 (1888), 491–546. 10.1098/rsta.1888.0016Search in Google Scholar

[26] R. D. Mindlin, Micro-structure in linear elasticity, Arch. Rational Mech. Anal., 16 (1964), 51–78. 10.1007/BF00248490Search in Google Scholar

[27] P. Monk and J. Sun, Finite element methods for Maxwell’s transmission eigenvalues, SIAM J. Sci. Comput. 34 (2012), no. 3, B247–B264. 10.1137/110839990Search in Google Scholar

[28] S. Nagaraj, S. Petrides and L. F. Demkowicz, Construction of DPG Fortin operators for second order problems, Comput. Math. Appl. 74 (2017), no. 8, 1964–1980. 10.1016/j.camwa.2017.05.030Search in Google Scholar

[29] E. Priest and T. Forbes, Magnetic Reconnection, Cambridge University, Cambridge, 2000. 10.1017/CBO9780511525087Search in Google Scholar

[30] N. V. Roberts, L. Demkowicz and R. Moser, A discontinuous Petrov–Galerkin methodology for adaptive solutions to the incompressible Navier–Stokes equations, J. Comput. Phys. 301 (2015), 456–483. 10.1016/j.jcp.2015.07.014Search in Google Scholar

[31] J. Sun, A mixed FEM for the quad-curl eigenvalue problem, Numer. Math. 132 (2016), no. 1, 185–200. 10.1007/s00211-015-0708-7Search in Google Scholar

Received: 2021-12-30

Accepted: 2022-03-10

Published Online: 2022-05-24

Published in Print: 2022-07-01

DPG Methods for a Fourth-Order div Problem

Abstract

1 Introduction

2 Trace Spaces and Norms

2.1 Canonical Trace Operators

Proposition 2 ([4, Lemma 2.1])

2.2 A grad-div Trace Operator

Proof

Proof

Proof

Proof

3 First-Order Variational Formulation and DPG Method

3.1 Fully Discrete Method

3.2 Proofs for the First-Order Formulation

3.2.1 Stability of the Adjoint Problem

Proof

Proof

3.2.2 Proofs of Theorems 7, 8

4 Second-Order Variational Formulation and DPG Method

4.1 Fully Discrete Method

4.2 Proofs for the Second-Order Formulation

4.2.1 Stability of the Adjoint Problem

Proof

Proof

4.2.2 Proofs of Theorems 12, 13

4.2.3 Fortin Operator

Proof

Proof

5 Numerical Results

Example with Smooth Solution

Example with Singular Solution

References

Journal and Issue

Articles in the same Issue