A priori error estimates for compatible spectral discretization of the Stokes problem

A priori error estimates for compatible spectral discretization of the Stokes problem for all admissible boundary conditions

Jasper Kreeft Delft University of Technology, Faculty of Aerospace Engineering,
Kluyverweg 2, 2629 HT Delft, The Netherlands.
 and  Marc Gerritsma
September 15, 2019

This paper describes the recently developed mixed mimetic spectral element method for the Stokes problem in the vorticity-velocity-pressure formulation. This compatible discretization method relies on the construction of a conforming discrete Hodge decomposition, that is based on a bounded projection operator that commutes with the exterior derivative. The projection operator is the composition of a reduction and a reconstruction step. The reconstruction in terms of mimetic spectral element basis-functions are tensor-based constructions and therefore hold for curvilinear quadrilateral and hexahedral meshes.

For compatible discretization methods that contain a conforming discrete Hodge decomposition, we derive optimal a priori error estimates which are valid for all admissible boundary conditions on both Cartesian and curvilinear meshes. These theoretical results are confirmed by numerical experiments. These clearly show that the mimetic spectral elements outperform the commonly used -compatible Raviart-Thomas elements.

Key words and phrases:
Stokes problem, mixed finite elements, mimetic/compatible discretization, error estimates
1991 Mathematics Subject Classification:
Primary 76D07, 65N30; Secondary 65M70, 12Y05, 13P20
Jasper Kreeft is funded by STW Grant 10113
This paper is in final form and no version of it will be submitted for publication elsewhere.

]J.J.kreeft@gmail.com, M.I.Gerritsma@TUDelft.nl

1. Introduction

Let , , be a bounded contractible domain with boundary . On this domain we consider the Stokes problem, consisting of the equations for conservation of momentum and for conservation of mass,


where the stress tensor is given by


with the velocity vector, the pressure, the forcing term, the mass source and the kinematic viscosity. For analysis purposes we choose .

This paper considers the recently developed mixed mimetic spectral element method (MMSEM) [40, 41]. This compatible finite/spectral element method is based on the compatible discretization of the exterior derivative from differential geometry, which represents the vector operators, grad, curl and div. The Stokes problem expressed in terms of these vector operations is known as the vorticity-velocity-pressure (VVP) formulation, [9, 23]. For the VVP formulation, the Laplace operator is split using the vector identity, , and by introducing vorticity as auxiliary variable, . The VVP formulation of the Stokes problem becomes


Following [9, 40] we make a distinction between the operators grad, curl and div, that correspond to the classical Newton-Leibnitz, Stokes circulation and Gauss divergence theorems, and the operators -grad, curl and -div that are their formal Hilbert adjoints,

The distinction between the two types of differential operators is made explicitly, because the construction of our conforming finite element spaces relies on the three mentioned integration theorems, while the mixed formulation relies on the formal Hilbert adjoint relations. While in vector calculus this distinction is not common, in differential geometry these structures naturally appear since they make a clear distinction between metric-free (topological) and metric-dependent operations.

The MMSEM is a compatible discretization method that relies on the construction of a conforming discrete Hodge-decomposition, which implies a discrete Poincaré inequality. It requires the development of a bounded projection operator that commutes with the exterior derivative. The bounded projection is a composition of a reduction by means of integration and mimetic spectral element basisfunctions as reconstruction.

The reduction onto -dimensional submanifolds result in the discrete unknowns representing integral quantities. This is one of the major differences with related methods as the Marker and Cell scheme [32] and the lowest-order Raviart-Thomas and Nédélec compatible finite elements [45, 48], where use is made of averaged quantities.

The basis functions, used for the reconstruction, are constructed using tensor products of one dimensional nodal and edge interpolation basis functions [28], and therefore hold for quadrilateral and hexahedral meshes. They belong to the class of compatible finite elements, and were constructed based on the mimetic framework first described in [38] and later extended in [11]. The mimetic framework, including the mimetic spectral elements, were extensively described in [41]. This mimetic framework relies on the languages of differential geometry instead of vector calculus, and algebraic topology as its discrete counterpart.

