GOAFEM for Nonsymmetric Problems

Convergence of Goal-Oriented Adaptive Finite Element Methods for Nonsymmetric Problems

Michael Holst mholst@math.ucsd.edu  and  Sara Pollock snpolloc@math.ucsd.edu Department of Mathematics
University of California San Diego
La Jolla CA 92093
July 15, 2019
Abstract.

In this article we develop convergence theory for a class of goal-oriented adaptive finite element algorithms for second order nonsymmetric linear elliptic equations. In particular, we establish contraction results for a method of this type for Dirichlet problems involving the elliptic operator with Lipschitz, almost-everywhere symmetric positive definite, with divergence-free, and with . We first describe the problem class and review some standard facts concerning conforming finite element discretization and error-estimate-driven adaptive finite element methods (AFEM). We then describe a goal-oriented variation of standard AFEM (GOAFEM). Following the recent work of Mommer and Stevenson for symmetric problems, we establish contraction of GOAFEM and convergence in the sense of the goal function. Our analysis approach is signficantly different from that of Mommer and Stevenson, combining the recent contraction frameworks developed by Cascon, Kreuzer, Nochetto and Siebert; by Nochetto, Siebert and Veeser; and by Holst, Tsogtgerel and Zhu. We include numerical results demonstrating performance of our method with standard goal-oriented strategies on a convection problem .

Key words and phrases:
Adaptive methods, elliptic equations, non-symmetric problems, quasi-orthogonality, duality, approximation theory, residual-based error estimator, convergence, contraction, optimality, a priori estimates, a posteriori estimates, goal oriented
MH was supported in part by NSF Awards 0715146 and 0915220, and by DOD/DTRA Award HDTRA-09-1-0036.
SP was supported in part by NSF Award 0715146.

1. Introduction

In this article we develop convergence theory for a class of goal-oriented adaptive finite element methods for second order nonsymmetric linear elliptic equations. In particular, we report contraction results for a method of this type for the problem

(1.1)
(1.2)

with a polyhedral domain, , with Lipschitz, almost-everywhere symmetric positive definite (SPD), with divergence-free, and with . The standard weak formulation of this problem reads: Find such that

(1.3)

where

(1.4)

Our approach is to first describe the problem class in some detail, and review some standard facts concerning conforming finite element discretization and error-estimate-driven adaptive finite element methods (AFEM). We will then describe a goal-oriented variation of standard AFEM (GOAFEM). Following the recent work of Mommer and Stevenson [18] for symmetric problems, we establish contraction of GOAFEM and convergence in the sense of the goal function. Our analysis approach is signficantly different from that of Mommer and Stevenson [18], combining the recent contraction frameworks of Cascon, Kreuzer, Nochetto and Siebert  [7], of Nochetto, Siebert and Veeser [19], and of Holst, Tsogtgerel and Zhu [16]. We also give some numerical results comparing our goal-oriented method both to the one presented in [18] and the dual weighted residual (DWR) method as in [2, 4, 9, 13, 14, 10], among others. Unlike the existing literature on the DWR method, we prove strong convergence of our goal-oriented method. We establish contraction of the goal error in terms of the energy norm errors and error estimators of the primal and dual problems, and indicate how this implies optimality in terms of the global error. Controlling this overestimate of the error shows convergence of the method to the goal, although not optimality in this sense. Our numerical results demonstrate, however, that the algorithm presented here performs at least comparably to and in some cases better than DWR and the method in [18] on a variety of convection dominated linear problems.

The goal-oriented problem concerns achieving a target quality in a given linear functional of the weak solution of the problem (1.3). For example, , the average value of over some subdomain . By writing down the adjoint operator, , we consider the adjoint or dual problem: find such that . It has been shown for the symmetric form () of problem (1.1)–(1.2) with piecewise constant SPD diffusion cofficient (and with ), that by solving the primal and dual problems simultaneously, one may converge to an approximation of faster than by approximating and then , when forcing contraction in only the primal problem [18]. We will follow the same general approach to establish similar goal-oriented AFEM results for nonsymmetric problems. In order to handle nonsymmetry, we will follow the technical approach in [17, 7, 16], and rely largely on establishing quasi-orthogonality. Contraction results are established in [17, 7] for (1.1)–(1.2) in the case that is SPD, Lipschitz or piecewise Lipschitz, is divergence-free, and . In [16], quasi-orthogonality is used as the basis for establishing contraction of AFEM for two classes of nonlinear problems. As in these earlier efforts, relying on quasi-orthogonality will require that we assume that the initial mesh is sufficiently fine, and that the solution to the dual problem is sufficiently smooth, e.g. in .

