Uniform Preconditioners of Linear Complexity for Problems of Negative Order

Rob Stevenson; Raymond van Venetië

doi:10.1515/cmam-2020-0052

Open Access Published by De Gruyter December 2, 2020

Uniform Preconditioners of Linear Complexity for Problems of Negative Order

Rob Stevenson and Raymond van Venetië

From the journal Computational Methods in Applied Mathematics

https://doi.org/10.1515/cmam-2020-0052

Abstract

We propose a multi-level type operator that can be used in the framework of operator (or Caldéron) preconditioning to construct uniform preconditioners for negative order operators discretized by piecewise polynomials on a family of possibly locally refined partitions. The cost of applying this multi-level operator scales linearly in the number of mesh cells. Therefore, it provides a uniform preconditioner that can be applied in linear complexity when used within the preconditioning framework from our earlier work [Uniform preconditioners for problems of negative order, Math. Comp. 89 (2020), 645–674].

Keywords: Operator Preconditioning; Uniform Preconditioners; Operators of Negative Order; Finite and Boundary Elements

MSC 2010: 65F08; 65N38; 65N30; 45Exx

1 Introduction

In this work, we construct a multi-level type preconditioner for operators of negative orders - 2 ⁢ s ∈ [ - 2 , 0 ] that can be applied in linear time and yields uniformly bounded condition numbers. The preconditioner will be constructed using the framework of “operator preconditioning” studied earlier in e.g. [12, 3, 6, 10]. The role of the “opposite order operator” will be fulfilled by a multi-level type operator, based on the work of Wu and Zheng in [14].

For some 𝑑-dimensional domain (or manifold) Ω, a measurable, closed, possibly empty γ ⊂ ∂ ⁡ Ω and an s ∈ [ 0 , 1 ] , we consider the Sobolev spaces

W := [ L 2 ⁢ ( Ω ) , H 0 , γ 1 ⁢ ( Ω ) ] s , 2 , V := W ′ ,

with H 0 , γ 1 ⁢ ( Ω ) being the closure in H 1 ⁢ ( Ω ) of the smooth functions on Ω that vanish at 𝛾. Let ( V T ) T ∈ T ⊂ V be a family of piecewise or continuous piecewise polynomials of some fixed degree w.r.t. uniformly shape regular, possibly locally refined partitions. With, for T ∈ T , A T : V T → V T ′ being some boundedly invertible linear operator, we are interested in constructing a preconditioner G T : V T ′ → V T such that the preconditioned operator G T ⁢ A T : V T → V T is uniformly boundedly invertible, and an application of G T can be evaluated in O ⁢ ( dim ⁡ V T ) arithmetic operations.

In order to create such a preconditioner, we will use the framework described in our earlier work [10]. Given V T , we constructed an auxiliary space W T ⊂ W with dim ⁡ W T = dim ⁡ V T such that, for D T defined by ( D T v ) ( w ) := ⟨ v , w ⟩ L 2 ⁢ ( Ω ) ( v ∈ V T , w ∈ W T ) and some suitable “opposite order” operator B T W : W T → W T ′ , a preconditioner G T of the form G T := D T - 1 ⁢ B T W ⁢ ( D T ′ ) - 1 is found. The space W T is equipped with a basis that, modulo a scaling, is biorthogonal to the canonical basis for V T so that the representation of D T is an invertible diagonal matrix.

With S T , 0 0 , 1 ⊂ W being the space of continuous piecewise linears w.r.t. 𝒯, zero on 𝛾, the above preconditioning approach hinges on the availability of a uniformly boundedly invertible operator B T S : S T , 0 0 , 1 → ( S T , 0 0 , 1 ) ′ , which is generally the most demanding requirement. For example, if s = 1 2 and γ = ∅ , a viable option is to take B T S as the discretized hypersingular operator. While this induces a uniform preconditioner, the application of B T S cannot be evaluated in linear complexity.

In this work, we construct a suitable multi-level type operator B T S that can be applied in linear complexity. For this construction, we require 𝕋 to be a family of conforming partitions created by newest vertex bisection [7, 13]. In the aforementioned setting of having an arbitrary s ∈ [ 0 , 1 ] , this multi-level operator B T S induces a uniform preconditioner G T , i.e., G T ⁢ A T is uniformly well-conditioned, where the cost of applying G T scales linearly in dim ⁡ V T . We also show that the preconditioner extends to the more general manifold case, where Ω is a 𝑑-dimensional (piecewise) smooth Lipschitz manifold, and the trial space V T is the parametric lift of a space of piecewise or continuous piecewise polynomials.

Finally, we remark that common multi-level preconditioners based on overlapping subspace decompositions are known not to work well for operators of negative order. A solution is provided by resorting to direct sum multi-level subspace decompositions. Examples are given by wavelet preconditioners, or closely related, the preconditioners from [2], for the latter assuming quasi-uniform partitions.

For - s = - 1 2 , an optimal multi-level preconditioner based on a non-overlapping subspace decomposition for operators defined on the boundary of a 2- or 3-dimensional Lipschitz polyhedron was recently introduced in [4].

1.1 Outline

In Section 2, we summarize the (operator) preconditioning framework from [10]. In Section 3, we provide the multi-level type operator that can be used as the “opposite order” operator inside the preconditioner framework. In Section 4, we comment on how to generalize the results to the case of piecewise smooth manifolds. In Section 5, we conclude with numerical results.

1.2 Notation

In this work, by λ ≲ μ , we mean that 𝜆 can be bounded by a multiple of 𝜇, independently of parameters which 𝜆 and 𝜇 may depend on, with the sole exception of the space dimension 𝑑, or in the manifold case, on the parametrization of the manifold that is used to define the finite element spaces on it. Obviously, λ ≳ μ is defined as μ ≲ λ , and λ ≂ μ as λ ≲ μ and λ ≳ μ .