The use of differential geometry and algebraic topology enjoys increasing popularity for the development of compatible schemes, [5, 6, 11, 12, 13, 21, 34, 35]. Compatible discretizations are often combined with mixed formulations. Mixed formulations are described extensively in among others [15, 30] and in terms of differential forms in [5, 6] for the Hodge-Laplacian and in [40] for the VVP formulation of the Stokes problem.

The MMSEM contains compatible finite elements that are compatible with all admissible types of boundary conditions for the Stokes problem in VVP formulation. We will show that the method obtains optimal rates of convergence for all variables on curvilinear meshes and for all admissible boundary conditions, i.e. standard and nonstandard. It is therefore extending the error estimates found in literature, which are often specifically constructed for certain types of boundary conditions, [1, 4, 10, 14, 24, 29]. To show optimal convergence a priori error estimates are derived.

This is an improvement with respect to the well-known Raviart-Thomas compatible finite elements. These are not compatible in case of Dirichlet boundary conditions and therefore lead to suboptimal convergence behavior, as was shown in [4, 24]. This non-compatibility results in a decrease in rate of convergence of maximal order.

From a physical/fluid dynamics point-of-view the new method is relevant because it combines optimal convergence with a pointwise divergence-free discretization (in absence of any mass source) of arbitrary order on curvilinear meshes, valid for all allowable types of boundary conditions, among which the no-slip condition.

The derived rates of convergence are confirmed using simple manufactured solution problems, discretized on both Cartesian and curvilinear meshes. The fact that the analysis holds for all admissible boundary conditions is also reflected in the numerical results.

This paper is organized as follows: First an introduction into differential geometry is given and the Stokes problem is reformulated in terms of differential forms. In Section 3 the mixed formulation is given and well-posedness is proven. In Section 4 the key properties of the mimetic discretization are explained that lead to compatible function spaces. This includes a discussion on the relevant properties of algebraic topology, the definitions op mimetic operators, the introduction of mimetic spectral element basisfunctions and finally the proof of discrete well-posedness. Having formulated the conforming/compatible finite element spaces, the error estimates are developed in Section 5 and the numerical results are shown in Section 6.

2. Notation and preliminaries

2.1. Differential forms

Differential forms offer significant benefits in the construction of structure-preserving spatial discretizations. For example, the coordinate-free action of the exterior derivative and generalized Stokes theorem give rise to commuting properties with respect to mappings between different manifolds. Acknowledging and respecting these kind of commuting properties are essential for the structure preserving behavior of the mimetic method.

Only those concepts from differential geometry which play a role in the remainder of this paper will be explained. More can be found in [2, 26, 27, 41].

Let denote a space of differential -forms or -forms, on a sufficiently smooth bounded -dimensional oriented manifold with boundary . Every element has a unique representation of the form


where with and where is a continuously differentiable scalar function, . Differential -forms are naturally integrated over -dimensional manifolds, i.e. for and , with ,


where indicates a duality pairing between the differential form and the geometry. Note that the -dimensional computational domain is indicated as , so without subscript. The differential forms live on manifolds and transform under the action of mappings. Let be a mapping between two manifolds. Then we can define the pullback operator, , expressing the -form on the -dimensional reference manifold, . The mapping, , and the pullback, , are each others formal adjoints with respect to a duality pairing (2.2),


where is an -dimensional submanifold of and a -dimensional submanifold of . A special case of the pullback operator is the trace operator. The trace of -forms to the boundary, , is the pullback of the inclusion of the boundary of a manifold, , see [41].

The wedge product, , of two differential forms and is a mapping: . The wedge product is a skew-symmetric operator, i.e. .

