Multilevel Solvers for Elliptic Problems on Domains

Peter Oswald peter.oswald@gmd.de    Institute for Algorithms and Scientific Computing, GMD - German National Research Center for Information Technology, D-53754 Sankt Augustin, Germany

$1 Introduction

For the numerical solution of elliptic boundary value problems for partial differential equations, multilevel methods have gained popularity over the last decade. This is mainly due to their nearly optimal complexity for a number of model problems. In many practical cases, they are based on multiresolution scale of nested finite-dimensional subspaces

0⊂V1⊂···⊂Vj⊂···

of a Hilbert space V serving as the energy space for the given variational problem. The scale (1.1) is used to produce stable subspace splittings

J=V0+V1+···+VJ

and to design fast iterative solvers related to such splittings for the discretized variational problem associated with a computational discretization space VJ (or a properly defined subspace J*⊂VJ). For 2mth order elliptic boundary value problems in Sobolev spaces, considerable progress has been achieved in the theoretical understanding of multilevel and multigrid methods as well as of other subspace correction methods for finite element discretizations (see [65, 67, 11, 60]). We survey some of these results in Section 2.

The underlying theory also applies to various wavelet discretizations; see, e.g., [26, 44, 45, 29] for some papers that deal with wavelet solvers for elliptic problems and are related to our approach. Roughly speaking, in these algorithms suitable “detail spaces” Wj ⊂ Vj are constructed, together with their algebraic bases, such that Vj = Vj–1 + Wj provides a stable splitting of Vj into “low frequency” (Vj–1) and “high frequency” (Wj) parts. Using this two-level decomposition recursively, (1.2) is replaced by

J=V0+˙W1+˙···+˙WJ.

One can use this splitting indirectly (i.e., the original problem is discretized with respect to a standard basis in VJ, and the wavelet decomposition behind (1.3) implicitly defines the structure of the multilevel preconditioner) or directly. In the latter case, the discretization is performed with respect to the wavelet basis and is automatically well conditioned. To achieve asymptotically optimal work estimates one has to use compression arguments. For elliptic problems of order 2 m, and Vj with locally supported basis functions, the first approach is often preferred since the VJ-discretization matrix is automatically sparse and available from standard engineering codes. The direct use of the wavelet decomposition (1.3) is promising for those situations in which the VJ-discretization is not a priori sparse, e.g., for integral equations.

In both cases, the explicit introduction of detail spaces Wj is the crucial step, and may add some theoretical and practical difficulties. For example, most of the popular examples of wavelet spaces (see [30, 16]) are derived in a one-dimensional, shift-invariant setting on . Multivariate examples on d are mostly obtained by tensor-product techniques. Adaptations to bounded intervals and domains have been studied in, e.g., [2, 17, 19, 18]. However, up to now there has been no comprehensive study of the practical potential of discretizations using multilevel structures based on shift-invariance and dyadic dilation (modulo boundary modifications) in the case of general, nonrectangular geometries. It is not completely clear to the author what will be left from the powerful wavelet machinery if the basic algebraic assumptions (invariance with respect to integer shifts and (dyadic) dilation) are significantly relaxed.

In our opinion, the departure from these assumptions is unavoidable for many engineering problems. A simple example is 2mth order elliptic PDEs with rapidly varying coefficients that exhibit a large ratio of ellipticity constants and are therefore far from a generic Hm-problem. Principal difficulties are to be expected if nonsymmetric equations with dominant low order hyperbolic parts, such as convection-diffusion equations, or nonlinear problems are to be studied. In our opinion, this robustness aspect is one of the target problems for future investigation, i.e., the adaptation of the multilevel concept to a class of operator equations (or even to an individual equation) still remains a decisive issue for practical implementations and engineering applications. An indicator of this tendency is analogous efforts within the FEM community and the renewed interest in the algebraic multigrid method and related algorithms.

There is another observation that dampens expectations concerning the practical use of wavelet solvers, even for standard symmetric elliptic boundary value problems. Numerical experiments [39, 49] show that for generic H1-problems, i.e., for second order elliptic equations, wavelet and prewavelet discretizations perform slightly worse than traditional finite element preconditioners associated with (1.2). While condition numbers of L2-problems (and sometimes also of Hs-problems with negative s) are usually improved and become uniformly bounded if J → ∞, the preconditioning effect in the H1-norm is reduced by a constant factor. This poses the problem of determining more carefully problem classes in which the use of wavelet preconditioners based on (1.3) is justified compared to simpler methods based only on the use of scaling functions and related to (1.2). Also, one may argue whether new wavelet families (Daubechies orthogonal wavelets [30], AFIF elements [51], etc.) are generally useful and are able to compete with traditional finite element and spline constructions from this more practical viewpoint.

The reader should not expect an answer to these more philosophical questions. In Sections 3–5, we concentrate on studying the influence of the domain geometry on the optimality of multilevel preconditioners resulting from a standard multiresolution analysis (these sections are sometimes rather technical and represent the original part of this paper). The concept is based on sequences {j} of subspaces of a fixed sequence {Vj} of subspaces in some Hs(d) which “live” on a uniform structure generated by shift-invariance principles and dyadic dilation as usual. These auxiliary subspaces j are spans of scaling functions (or, in the finite element terminology, nodal basis functions) of levels ≤ j such that a canonical Hs(d)-elliptic Galerkin discretization can be solved efficiently, e.g., by preconditioned iterative methods, with the preconditioner inherited from the generating system, or frame, consisting of scaling functions (see Section 4). Thus, we essentially stay with splittings of the type (1.2).

The connection of the auxiliary problems in j with the originally given Hs(Ω)-elliptic problem, with natural boundary conditions, on a generic bounded domain Ω ⊂ d is established in Section 3 where sufficiently rich subspaces Vj,Ω ⊂ Vj |Ω will be constructed. The construction consists of a local boundary modification which is similar to constructions outlined in other papers on wavelets on intervals and domains, too, but is relatively simple. A drawback is that, in contrast to {j}, the sequence Vj,Ω is not monotone, which requires additional considerations when adaptivity is an issue. Information between Vj,Ω and j is exchanged by appropriate restriction (Rj) and extension (Ej) operators. To obtain uniform condition number and work estimates, certain geometric conditions on Ω arise in a natural way. Asymptotically, they hold for...