For normed linear spaces Y and Z , in this paper, for convenience over ℝ, L ⁢ ( Y , Z ) will denote the space of bounded linear mappings Y → Z endowed with the operator norm ∥ ⋅ ∥ L ⁢ ( Y , Z ) . The subset of invertible operators in L ⁢ ( Y , Z ) with inverses in L ⁢ ( Z , Y ) will be denoted as L ⁢ is ⁢ ( Y , Z ) . The condition number of a C ∈ L ⁢ is ⁢ ( Y , Z ) is defined as κ Y , Z ⁢ ( C ) := ∥ C ∥ L ⁢ ( Y , Z ) ⁢ ∥ C - 1 ∥ L ⁢ ( Z , Y ) .

For Y a reflexive Banach space and C ∈ L ⁢ ( Y , Y ′ ) being coercive, i.e.,

inf 0 ≠ y ∈ Y ⁡ ( C ⁢ y ) ⁢ ( y ) ∥ y ∥ Y 2 > 0 ,

both 𝐶 and ℜ ⁡ ( C ) := 1 2 ⁢ ( C + C ′ ) are in L ⁢ is ⁢ ( Y , Y ′ ) with

∥ ℜ ( C ) ∥ L ⁢ ( Y , Y ′ ) ≤ ∥ C ∥ L ⁢ ( Y , Y ′ ) , ∥ C - 1 ∥ L ⁢ ( Y ′ , Y ) ≤ ∥ ℜ ( C ) - 1 ∥ L ⁢ ( Y ′ , Y ) = ( inf 0 ≠ y ∈ Y ( C ⁢ y ) ⁢ ( y ) ∥ y ∥ Y 2 ) - 1 .

The set of coercive C ∈ L ⁢ is ⁢ ( Y , Y ′ ) is denoted as L ⁢ is c ⁢ ( Y , Y ′ ) . If C ∈ L ⁢ is c ⁢ ( Y , Y ′ ) , then C - 1 ∈ L ⁢ is c ⁢ ( Y ′ , Y ) and ∥ ℜ ( C - 1 ) - 1 ∥ L ⁢ ( Y , Y ′ ) ≤ ∥ C ∥ L ⁢ ( Y , Y ′ ) 2 ∥ ℜ ( C ) - 1 ∥ L ⁢ ( Y ′ , Y ) .

Given a family of operators C i ∈ L ⁢ is ⁢ ( Y i , Z i ) ( L ⁢ is c ⁢ ( Y i , Z i ) ), we will write C i ∈ L ⁢ is ⁢ ( Y i , Z i ) ( L ⁢ is c ⁢ ( Y i , Z i ) ) uniformly in 𝑖, or simply “uniform”, when

sup i max ( ∥ C i ∥ L ⁢ ( Y i , Z i ) , ∥ C i - 1 ∥ L ⁢ ( Z i , Y i ) ) < ∞ or sup i max ( ∥ C i ∥ L ⁢ ( Y i , Z i ) , ∥ ℜ ( C i ) - 1 ∥ L ⁢ ( Z i , Y i ) ) < ∞ .

2 Preconditioning

Let ( T ) T ∈ T be a family of conforming partitions of a domain Ω ⊂ R d into (open) uniformly shape regular𝑑-simplices, where we assume that 𝛾 is the (possibly empty) union of ( d - 1 ) -faces of T ∈ T . For d ≥ 2 , such partitions automatically satisfy a uniform 𝐾-mesh property, and for d = 1 , we impose this as an additional condition. The discussion of the manifold case is postponed to Section 4.

Recalling that V T ⊂ V is a family of piecewise or continuous piecewise polynomials of some fixed degree w.r.t. 𝒯, let A T ∈ L ⁢ is ⁢ ( V T , V T ′ ) uniformly in T ∈ T . A common setting is that ( A T ⁢ v ) ⁢ ( v ~ ) := ( A ⁢ v ) ⁢ ( v ~ ) ( v , v ~ ∈ V T ) for some A ∈ L ⁢ is c ⁢ ( V , V ′ ) . We are interested in finding optimal preconditioners G T for A T , i.e., G T ∈ L ⁢ is ⁢ ( V T ′ , V T ) uniformly in T ∈ T , whose application moreover requires O ⁢ ( dim ⁡ V T ) arithmetic operations.

Recall the space

S T , 0 0 , 1 := { u ∈ H 0 , γ 1 ( Ω ) : u | T ∈ P 1 ( T ∈ T ) } ⊂ W

(thus equipped with ∥ ⋅ ∥ W ). In [10], using operator preconditioning, we reduced the issue of constructing such preconditioners G T to the issue of constructing B T S ∈ L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) uniformly. In the next section, we summarize this reduction.

2.1 Construction of Optimal Preconditioners

For the moment, consider the lowest order case of V T being either the space of piecewise constants or continuous piecewise linears. In [10], a space W T ⊂ W was constructed with dim ⁡ W T = dim ⁡ V T and

(2.1) inf T ∈ T ⁡ inf 0 ≠ v ∈ V T ⁡ sup 0 ≠ w ∈ W T ⁡ ⟨ v , w ⟩ L 2 ⁢ ( Ω ) ∥ v ∥ V ⁢ ∥ w ∥ W > 0 .

Moreover, W T ⊂ W was equipped with a locally supported basis Ψ T that, modulo a scaling, is L 2 ⁢ ( Ω ) -biorthogonal to the canonical basis Ξ T of V T .

As a consequence of (2.1), D T defined by ( D T ⁢ v ) ⁢ ( w ) := ⟨ v , w ⟩ L 2 ⁢ ( Ω ) ( v ∈ V T , w ∈ W T ) is in L ⁢ is ⁢ ( V T , W T ′ ) uniformly. We infer that, once we have constructed B T W ∈ L ⁢ is ⁢ ( W T , W T ′ ) uniformly, then, by taking

(2.2) G T := D T - 1 ⁢ B T W ⁢ ( D T ′ ) - 1 ,