An important operator in differential geometry is the exterior derivative, . It is induced by the generalized Stokes’ theorem, combining the classical Newton-Leibnitz, Stokes circulation and Gauss divergence theorems. Let be a -dimensional manifold and , then


where is a -dimensional manifold being the boundary of . The duality pairing in (2.4) shows that the exterior derivative is the formal adjoint of the boundary operator . The exterior derivative is independent of any metric and coordinate system. Applying the exterior derivative twice always leads to the null -form, for all . As a consequence, on contractible domains the exterior derivative gives rise to an exact sequence, called De Rham complex [27], and indicated by ,


In vector calculus a similar sequence exists, where, from left to right for , the ’s denote the vector operators grad, curl and div. The exterior derivative and wedge product are related according to Leibnitz’s rule as: for all and ,


The pullback operator and exterior derivative possess the following commuting property,


In this paper we will consider Hilbert spaces , where in (2.1) the functions . The pointwise inner-product of -forms, , is constructed using inner products of one-forms, that is based on the inner product on vector spaces, see [26, 27]. The wedge product and inner product induce the Hodge- operator, , a metric operator that includes orientation. Let , then


where is a unit volume form, . In geometric physics the Hodge- switches between an inner-oriented description of physical variables and an outer-oriented description. See [40, 41, 43, 52] for a thorough discussion on the concepts of inner and outer orientation. The space of square integrable -forms on can be equipped with a inner product, , given by,


The norm corresponding to the space is . Higher degree Sobolev spaces, , consists of all -forms as in (2.1) where , with corresponding norms and . The Hilbert space associated to the exterior derivative is defined as


and the norm corresponding to is defined as . The -semi-norm is the -norm of the exterior derivative, . Note that , where the left equality holds for and the right for . The -de Rham complex, also called Hilbert complex [16], , is the exact sequence of maps and spaces given by


In terms of vector operations the Hilbert complex becomes for ,

and for , either

The two are related by the Hodge- operator (2.8), see [46],

Remark 1.

The upper complex is associated with outer-oriented -forms, i.e. -forms that are associated with outer-oriented manifolds, and the lower complex is associated with inner-oriented -forms. In this paper we mainly consider the upper complex and circumvent the lower complex by means of integration by parts. Only the pressure and tangential velocity boundary conditions are given on the lower complex, as we will see in the following sections.

A similar double Hilbert complex can be constructed in . Since the exterior derivative is nilpotent, it ensures that the range, , of the exterior derivative on -forms is contained in the nullspace, , of the exterior derivative on -forms, .

Every space of -forms in the complex can be decomposed into the nullspace of , and its orthogonal complement, . This is the Hodge decomposition, where on contractible domains . By the Hodge decomposition it follows that the exterior derivative is an isomorphism .

The inner product gives rise to the formal Hilbert adjoint of the exterior derivative, the codifferential operator, . Let , then


In case of non-zero trace, and by combining (2.9), (2.4) and (2.6), we obtain integration by parts,


Also the codifferential operator is nilpotent, , i.e., its range is contained in its nullspace, , where and . In fact the codifferential is an isomorphism , where follows from the following Hodge decomposition, . On contractible manifolds this gives rise to the following exact sequence,


In vector notation from right to left the ’s denote the -grad, curl and -div operators in , as were also mentioned in the introduction. However, whereas the exterior derivative is a metric-free operator, the codifferential operator is metric-dependent. The Hodge-Laplace operator, , is constructed as a composition of the exterior derivative and the codifferential operator,


An important inequality in stability analysis, relating the -norm and the -norm, is Poincaré inequality.

Lemma 1 (Poincaré inequality).

[6] Consider the Hilbert complex , then the exterior derivative is a bounded bijection from to , and hence, by Banach’s bounded inverse theorem, there exists a constant such that


Finally, for Hilbert spaces with essential boundary conditions we write, , and for natural boundary conditions we consider the following trace map, .

2.2. Stokes problem in differential form notation