Following [16], the contraction argument developed in this paper will follow from first establishing three preliminary results for two successive AFEM approximations and , and then applying the Dörfler marking strategy:

  • Quasi-orthogonality (§3.1): There exists such that

  • Error estimator as upper bound on error (§3.2): There exists such that

  • Estimator reduction (§3.4): For the marked set that takes refinement , for positive constants and any

The marking strategy used is the original Dörfler strategy; elements are marked for refinement based on indicators alone. The marked set must satisfy

In this goal-oriented method, a second marked set is chosen based on an error indicator for the dual problem associated with the given goal functional, and the union of the two marked sets is then used for refinement. A main advantage of the approach in [7, 16] is that it does not require an interior node property. This allows us to establish the necessary results for contraction without taking full refinements of the mesh at each iteration. This improvement follows from the use of the local perturbation estimate or local Lipschitz property rather than the estimator as lower bound on error. We use the standard lower bound estimate as found in [17] for optimality arguments in the second part of the paper concerning quasi-optimality of the method.

There are three main notions of error used throughout this paper. The energy error , the quasi-error, and the total-error. The energy error is defined by the symmetric part of the bilinear form that arises from the given differential operator in (1.3). The quasi-error is the sum of the energy-error and scaled error estimator

and this is the quantity that is reduced at each iteration of the algorithm. In §3 the quasi-error is shown to satisfy

The total error includes the oscillation term rather than the estimator

The oscillation term captures the higher-frequency oscillations in the residual missed by the averaging of the finite element method. While the quasi-error is the focus of the contraction arguments, the total error is used in our discussion of complexity analysis.

Throughout this paper, the constant will denote a generic but global constant that may depend on the data and the condition of the initial mesh , and may change as an argument proceeds, without danger of confusion.

Outline of the paper. The remainder of the paper is structured as follows. In §2, we first describe the problem class and review some standard facts concerning conforming finite element discretization and error-estimate-driven adaptive finite element methods (AFEM). In §2.3, we then describe a goal-oriented variation of the standard approach to AFEM (GOAFEM). Following the recent work of Mommer and Stevenson for symmetric problems, in §3 we establish contraction of goal-oriented AFEM. We also then show convergence in §3.6 in the sense of the goal function. Our analysis combines the recent contraction frameworks developed in [7, 19, 16], applied now to the goal oriented problem. In §5, we present some numerical experiments comparing the method presented here with two standard goal oriented strategies. We recap the results in §6, and point out some remaining open problems.

2. Problem class, discretization, goal-oriented AFEM

2.1. Problem class, weak formulation, spaces and norms

Consider the nonsymmetric problem (1.3), where as in (1.4) we have

Here we have introduced the notation for the inner-product over . The adjoint or dual problem is: Find

(2.1)

where is the formal adjoint of , and where the functional is defined through

(2.2)

for some given . We will make the following assumptions on the data:

Assumption 2.1 (Problem data).

The problem data and dual problem data satisfy

  • , Lipschitz, and a.e. symmetric positive-definite:

    (2.3)
    (2.4)
  • , with , and divergence-free.

  • , with , and .

  • .

The native norm is the Sobolev norm given by

(2.5)

The norm of a vector valued function over domain is defined here as the norm of the norm of each component

(2.6)

Similarly, the norm of a matrix valued function over domain is defined as the Frobenius norm of the norm of each component

(2.7)

We note that one could employ other equivalent discrete norms in the definitions (2.1) and (2.1), however this choice simplifies the analysis.

Continuity of follows from the Hölder inequality, and bounding the norm of the function and its gradient by the norm

(2.8)

Coercivity follows from the Poincaré inequality with constant and the divergence-free condition

(2.9)

where the coercivity constant . Continuity and coercivity imply existence and uniqueness of the solution by the Lax-Milgram Theorem [12]. The adjoint operator is given by

Integration by parts on the convection term and the divergence-free condition imply

(2.10)

Define the energy semi-norm by

(2.11)

Non-negativity follows directly from the coercivity estimate (2.9)

(2.12)

which establishes the energy semi-norm as a norm. Putting this together with the reverse inequality

(2.13)

establishes the equivalence between the native and energy norms with the constant

2.2. Finite element approximation

We employ a standard conforming piecewise polynomial finite element approximation below.

Assumption 2.2 (Finite element mesh).

We make the following assumptions on the underlying simplex mesh:

  • The initial mesh is conforming.

  • The mesh is refined by newest vertex bisection [5], [18] at each iteration.

  • The initial mesh is sufficiently fine. In particular, it satisfies (3.6).