we have G T ∈ L ⁢ is ⁢ ( V T ′ , V T ) uniformly. Biorthogonality, modulo a scaling, of the bases Ψ T and Ξ T implies that the matrix representation of D T is diagonal so that D T - 1 and its adjoint can be applied in linear complexity.

The aforementioned space W T is a subspace of S T , 0 0 , 1 ⊕ B T ⊂ W , where B T is a “bubble space” with dim ⁡ B T = O ⁢ ( # ⁢ T ) , such that the projector I T on S T , 0 0 , 1 ⊕ B T defined by ran ⁡ I T = S T , 0 0 , 1 and ran ⁡ ( Id - I T ) = B T is “local” and uniformly bounded, and the canonical basis Θ T of “bubbles” for B T is, when normalized in ∥ ⋅ ∥ W , a uniformly Riesz basis for B T . Because of the latter, B T B defined by

( B T B ⁢ c ⊤ ⁢ Θ T ) ⁢ ( d ⊤ ⁢ Θ T ) := β ⁢ ( Δ T ⁢ c ) ⊤ ⁢ d

for some diagonal Δ T ≂ diag ⁡ ( ⟨ Θ T , Θ T ⟩ W ) and constant β > 0 is in L ⁢ is c ⁢ ( B T , B T ′ ) uniformly.

Given some “opposite order” operator B T S ∈ L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) , by taking

(2.3) B T W := I T ′ ⁢ B T S ⁢ I T + ( Id - I T ) ′ ⁢ B T B ⁢ ( Id - I T ) ,

it holds that B T W ∈ L ⁢ is c ⁢ ( W T , W T ′ ) uniformly [11, Proposition 5.1], which makes G T a uniform preconditioner.

2.2 Implementation of G T

Recalling the aforementioned bases Ξ T , Ψ T and Θ T for V T , W T and B T , respectively, equipping S T , 0 0 , 1 with the nodal basis Φ T and equipping V T ′ , W T ′ , B T ′ and ( S T , 0 0 , 1 ) ′ with the dual bases Ξ T ′ , Ψ T ′ , Θ T ′ and Φ T ′ , respectively, the representation of A T ∈ L ⁢ ( V T , V T ′ ) is the stiffness matrix A T := ( A T ⁢ Ξ T ) ⁢ ( Ξ T ) := [ ( A T ⁢ η ) ⁢ ( ξ ) ] ( ξ , η ) ∈ Ξ T , and the representation of G T ∈ L ⁢ ( V T ′ , V T ) is the matrix G T := ( G ⁢ Ξ T ′ ) ⁢ ( Ξ T ′ ) . It is given by

(2.4) G T = D T - 1 ⁢ ( p T ⊤ ⁢ B T S ⁢ p T + q T ⊤ ⁢ B T B ⁢ q T ) ⁢ D T - ⊤ ,

where both

D T := ( D T ⁢ Ξ T ) ⁢ ( Ψ T ) , B T B := ( B T B ⁢ Θ T ) ⁢ ( Θ T )

are diagonal, both

p T := ( I T ⁢ Ψ T ) ⁢ ( Φ T ′ ) , q T := ( ( Id - I T ) ⁢ Ψ T ) ⁢ ( Θ T ′ )

are uniformly sparse and

(2.5) B T S := ( B T S ⁢ Φ T ) ⁢ ( Φ T ) .

Note that the cost of the application of G T scales linearly in # ⁢ T as soon as this holds true for the application of B T S .

The above preconditioning approach is summarized in the following theorem.

Theorem 2.1

Theorem 2.1 ([10, Section 3])

Given a family B T S ∈ L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) uniformly in T ∈ T , then for B T W as described in (2.3), the operator G T from (2.2) is a uniform preconditioner. Furthermore, if the matrix representation B T S , cf. (2.5), can be applied in O ⁢ ( # ⁢ T ) operations, then the matrix representation of the preconditioner G T , cf. (2.4), can be applied in O ⁢ ( # ⁢ T ) operations.

Because B T W in (2.3) is given as the sum of two operators that “act” on different subspaces of W T , the condition number of the preconditioned system depends on the relative scaling of both these operators which can be steered by selecting the parameter 𝛽. A suitable 𝛽 will be selected experimentally.

Alternatively, [11, Proposition 5.1] shows that a value of 𝛽 is reasonable if it is chosen such that the interval bounded by the coercivity and boundedness constants of B T S is included in that interval corresponding to B T B or vice versa. Also these coercivity and boundedness constants can be approximated experimentally or by making some theoretical estimates.

Constructions of Ψ T , Θ T and Δ T , and resulting explicit formulas for matrices D T , B T B , p T , q T are derived in [10]. For ease of reading, we recall these formulas below for the case that V T is the space of piecewise constants. For the continuous piecewise linear case, we refer to [10, Section 4.2].

2.2.1 Piecewise Constant Trial Space V T

For T ∈ T , we define N T as the set of vertices of 𝒯, and N T 0 as the set of vertices of 𝒯 that are not on 𝛾. For ν ∈ N T , we set its valence d T , ν := # ⁢ { T ∈ T : ν ∈ T ¯ } . For T ∈ T , and with N T denoting the set of its vertices, we set N T , T 0 := N T 0 ∩ N T .

If one considers V T as the space of discontinuous piecewise constants, i.e.

V T = S T - 1 , 0 := { u ∈ L 2 ( Ω ) : u | T ∈ P 0 ( T ∈ T ) } ⊂ V ,

equipped with the canonical basis Ξ T := { 1 T : T ∈ T } , then we find, for arbitrary constant β > 0 ,

D T = diag ⁡ { | T | : T ∈ T } , ( p T ) ν ⁢ T = { d T , ν - 1 if ⁢ ν ∈ N T , T 0 , 0 if ⁢ ν ∉ N T , T 0 ,

B T B = β ⁢ D T 1 - 2 ⁢ s d , ( q T ) T ′ ⁢ T = δ T ′ ⁢ T - 1 d + 1 ⁢ ∑ ν ∈ N T , T 0 ∩ N T , T ′ 0 d T , ν - 1 .

2.3 Higher Order Case

For higher order discontinuous or continuous finite element spaces V T , suitable preconditioners G T can be built either from the current preconditioner G T for the lowest order case by application of a subspace correction method (most conveniently in the discontinuous case where, on each element, the space of polynomials of some fixed degree is split into the space of constants and its orthogonal complement), or by expanding W T by enlarging the bubble space B T . While referring to [10] for details, we recall that, with either option, the construction of an optimal preconditioner G T that can be applied in linear complexity hinges on the availability of an operator B T S ∈ L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) uniformly in T ∈ T , that can be applied in linear complexity.

3 An Operator B T S of Multi-level Type

In this section, we will introduce an operator B T S ∈ L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) of multi-level type. The operator B T S is based on a stable multi-level decomposition of S T , 0 0 , 1 given by Wu and Zheng [14]. Usually, such a stable multi-level decomposition is used as a theoretical tool for proving optimality of an additive (or multiplicative) Schwarz type preconditioner for an operator in L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) . In this work, however, we are going to use their results for the construction of the operator B T S ∈ L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) for which it is crucial that its application can be implemented in linear complexity.