Consider again a bounded contractible domain . Because we require exact conservation of mass and because we can perform exact discretization of the exterior derivative, see Section 4.2, we use the following formulation for the Stokes problem: let , then the VVP formulation is given by


In the VVP formulation the pressure in (2.18b) acts as a Lagrange multiplier for the constraint on velocity, (2.18c), whereas velocity in (2.18a) acts as a Lagrange multiplier for the constraint on vorticity in (2.18b).

Let be the boundary of , where

We will impose the tangential vorticity and normal velocity as essential boundary conditions, and the tangential velocity and the pressure plus divergence of velocity as the natural boundary conditions:


Then the boundary can be partitioned into four sections, , with for , where


This decomposition, introduced before in [23, 36, 40], shows all admissible boundary conditions. It will also follow directly from the mixed formulation, see (3.7), Section 3.

In case, , , no pressure boundary conditions are prescribed, and so the pressure is only determined up to an element , i.e. up to a constant. As a post processing step either the pressure in a point in can be set, or a zero average pressure can be imposed; i.e. . In case , no velocity boundary conditions are prescribed, and so the solution of velocity is determined modulo a curl-free element, i.e. modulo .

3. Mixed formulation

3.1. Mixed formulation of Stokes problem

The use of a mixed formulation is based on the following reasoning; We know how to discretize exactly the metric-free exterior derivative , but it is less obvious how to treat the codifferential operator .

3.1.1. Generalized Poisson problem

Take for example the generalized Poisson problem using the Hodge-Laplacian acting on -forms, , on with boundary . A standard Galerkin approach, using integration by parts (2.14), would give; find with and , given , such that


It has a corresponding minimization problem for an energy functional over the space . The standard Galerkin formulation is coercive, which immediately implies stability. Corresponding to this standard Galerkin formulation one usually chooses a -conforming approximation space. This could be a standard continuous piecewise polynomial vector space based on nodal interpolation.

However, in case of a nonconvex polyhedral or curvilinear or noncontractible domain , for allmost all , . Consequently, the solution will be stable but inconsistent in general, [20]. In other words, the solution converges to the wrong solution. Unfortunately, it seems not possible to construct conforming finite element spaces. Alternatively, one proposed to use mixed formulations, [15]. In contrast to standard Galerkin, the mixed formulation uses integration by parts (2.14) to express each codifferential in terms of an exterior derivative and suitable boundary conditions.

Consequently, mixed formulations require only -conforming finite element spaces, which are much easier to construct. Therefore, in all cases mixed formulations do converge to the true solution. Mixed formulations correspond to saddle point problems instead of minimization problems.

The derivation of the mixed formulation of the Poisson problem consists of three steps:

  1. Introduce an auxiliary variable in ,

  2. multiply both equations by test functions using -inner products,

  3. use integration by parts, as in (2.14), to express the remaining codifferentials in terms of the exterior derivatives and boundary integrals.

Again the boundary may constitute up to four different types of boundary conditions,


Then also for the generalized Poisson problem the boundary consists up to four sections as defined in (2.20). To obtain a unique solution for the corresponding mized formulation, we define the following Hilbert spaces,


with corresponding norms, , respectively. The resulting mixed formulation for the Poisson problem for all becomes: find , given , for all , such that


Note that, for a scalar Poisson, it is not a choice whether to use Galerkin or mixed formulation, but it depends on whether the scalar is a 0-form or an -form. This is determined by the physics.

3.1.2. Stokes problem

In a similar way the mixed formulation of the VVP formulation of the Stokes problem is obtained. Consider the Hilbert spaces and defined in the previous section, where , and define the following Hilbert space


with corresponding norm and where . Then the mixed formulation of the VVP formulation reads: find , for the given data , and natural boundary conditions , , for all , such that


Again use is made of integration by parts, (2.14).

Proposition 1.