Based on assumptions 2.2 we have the following mesh constants.

  • Define

    (2.14)

    In particular, is the initial mesh diameter.

  • Define the mesh constant where and then for any two elements in the same generation and as neighboring elements may differ by at most one generation for any two neighboring elements

    (2.15)
  • The minimal angle condition satisfied by newest vertex bisection implies the meshsize is comparable to , the size of any true-hyperface of . In particular, there is a constant

    (2.16)

Let the set of conforming meshes derived from the initial mesh . Define by For a conforming mesh with a conforming refinement we say . The set of refined elements is given by

(2.17)

Define the finite element space

(2.18)

For subsets ,

(2.19)

where is the space of polynomials degree degree over . Denote the patch about

(2.20)

For a -simplex , an true-hyperface is a dimensional face of , e.g., a face in 3D or an edge in 2D. Define the discrete primal problem: Find

(2.21)

and the discrete dual problem

(2.22)

2.3. Goal oriented AFEM (GOAFEM)

As in [18] the goal oriented adaptive finite element method (GOAFEM) is based on the standard AFEM algorithm:

(2.23)

In the goal oriented method, one enforces contraction of the quasi-error in both the primal problem and an associated dual problem. As shown in section §3.6, the error in the goal-function satisfies the bound

This motivates driving down the energy-error in both the primal and dual problems at each iteration. As noted in [7] the residual-based error estimator does not exhibit monotone behavior in general, although it is monotone non-increasing with respect to nested mesh refinement when applied to the same (coarse) function. The quasi-error is shown to contract for each problem for which mesh refinement satisfies the Dörfler property. However, refining the mesh with respect to the primal problem does not guarantee the quasi-error in the dual problem will be non-increasing, and vice-versa. As such, the procedures SOLVE and ESTIMATE are performed for each of the primal and dual problems. The marked set is taken to be the union of marked sets from the primal and dual problems, each chosen to satisfy the Dörfler property. This method produces a sequence of refinements for which the quasi-error in the both the primal and dual problems contract at each step. The requirement to reduce the quasi-error rather than the energy error as in [18] is why the marking strategy in this method differs from the one shown effective for the Laplacian. Our numerical results demonstrate similar behavior of both methods, although the method presented here has the advantage that the code takes fewer iterations of (2.23) to achieve similar results.

Procedure SOLVE. The contraction result supposes the exact Galerkin solution is found on each mesh refinement.

Procedure ESTIMATE. The estimation of the error on each element is determined by a standard residual-based estimator. The residuals over element interiors and jump-residuals over the boundaries are based on the local strong forms of the elliptic operator and its adjoint as follows.

(2.24)

The residuals for the primal and dual problems using the sign convention in [7] are:

(2.25)

While the primal and dual solutions of (1.3) and (2.1) respectively satisfy

the residuals for the primal and dual problems are in general different. The jump residual for the primal and dual problems is

(2.26)

where jump operator is given by

(2.27)

and is taken to be the appropriate outward normal defined piecewise on . On boundary edges we have

so that . For clarity, we will also employ the notation

and similarly for the other strong form operators. The error indicator is given as

(2.28)

The dual error-indicator is then given by

(2.29)

The error estimators are given by the sum of error indicators over elements in the space where or .

(2.30)

The dual energy estimator is:

(2.31)

The contraction results for the quasi-error presented below will be shown to hold for where the error estimator and oscillation are defined in terms of the norm. While complexity results are shown only for , the contraction results for are useful for nonlinear problems; see [16].

For analyzing oscillation, for let the orthogonal projector defined by the best approximation in over mesh and . Define now the oscillation on the elements for the primal problem by

(2.32)

and analogously for the dual problem. For subsets set

(2.33)

The data estimator and data oscillation, identical for both the primal and dual problems, are given by

(2.34)
(2.35)

The data estimator and oscillation over the mesh or a subset are given by the maximum data estimator (oscillation) over elements in the mesh or subset: For

The data estimator and data oscillation on the initial mesh

As the grid is refined, the data estimator and data oscillation terms satisfy the monotonicity property [7] for refinements

(2.36)

Procedure MARK. The Dörfler marking strategy for the goal-oriented problem is based on the following steps as in [18]:

  • Given , mark sets for each of the primal and dual problems:

    • Mark a set such that,

      (2.37)
    • Mark a set such that,

      (2.38)
  • Let the union of sets found for the primal and dual problems respectively.

The set differs from that in [18], where the set of lesser cardinality between is used. In the case of the nonsymmetric problem the error reduced at each iteration is the quasi-error rather than the energy error as in the symmetric problem [18]. This error for each problem is guaranteed to contract based on the refinement satisfying the Dörfler property. As such, refining the mesh with respect to one problem does not guarantee the quasi-error in the other problem is nonincreasing. Sets with optimal cardinality (up to a factor of 2) can be chosen in linear time  by binning the elements rather than performing a full sort [18].