3.1 Definition and Analysis of B T S

For d ≥ 2 , let 𝕋 be the family of all conforming partitions of Ω into 𝑑-simplices that can be created by newest vertex bisection starting from some given conforming initial partition T ⊥ that satisfies a matching condition [9].

With T := ⋃ T ∈ T { T : T ∈ T } and N := ⋃ T ∈ T N T , for T ∈ T , let gen ⁡ ( T ) be the number of bisections needed to create 𝑇 from its ancestor T ′ ∈ T ⊥ , and for ν ∈ N , let gen ⁡ ( ν ) := min ⁡ { gen ⁡ ( T ) : T ∈ T , ν ∈ N T } . Notice that | T | ≂ 2 - gen ⁡ ( T ) . For T ∈ T , let Q T denote the L 2 ⁢ ( T ) -orthogonal projector onto P 1 ⁢ ( T ) .

The case d = 1 can be included by letting 𝕋 be the family of a partitions of Ω that can be constructed by bisections from T ⊥ = { Ω } such that the generations of any two neighboring subintervals in any T ∈ T differ by not more than one.

For T ∈ T , set L = L ⁢ ( T ) := max T ∈ T ⁡ gen ⁡ ( T ) , and define

T ⊥ = T 0 ≺ T 1 ≺ ⋯ ≺ T L = T ⊂ T

by constructing T j - 1 from T j by removing all ν ∈ N j := N T j from the latter for which gen ⁡ ( ν ) = j . For ν ∈ N j 0 := N T j 0 , we define ω j ⁢ ( ν ) = ⋃ { T ∈ T j : ν ∈ N T } .

With this hierarchy of partitions, we define an averaging quasi-interpolator Π j ∈ L ⁢ ( S T , 0 0 , 1 , S T j , 0 0 , 1 ) by

(3.1) ( Π j u ) ( ν ) := ∑ { T ∈ T j : ν ∈ N T } | T | ⁢ ( Q T ⁢ u ) ⁢ ( ν ) ∑ { T ∈ T j : ν ∈ N T } | T | ( u ∈ S T , 0 0 , 1 , ν ∈ N j 0 ) .

Since S T j , 0 0 , 1 is a space of continuous piecewise linears, it indeed suffices to define Π j ⁢ u at the vertices N j 0 . Recall that S T , 0 0 , 1 ⊂ W := [ L 2 ⁢ ( Ω ) , H 0 , γ 1 ⁢ ( Ω ) ] s , 2 for some s ∈ [ 0 , 1 ] . The next theorem shows that Π j induces a stable multi-level decomposition of S T , 0 0 , 1 .

Theorem 3.1

Theorem 3.1 ([14, Lemma 3.7])

For the averaging quasi-interpolator Π j from (3.1), and Π - 1 := 0 , it holds that

∥ u ∥ W 2 ≂ ∑ j = 0 L 4 j ⁢ s / d ∥ ( Π j - Π j - 1 ) u ∥ L 2 ⁢ ( Ω ) 2 ( u ∈ S T , 0 0 , 1 ) .

Proof

In [14], the inequality “ ≳ ” was proven for the case s = 1 , d ∈ { 2 , 3 } and γ = ∂ ⁡ Ω . The arguments, however, immediately extend to s ∈ [ 0 , 1 ] , d ≥ 1 and γ ⊊ ∂ ⁡ Ω .

The proof of the other inequality “ ≲ ” follows from well-known arguments. For some t ∈ ( 1 , 3 2 ) , let H r ⁢ t := [ L 2 ⁢ ( Ω ) , H 0 , γ 1 ⁢ ( Ω ) ∩ H t ⁢ ( Ω ) ] r , 2 for r ∈ [ 0 , 1 ] . Then H s ≃ W by the reiteration theorem, and for r ∈ [ 0 , 1 ] , ∥ ⋅ ∥ H r ⁢ t ≲ 2 j ⁢ r ⁢ t / d ⁢ ∥ ⋅ ∥ L 2 ⁢ ( Ω ) on S T j , 0 0 , 1 .

Let u ∈ S T , 0 0 , 1 be written as ∑ j = 0 L u j with u j ∈ S T j , 0 0 , 1 . Then, for ε ∈ ( 0 , s ) , ε ≤ t - s , we have

∥ u ∥ H s ⁢ ( Ω ) 2 ≲ ∑ j = 0 L ∑ i = j L ∥ u j ∥ H s + ε ⁢ ( Ω ) ⁢ ∥ u i ∥ H s - ε ⁢ ( Ω ) ≲ ∑ j = 0 L ∑ i = j L 2 j ⁢ ( s + ε ) / d ⁢ 2 i ⁢ ( s - ε ) / d ⁢ ∥ u j ∥ L 2 ⁢ ( Ω ) ⁢ ∥ u i ∥ L 2 ⁢ ( Ω ) ≲ ∑ j = 0 L 4 j ⁢ s / d ⁢ ∥ u j ∥ L 2 ⁢ ( Ω ) 2 . ∎