[7] Problems (2.18)-(2.19) and (3.7) are equivalent, in the sense that any triple is a solution of problem (2.18)-(2.19) if and only if it is a solution of problem (3.7).

3.2. Well-posedness of mixed formulation

Before we continue first define the following nullspaces of ,

and consider the following decompositions, and . Since vorticity is defined as , we have , and because we consider contractible domains only, it follows that . Note that for , . A similar decomposition can be made for . Define

then . The velocity is decomposed as , where and .

We can write the mixed formulation of (3.7) in a more general representation, using four continuous bilinear forms,

and three continuous linear forms

The mixed formulation becomes


There exists continuity constants such that


By Cauchy-Schwarz we know that , however we write for generality purpose. The continuous linear forms are bounded such that


At first restrict to all . This gives the vorticity-velocity subproblem, which is a saddle point problem:

Proposition 2.

[23] System (3.12) has a unique solution if there exists positive constants , such that we have coercivity in the kernel of ,


and satisfies the following inf-sup condition for ,


The proof of (3.13) is straightforward, see e.g. [15]. For (3.14), we have , where . Thus, given there exists a unique such that and by Lemma 1. Therefore

Proposition 3.

The full problem (3.9) has a unique solution if conditions (3.13) and (3.14) from Proposition 2 are satisfied and additionally is there exists positive constants , such that we have coercivity in the range of ,


and satisfies the following inf-sup condition for ,


The proof is similar to that of Proposition 2. See also [8], Section 7.1. ∎

So well-posedness of the Stokes problem (3.7) relies only on the Hodge decomposition and the Poincaré inequality.

Corollary 1.

[7, 23] Problem (3.7) is well-posed according to Propositions 2 and 3. That is, for any given data and and natural boundary conditions and , there exists a unique solution satisfying (3.7). Moreover, this solution satisfies:


where is a constant depending only on the Poincaré constant and the continuity constants.

4. Compatible spectral discretization

Well-posedness of the Stokes problem in VVP formulation relies solely on the Hodge decomposition and the Poincaré inequality. For a compatible discretization, these properties need to be respected as well in the finite dimensional spaces. Key ingredient to obtain a discrete Hodge decomposition and discrete Poincaré inequality is the construction of a bounded projection operator that commutes with the exterior derivative.

The compatible spectral discretization consists of three parts. First, the discrete structure is described in terms of chains and cochains from algebraic topology, the discrete counterpart of differential geometry. This discrete structure mimics many of the properties from differential geometry. Secondly, mimetic operators are introduced that relate the continuous formulation in terms of differential forms to the discrete representation based on cochains and finite dimensional differential forms. Thirdly, mimetic spectral element basis functions are described following the definitions of the mimetic operators. In this paper we address these topics only briefly. More details of the mimetic spectral element method can be found in [40, 41]. Finally, well-posedness of the discrete numerical formulation is proven and interpolation error estimates are given.

4.1. Algebraic Topology

Let be an oriented cell-complex covering the manifold , describing the topology of the mesh, and consisting of -cells , . The two most popular classes of -cells in literature to describe the topology of a manifold are either in terms of simplices, see for instance [44, 51, 53], or in terms of cubes, see [42, 52]. From a topological point of view both descriptions are equivalent, see [22]. Despite this equivalence of simplicial complexes and cubical complexes, the reconstruction maps in terms of basis functions, to be discussed in Section 4.2, differ significantly. For mimetic methods based on simplices see [5, 6, 21, 47], whereas for mimetic methods based on singular cubes see [3, 37, 39, 50]. We restrict ourselves to -cubes, although we will keep calling them -cells.

The ordered collection of all -cells in generate a basis for the space of -chains, . Then a -chain, , is a formal linear combination of -cells, ,


The boundary operator on -chains, , is an homomorphism defined by [33, 44],


The boundary of a -cell will then be a -chain formed by the oriented faces of . Like the exterior derivative, applying the boundary operator twice on a -chain gives the null -chain,