Abstract
We propose a multilevel 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 multilevel 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].
1 Introduction
In this work, we construct a multilevel 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 multilevel 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 multilevel 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 multilevel operator
B
T
S
induces a uniform preconditioner
G
T
, i.e.,
G
T
A
T
is uniformly wellconditioned, 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 multilevel 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 multilevel subspace decompositions. Examples are given by wavelet preconditioners, or closely related, the preconditioners from [2], for the latter assuming quasiuniform partitions.
For

s
=

1
2
, an optimal multilevel preconditioner based on a nonoverlapping subspace decomposition for operators defined on the boundary of a 2 or 3dimensional 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 multilevel 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 Multilevel 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 multilevel type. The operator
B
T
S
is based on a stable multilevel decomposition of
S
T
,
0
0
,
1
given by Wu and Zheng [14]. Usually, such a stable multilevel 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 quasiinterpolator
Π
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 multilevel decomposition of
S
T
,
0
0
,
1
.
Theorem 3.1
Theorem 3.1 ([14, Lemma 3.7])
For the averaging quasiinterpolator
Π
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 wellknown 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 multilevel 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 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.
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 elementwise 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.
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 nonzero 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 nonzero 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
or
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 meshdependent 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.
It remains to discuss the construction of an operator
B
T
S
of multilevel 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 panelwise
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 multilevel 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 panelwise 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 twodimensional 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 multilevel 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 multilevel 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 “energynorms”
(
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 multilevel 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 multilevel 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
References
[1] M. Bebendorf, Approximation of boundary element matrices, Numer. Math. 86 (2000), no. 4, 565–589. Search 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. Search 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. Search 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. Search 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. Search in Google Scholar
[6] R. Hiptmair, Operator preconditioning, Comput. Math. Appl. 52 (2006), no. 5, 699–706. Search 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. Search 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. Search in Google Scholar
[9] R. P. Stevenson, The completion of locally refined simplicial partitions created by bisection, Math. Comp. 77 (2008), 227–241. Search in Google Scholar
[10] R. P. Stevenson and R. van Venetië, Uniform preconditioners for problems of negative order, Math. Comp. 89 (2020), 645–674. Search 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. Search 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. Search in Google Scholar
[13] C. T. Traxler, An algorithm for adaptive mesh refinement in 𝑛 dimensions, Computing 59 (1997), no. 2, 115–137. Search in Google Scholar
[14] J. Wu and H. Zheng, Uniform convergence of multigrid methods for adaptive meshes, Appl. Numer. Math. 113 (2017), 109–123. Search in Google Scholar
Received: 20200410
Revised: 20200825
Accepted: 20201116
Published Online: 20201202
Published in Print: 20210401
© 2020 Stevenson and van Venetië, published by De Gruyter
This work is licensed under the Creative Commons Attribution 4.0 International License.