The relevance of the multi-level decomposition from Theorem 3.1 by Wu and Zheng lies in the fact that ( Π j ⁢ u ) ⁢ ( ν ) can only differ from ( Π j - 1 ⁢ u ) ⁢ ( ν ) in any ν ∈ N j 0 ∖ N j - 1 0 as well as in only two ^[1] of its neighbors in N j - 1 0 (the endpoints of the edge on which 𝜈 was inserted).

Proposition 3.2

Proposition 3.2 ([14, Lemma 3.1])

With, for j ≥ 1 , M j 0 := { ν ∈ N j - 1 0 : ω j ⁢ ( ν ) = ω j - 1 ⁢ ( ν ) } , it holds that, for ν ∈ M j 0 , ( ( Π j - Π j - 1 ) ⁢ u ) ⁢ ( ν ) = 0 ; see Figure 1.

$Figure 1 For d = 3 d=3 , a tetrahedron T ∈ T j - 1 T\in\mathcal{T}_{j-1} and its bisection. The dots indicate all vertices in N j 0 ∖ M j 0 N_{j}^{0}\setminus M^{0}_{j} .$

Figure 1

For d = 3 , a tetrahedron T ∈ T j - 1 and its bisection. The dots indicate all vertices in N j 0 ∖ M j 0 .

Remark 3.3

The proof from [14] given for d ∈ { 2 , 3 } generalizes to d ≥ 1 . Indeed, the arguments that are used are based on the fact that the basis for S 1 ⁢ ( T ) that is dual to the nodal basis takes equal values in all but one nodal point. This is a consequence of the fact that the mass matrix of the nodal basis for S 1 ⁢ ( T ) , and so its inverse, is invariant under permutations of the barycentric coordinates, which holds true in any dimension.

As a consequence of Proposition 3.2, we have

∥ ( Π j - Π j - 1 ) ⁢ u ∥ L 2 ⁢ ( Ω ) 2 ≂ 2 - j ⁢ ∑ ν ∈ N j 0 ∖ M j 0 | ( ( Π j - Π j - 1 ) ⁢ u ) ⁢ ( ν ) | 2 .

From Theorem 3.1, we conclude that B T S = ( B T S ) ′ ∈ L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) defined by

(3.2) ( B T S ⁢ u ) ⁢ ( v ) := ∑ j = 0 L 2 j ⁢ ( 2 ⁢ s / d - 1 ) ⁢ ∑ ν ∈ N j 0 ∖ M j 0 ( ( Π j - Π j - 1 ) ⁢ u ) ⁢ ( ν ) ⁢ ( ( Π j - Π j - 1 ) ⁢ v ) ⁢ ( ν )

is uniform, i.e.

sup T ∈ T ⁡ max ⁡ ( ∥ B T S ∥ L ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) , ∥ ( B T S ) - 1 ∥ L ⁢ ( ( S T , 0 0 , 1 ) ′ , S T , 0 0 , 1 ) ) < ∞ .

3.2 Implementation of B T S

Since the operator Π j is a weighted local L 2 ⁢ ( Ω ) projection, it allows for a natural implementation by considering S T - 1 , 1 , the space of discontinuous piecewise linears w.r.t. 𝒯. Recall the nodal basis Φ T for S T , 0 0 , 1 , and equip S T - 1 , 1 with the element-wise nodal basis.

Denote E T for the representation of the embedding S T , 0 0 , 1 into S T - 1 , 1 . For 0 ≤ j ≤ L , let R j be the representation of the L 2 ⁢ ( Ω ) -orthogonal projector of S T - 1 , 1 onto S T j - 1 , 1 , and let R - 1 := 0 . For 0 ≤ j ≤ L , let H j be the representation of the averaging operator H j : S T j - 1 , 1 → S T j , 0 0 , 1 defined by

(3.3) ( H j u ) ( ν ) = ∑ { T ∈ T j : ν ∈ N T } | T | ⁢ u | T ⁢ ( ν ) ∑ { T ∈ T j : ν ∈ N T } | T | ( ν ∈ N T 0 ) ,

and let H - 1 := 0 . For 1 ≤ j ≤ L , let P j be the representation of the embedding S T j - 1 , 0 0 , 1 → S T j , 0 0 , 1 (often called prolongation), and let P 0 := 0 .

Then the representation B T S of B T S from (3.2) is given by

B T S = E T ⊤ ⁢ ( ∑ j = 0 L ( H j ⁢ R j - P j ⁢ H j - 1 ⁢ R j - 1 ) ⊤ ⁢ 2 j ⁢ ( 2 ⁢ s / d - 1 ) ⁢ ( H j ⁢ R j - P j ⁢ H j - 1 ⁢ R j - 1 ) ) ⁢ E T .

Applying E T amounts to duplicating values at any internal node with a number equal to the valence of that node.

By representing 𝒯 as the leaves of a binary tree with roots being the simplices of T ⊥ , computing for x → ∈ ran ⁡ E T the sequence ( R j ⁢ x → ) 0 ≤ j ≤ L amounts to computing, while traversing from the leaves to the root, for any parent and both its children the orthogonal projection of a piecewise linear function on the children to the space of linears on the parent. For d = 2 , the matrix representation of the latter projection is given in Figure 2.

Figure 2

Numbering of the vertices of the parent and that of both children for d = 2 , and the resulting matrix representation of the orthogonal projection of the space of piecewise linears on the children to the space of linears on the parent.

Proposition 3.4