Procedure REFINE. The refinement (including the completion) is performed according to newest vertex bisection [5]. The complexity and other properties of this procedure are now well-understood, and will simply be exploited here.

3. Contraction and convergence theorems

The key elements of the main contraction argument constructed below are quasi-orthogonality 3.1, error estimator as upper-bound on energy-norm error 3.2 and estimator reduction 3.4. Estimator-reduction is shown via the local-perturbation estimate 3.3. The local perturbation of the oscillation is presented here and used in §4. Mesh refinements (respectively ) are assumed conforming, and is assumed the Galerkin solution on refinement . The following results hold for both the primal and dual problems which differ by the sign of the convection term; therefore, they are established here only for the primal problem.

3.1. Quasi-orthogonality

Orthogonality in the energy-norm does not generally hold in the nonsymmetric problem. We use the weaker quasi-orthogonality result to establish contraction of AFEM (GOAFEM). The following is a variation on Lemma 2.1 in [17] (see also [16]).

Lemma 3.1 (Quasi-orthogonality).

Let the problem data satisfy Assumption 2.1 and the mesh satisfy conditions (1) and (2) of Assumption 2.2. Let with . Let the solution to (2.21), . There exists a constant depending on the problem data and initial mesh , and a number dictated only by the angles of , such that if the meshsize of the initial mesh satisfies , then

(3.1)

where

Equality holds (usual orthogonality) when in , in which case the problem is symmetric.

Proof.

The proof follows close that of Lemma 2.1 in [17]. Let

By Galerkin orthogonality

(3.2)

Rearranging and applying the divergence-free condition on the convection term

Applying Hölder’s inequality and coercivity (2.9) followed by Young’s inequality with constant to be determined,

(3.3)

By a duality argument for some assuming for some depending on the angles of

(3.4)

The details of this argument as described in the appendix §7 may also be found in [1] and [8]. Applying (3.4) and (3.3) to (3.2),

(3.5)

Choose to equate coefficients

then

Assuming the initial mesh as characterized by satisfies

(3.6)

the quasi-orthogonality result holds. ∎

Note that by (3.2) we also have

(3.7)

Similarly to (3.3)

(3.8)

which under the same assumptions yields the estimate

(3.9)

where .

3.2. Error estimator as global upper-bound

We now recall the property that the error estimator is a global upper bound on the error. The proof is fairly standard; see e.g. [18] (Proposition 4.1), [17] (3.6), and [16].

Lemma 3.2 (Error estimator as global upper-bound).

Let the problem data satisfy Assumption 2.1 and the mesh satisfy conditions (1) and (2) of Assumption 2.2. Let with . Let the solution to (2.21), and the solution to (1.3). Let

Then for global constant depending on the problem data and initial mesh

(3.10)

and in particular

(3.11)

3.3. Local perturbation

The local perturbation property established in [7], analogous to the local Lipshitz property in [16], is a key step in establishing the contraction result. This is a minor variation on Proposition 3.3 in [7] which deals with a symmetric problem. Here, we include a convection term in the estimate. In particular, (3.12) shows that the difference in the error indicators over an element between two functions in a given finite element space may be bounded by a fixed factor of the native norm over the patch of the difference in functions. In contrast with the analogous result in [7] the estimate (3.13) involves a fixed factor of the native norm over an individual element rather than a patch as by the continuity of the oscillation term does not involve the jump residual.

We include the proof of  (3.12) for completeness. The proof of (3.13) may be found in [7] with the final result inferred by the absence of the jump residual in the oscillation term.

Lemma 3.3 (Local perturbation).

Let the problem data satisfy Assumption 2.1 and the mesh satisfy condition (1) of Assumption 2.2. Let . For all and for any

(3.12)
(3.13)

where recalling (2.20) is the union of with elements in sharing a true-hyperface with . The constants depend on the initial mesh , the dimension and the polynomial degree .

Proof of (3.12).

From (2.28)

(3.14)

Denote by . Set . By linearity

and

For by the triangle inequality

For using the generalized triangle-inequality

(3.15)

we have

Consider the second term on the RHS . By definition (2.24) of , the product rule applied to the diffusion term and the triangle-inequality

where is the Hessian of . Consider each term. The first diffusion term

(3.16)

by the inequality

(3.17)

Applying (3.17) and inverse-estimate [6] to the second diffusion term

(3.18)

For the reaction term

(3.19)

For the convection term applying (3.17)

(3.20)

Consider the the jump-residual term . For each interior true-hyperface