The application of B T S can be computed in O ⁢ ( # ⁢ T ) operations.

Proof

Because the number of nodes in a binary tree is less than 2 times the number of its leaves, for x → ∈ R dim ⁡ S T , 0 0 , 1 , the computation of the sequence ( R j ⁢ E T ⁢ x → ) 0 ≤ j ≤ L takes O ⁢ ( # ⁢ T ) operations. From Proposition 3.2, recall that any vector in ran ⁡ H j ⁢ R j - P j ⁢ H j - 1 ⁢ R j - 1 vanishes at M j 0 so that the number of its non-zero entries is bounded by # ⁢ ( N j 0 ∖ M j 0 ) ≤ 3 ⁢ # ⁢ ( N j 0 ∖ N j - 1 0 ) . Knowing already R j ⁢ E T ⁢ x → and R j - 1 ⁢ E T ⁢ x → , computing any non-zero entry of ( H j ⁢ R j - P j ⁢ H j - 1 ⁢ R j - 1 ) ⁢ E T ⁢ x → requires O ⁢ ( 1 ) operations. ∎

We conclude that the operator B T S , with above matrix representation B T S , satisfies the requirements of Theorem 2.1.

4 Manifold Case

Let Γ be a compact 𝑑-dimensional Lipschitz, piecewise smooth manifold in R d ′ for some d ′ ≥ d with or without boundary ∂ ⁡ Γ . For some closed measurable γ ⊂ ∂ ⁡ Γ and s ∈ [ 0 , 1 ] , let

W := [ L 2 ⁢ ( Γ ) , H 0 , γ 1 ⁢ ( Γ ) ] s , 2 , V := W ′ .

We assume that Γ is given as the essentially disjoint union of ⋃ i = 1 p χ i ⁢ ( Ω i ) ¯ , with, for 1 ≤ i ≤ p , χ i : R d → R d ′ being some smooth regular parametrization and Ω i ⊂ R d an open polytope. Without loss of generality, assuming that, for i ≠ j , Ω ¯ i ∩ Ω ¯ j = ∅ , we define

χ : Ω := ⋃ i = 1 p Ω i → ⋃ i = 1 p χ i ( Ω i ) by χ | Ω i = χ i .

Let 𝕋 be a family of conforming partitions 𝒯 of Γ into “panels” such that, for 1 ≤ i ≤ p , χ - 1 ⁢ ( T ) ∩ Ω i is a uniformly shape regular conforming partition of Ω i into 𝑑-simplices (that, for d = 1 , satisfies a uniform 𝐾-mesh property). We assume that 𝛾 is a (possibly empty) union of “faces” of T ∈ T (i.e., sets of type χ i ⁢ ( e ) , where 𝑒 is a ( d - 1 ) -dimensional face of χ i - 1 ⁢ ( T ) ).

We set

V T := { u ∈ L 2 ( Γ ) : u ∘ χ | χ - 1 ⁢ ( T ) ∈ P 0 ( T ∈ T ) } ⊂ V

V T := { u ∈ C ( Γ ) : u ∘ χ | χ - 1 ⁢ ( T ) ∈ P 1 ( T ∈ T ) } ⊂ V ,

equipped with canonical basis Ξ T , and, for the construction of a preconditioner,

S T , 0 0 , 1 := { u ∈ H 0 , γ 1 ( Γ ) : u ∘ χ | χ - 1 ⁢ ( T ) ∈ P 1 ( T ∈ T ) } ⊂ W ,

equipped with canonical basis Φ T .

As in the domain case, a space W T ⊂ W can be constructed with

dim ⁡ W T = dim ⁡ V T and inf T ∈ T ⁡ inf 0 ≠ v ∈ V T ⁡ sup 0 ≠ w ∈ W T ⁡ ⟨ v , w ⟩ L 2 ⁢ ( Γ ) ∥ v ∥ V ⁢ ∥ w ∥ W > 0 ,

which can be equipped with a locally supported basis Ψ T that, modulo a scaling, is L 2 ⁢ ( Γ ) -biorthogonal to Ξ T . Now assuming that a family of B T S ∈ L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) uniformly is available, the construction of an optimal preconditioner G T follows exactly the same lines as outlined in Section 2 for the domain case.

For the case that Γ is not piecewise polytopal, a hidden problem is, however, that above construction of Ψ T requires exact integration of lifted polynomials over the manifold. To circumvent this problem, in [10], we have relaxed the condition of L 2 ⁢ ( Γ ) -biorthogonality of Ξ T and Ψ T to biorthogonality w.r.t. to a mesh-dependent scalar product obtained from the L 2 ⁢ ( Γ ) -scalar product by replacing the Jacobian on the pull back of each panel by its mean. It was shown that the resulting preconditioner is still optimal and that the expression for its matrix representation (for the moment without the representation of B T S ), that was recalled in Section 2.2.1 for the piecewise constant case, applies verbatim by only reading | T | as the volume of the panel.^[2]

It remains to discuss the construction of an operator B T S of multi-level type, where it is now assumed that 𝕋 is a family corresponding to newest vertex bisection. An exact copy of the construction of B T S given in the domain case would require the application of the panel-wise L 2 ⁢ ( T ) -orthogonal projector Q T , cf. (3.1), which generally poses a quadrature problem. Reconsidering the domain case, the proof of [14, Lemma 3.7] (which provides the proof of the inequality “ ≳ ” in our Theorem 3.1) builds on the fact that, for T 0 ≺ T 1 ≺ ⋯ being a sequence of uniformly refined partitions, the decomposition S T L , 0 0 , 1 = ∑ j = 0 L S T j , 0 0 , 1 ∩ ( S T j - 1 , 0 0 , 1 ) ⟂ L 2 ⁢ ( Ω ) , where S T - 1 , 0 0 , 1 := { 0 } , is stable, uniformly in 𝐿, w.r.t. the norm on W . This stability holds also true when the orthogonal complements are taken w.r.t. a weighted L 2 ⁢ ( Ω ) -scalar product for any weight 𝑤 with w , 1 / w ∈ L ∞ ⁢ ( Ω ) .

This has the consequence that, for the construction of the multi-level operator B T S in the manifold case, we may equip L 2 ⁢ ( Γ ) with scalar product

∑ i = 1 p ∫ Ω i u ⁢ ( χ i ⁢ ( x ) ) ⁢ v ⁢ ( χ i ⁢ ( x ) ) ⁢ d x ,

which is constructed from the canonical L 2 ⁢ ( Γ ) -scalar product by simply omitting the Jacobians | ∂ ⁡ χ i ⁢ ( x ) | . With this modified scalar product, the panel-wise orthogonal projector Q T is the same as in the domain case. We conclude that the resulting B T S as in (3.2) is in L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) uniformly and that its application can be performed in linear complexity. Indeed, its implementation is equal as in the domain case as described in Section 3.2 when | T | in (3.3) is read as | χ - 1 ⁢ ( T ) | .

5 Numerical Experiments

Let Γ = ∂ [ 0 , 1 ] 3 ⊂ R 3 be the two-dimensional manifold without boundary given as the boundary of the unit cube, W := H 1 / 2 ⁢ ( Γ ) , V := H - 1 / 2 ⁢ ( Γ ) . We consider the trial space V T = S T - 1 , 0 ⊂ V of discontinuous piecewise constants. We will evaluate preconditioning of the discretized single layer operator A T ∈ L ⁢ is c ⁢ ( V T , V T ′ ) .

The role of the opposite order operator in L ⁢ is c ⁢ ( S T , 0 0 , 1 , ( S T , 0 0 , 1 ) ′ ) from Section 2.1 will be fulfilled by the multi-level operator B T S from (3.2). Equipping S T , 0 0 , 1 with the nodal basis Φ T , the matrix representation of the preconditioner G T from Section 2.1 reads as

G T = D T - 1 ⁢ ( p T ⊤ ⁢ B T S ⁢ p T + β ⁢ q T ⊤ ⁢ D T 1 / 2 ⁢ q T ) ⁢ D T - 1

for D T = diag ⁡ { | T | : T ∈ T } , uniformly sparse p T and q T as given in Section 2.1, and with the representation of the multi-level operator B T S given by

B T S = E T ⊤ ⁢ ( ∑ j = 0 L ( H j ⁢ R j - P j ⁢ H j - 1 ⁢ R j - 1 ) ⊤ ⁢ 2 - j / 2 ⁢ ( H j ⁢ R j - P j ⁢ H j - 1 ⁢ R j - 1 ) ) ⁢ E T

for the representations E T , H j , R j and P j as provided in Section 3.2 (the minor adaptations in the manifold case described in Section 4 to the matrix representations from Sections 2.1 and 3.2 vanish in the current simple case).

The BEM++ software package [8] is used to approximate the matrix representation of the discretized single layer operator A T by hierarchical matrices based on adaptive cross approximation [5, 1].

Equipping V T and R dim ⁡ V T with “energy-norms” ( A T ⋅ ( ⋅ ) or ∥ A T 1 / 2 ⋅ ∥ , respectively, we calculated the (spectral) condition numbers κ L ⁢ ( V T , V T ) ⁢ ( G T ⁢ A T ) = κ L ⁢ ( R dim ⁡ V T , R dim ⁡ V T ) ⁢ ( G T ⁢ A T ) = ρ ⁢ ( G T ⁢ A T ) ⁢ ρ ⁢ ( ( G T ⁢ A T ) - 1 ) , where ρ ⁢ ( ⋅ ) is the spectral radius, using the Lanczos method.

As initial partition T ⊥ = T 1 of Γ, we take a conforming partition consisting of 2 triangles per side, so 12 triangles in total, with an assignment of the newest vertices that satisfies the matching condition. We fixed β = 5.3 , being the value for which, for a relative small uniform refinement 𝒯 of T ⊥ , we found ρ ⁢ ( D T - 1 ⁢ p T ⊤ ⁢ B T S ⁢ p T ⁢ D T - 1 ⁢ A T ) = ρ ⁢ ( D T - 1 ⁢ β ⁢ q T ⊤ ⁢ D T 1 / 2 ⁢ q T ⁢ D T - 1 ⁢ A T ) .

5.1 Uniform Refinements

Here we let 𝕋 be the sequence { T k } k ≥ 1 of (conforming) uniform refinements, that is, T k ≻ T k - 1 is found by bisecting each triangle from T k - 1 into 2 subtriangles using newest vertex bisection.

Table 1

Spectral condition numbers of the preconditioned single layer system discretized by piecewise constants S T - 1 , 0 using uniform refinements. Preconditioner G T is constructed using the multi-level operator with β = 5.3 . The last column indicates the number of seconds per degree of freedom per application of G T .

dofs	κ S ⁢ ( A T )	κ S ⁢ ( G T ⁢ A T )	sec / dof
12	14.5	2.6	2.6 ⋅ 10 - 5
48	31.0	2.7	1.4 ⋅ 10 - 5
192	59.9	2.8	4.9 ⋅ 10 - 6
768	118.7	3.3	1.4 ⋅ 10 - 6
3072	234.6	3.8	6.3 ⋅ 10 - 7
12288	450.4	4.1	3.3 ⋅ 10 - 6
49152	852.5	4.3	6.5 ⋅ 10 - 7
196608	1566.4	4.5	7.3 ⋅ 10 - 7
786432	2730.5	4.6	7.8 ⋅ 10 - 7

Table 1 shows the condition numbers of the preconditioned system in this situation. The condition numbers are relatively small, and the timing results show that the implementation of the preconditioner is indeed linear.

5.2 Local Refinements

Here we take 𝕋 as a sequence { T k } k ≥ 1 of (conforming) locally refined partitions, where T k ≻ T k - 1 is constructed by applying newest vertex bisection to all triangles in T k - 1 that touch a corner of the cube.

Table 2 contains results for the preconditioned single layer operator discretized by piecewise constants S T - 1 , 0 . The preconditioned condition numbers are nicely bounded, and the timing results confirm that our implementation of the preconditioner is of linear complexity, also in the case of locally refined partitions.

Table 2

Spectral condition numbers of the preconditioned single layer system discretized by piecewise constants S T - 1 , 0 using local refinements at each of the eight cube corners. Operator G T is applied using the multi-level operator with β = 5.3 . The second column is defined by h T , min := min T ∈ T ⁡ | T | . The last column indicates the number of seconds per degree of freedom per application of G T .

dofs	h T , min	κ S ⁢ ( G T ⁢ A T )	sec / dof
12	1.4 ⋅ 10 0	2.63	2.5 ⋅ 10 - 5
336	8.8 ⋅ 10 - 2	2.73	2.4 ⋅ 10 - 6
720	5.5 ⋅ 10 - 3	2.91	1.8 ⋅ 10 - 6
1104	3.4 ⋅ 10 - 4	2.96	1.8 ⋅ 10 - 6
1488	2.1 ⋅ 10 - 5	2.99	2.2 ⋅ 10 - 6
1872	1.3 ⋅ 10 - 6	2.98	2.0 ⋅ 10 - 6
2256	8.4 ⋅ 10 - 8	3.00	2.3 ⋅ 10 - 6
2640	5.2 ⋅ 10 - 9	3.00	2.0 ⋅ 10 - 6
3024	3.2 ⋅ 10 - 10	3.01	2.3 ⋅ 10 - 6
3408	2.0 ⋅ 10 - 11	3.01	2.5 ⋅ 10 - 6
3696	2.5 ⋅ 10 - 12	3.01	2.6 ⋅ 10 - 6

Funding source: Nederlandse Organisatie voor Wetenschappelijk Onderzoek

Award Identifier / Grant number: 613.001.652

Funding statement: The second author has been supported by the Netherlands Organization for Scientific Research (NWO) under contract no. 613.001.652.

References

[1] M. Bebendorf, Approximation of boundary element matrices, Numer. Math. 86 (2000), no. 4, 565–589. 10.1007/PL00005410Search in Google Scholar

[2] J. H. Bramble, J. E. Pasciak and P. S. Vassilevski, Computational scales of Sobolev norms with application to preconditioning, Math. Comp. 69 (2000), no. 230, 463–480. 10.1090/S0025-5718-99-01106-0Search in Google Scholar

[3] S. H. Christiansen and J.-C. Nédélec, Des préconditionneurs pour la résolution numérique des équations intégrales de frontière de l’acoustique, C. R. Acad. Sci. Paris Sér. I Math. 330 (2000), no. 7, 617–622. 10.1016/S0764-4442(00)00225-1Search in Google Scholar

[4] T. Führer, A. Haberl, D. Praetorius and S. Schimanko, Adaptive bem with inexact pcg solver yields almost optimal computational costs, Numer. Math. 141 (2019), no. 4, 967–1008. 10.1007/s00211-018-1011-1Search in Google Scholar

[5] W. Hackbusch, A Sparse Matrix Arithmetic Based on ℋ-Matrices. Part I: Introduction to ℋ-Matrices, Computing 62 (1999), no. 2, 89–108. 10.1007/s006070050015Search in Google Scholar

[6] R. Hiptmair, Operator preconditioning, Comput. Math. Appl. 52 (2006), no. 5, 699–706. 10.1016/j.camwa.2006.10.008Search in Google Scholar

[7] J. Maubach, Local bisection refinement for 𝑛-simplicial grids generated by reflection, SIAM J. Sci. Comput. 16 (1995), no. 1, 210–227. 10.1137/0916014Search in Google Scholar

[8] W. Śmigaj, T. Betcke, S. Arridge, J. Phillips and M. Schweiger, Solving boundary integral problems with BEM++, ACM Trans. Math. Softw. 41 (2015), 10.1145/2590830. 10.1145/2590830Search in Google Scholar

[9] R. P. Stevenson, The completion of locally refined simplicial partitions created by bisection, Math. Comp. 77 (2008), 227–241. 10.1090/S0025-5718-07-01959-XSearch in Google Scholar

[10] R. P. Stevenson and R. van Venetië, Uniform preconditioners for problems of negative order, Math. Comp. 89 (2020), 645–674. 10.1090/mcom/3481Search in Google Scholar

[11] R. P. Stevenson and R. van Venetië, Uniform preconditioners for problems of positive order, Comput. Math. Appl. 79 (2020), no. 12, 3516–3530. 10.1016/j.camwa.2020.02.009Search in Google Scholar

[12] O. Steinbach and W. L. Wendland, The construction of some efficient preconditioners in the boundary element method, Adv. Comput. Math. 9 (1998), no. 1–2, 191–216. 10.1023/A:1018937506719Search in Google Scholar

[13] C. T. Traxler, An algorithm for adaptive mesh refinement in 𝑛 dimensions, Computing 59 (1997), no. 2, 115–137. 10.1007/BF02684475Search in Google Scholar

[14] J. Wu and H. Zheng, Uniform convergence of multigrid methods for adaptive meshes, Appl. Numer. Math. 113 (2017), 109–123. 10.1016/j.apnum.2016.11.005Search in Google Scholar

Received: 2020-04-10

Revised: 2020-08-25

Accepted: 2020-11-16

Published Online: 2020-12-02

Published in Print: 2021-04-01

This work is licensed under the Creative Commons Attribution 4.0 International License.

Uniform Preconditioners of Linear Complexity for Problems of Negative Order

Abstract

1 Introduction

1.1 Outline

1.2 Notation

2 Preconditioning

2.1 Construction of Optimal Preconditioners

2.2 Implementation of G T

Theorem 2.1 ([10, Section 3])

2.2.1 Piecewise Constant Trial Space V T

2.3 Higher Order Case

3 An Operator B T S of Multi-level Type

3.1 Definition and Analysis of B T S

Theorem 3.1 ([14, Lemma 3.7])

Proof

Proposition 3.2 ([14, Lemma 3.1])

3.2 Implementation of B T S

Proof

4 Manifold Case

5 Numerical Experiments

5.1 Uniform Refinements

5.2 Local Refinements

References

Journal and Issue

Articles in the same Issue