A posteriori error analysis of an enriched Galerkin method of order one for the Stokes problem

A Approximation of the boundary data

In practical situations, the boundary data g_D is given by measurements or by the results of a previous computation; in both cases, it is a discrete set of values. To simplify the a posteriori error analysis, it is convenient to assume that this discrete set is interpolated so as to be the trace of a function of M0kThd on Γ, say Π_h(g_D), ^0kand the condition (2.1) is enforced by adding a correction of the form

This is adequate since the function x is a polynomial of degree one, thus is in the discrete space, and its divergence is the constant d, the dimension; clearly,

Then the discrete set g_D is approximated by

By construction, g_h,D satisfies (2.1), and the discrete problem (3.5) is solved with the data g_h,D. To simplify notation, the index h in g_D was dropped in the previous analysis.

The same approach is used in the case of a manufactured solution, but in this case we can measure the error made when replacing g_D with g_h,D. For Π_h, we use the Scott–Zhang interpolation operator of degree k in each element, globally continuous. First,

Now,

which is a bounded quantity and of course, ‖x‖_[L²_(Γ)]d is also a bounded quantity. Therefore the approximation error of g_D is bounded as follows:

and

The a posteriori error analysis developed in the previous sections holds with g_h,D instead of g_D and (u, p) replaced by (u(h), p(h)) ∈ H¹(Ω)^d × L(Ω), solution of

As g_h,D is a continuous piecewise polynomial, it belongs to [H¹(Γ)]^d and so there exists some s, 0 < s < 1/2, such that u(h) ∈ [H^1+s(Ω)]^d and p(h) ∈ H^s(Ω). Hence the error equalities of Section 4.1 are well-defined and all reliability and efficiency bounds derived above hold with u and p replaced by u(h) and p(h). To express the reliability bounds in terms of u and p, we use the following:

Regarding efficiency, we must add to the right-hand side of (4.51) the terms

and to that of (4.56),

These quantities can be bounded by the global estimates in (A.4), but this is not optimal. Sharper bounds can be derived by applying the localizing techniques of 34, but this is outside the scope of this work.

B Choice of penalty parameters

The penalty parameters σ_e are chosen to guarantee that the form a_ε is elliptic in Vh,kEGd. In the symmetric or incomplete case, we consider an interior edge or face e (the case of a boundary e is simpler) and evaluate (∇ u_h ⋅ n_e , [v_h])_e, where u_h denotes the value of u_h restricted to one of the elements E adjacent to e,

Thus,we must estimate ∇uhEL2(e)d×d. In the case of general k, or in the case of quadrilaterals or hexahedra, we use the fact that, after transformation, the function ∇ u_h belongs to a finite dimensional space in the reference element and apply an equivalence of norms to switch from the surface norm to the volume norm; this brings a constant, say ̂D, that depends only on d and the degree of the polynomials in the reference element, hence is independent of h. But in the simplex case when k = 1, this is trivial because ∇ u_h is a constant tensor. To simplify, we study this case here but keep in mind that there is an additional constant in the quadrilateral or hexahedral cases. We have

Therefore, when E₁ and E₂ are the two elements sharing e, we have

Now, we sum over all interior e:

Let us take the maximum of the factor

and observe that the second factor is of the order of h⁻¹ , whatever the dimension. Therefore, _e

where depends only on the regularity of the mesh. It remains to count the number of occurrences of an element E in the sum over e. Each face involves its two adjacent elements. Since a simplex has d + 1 faces, it is counted d + 1 times. Therefore,

In the case of a quadrilateral or hexahedral element, d + 1 is replaced by 2d, the number of faces, and the constantis multiplied by ̂D. Of course, this formula is also valid for boundary faces, except that the factor 2(d+1)/2 is replaced by 1 and the sums are taken over εh∂. It is used with Young’s inequality to deduce the value of min σ_e that guarantees ellipticity of a_ε when ε = 0 or ε = −1. When ε = 1 (non-symmetric case) strictly speaking σ_e could be zero but experiments have shown that the choice σ_e = 1 brings more stability.

C The case of hexahedra with curved faces

From the computer implementation point of view, handling curved faces is more difficult, in particular because the normal vector to each face n_e is now a function of its position x. But when the family of meshes is regular, this function does not vary much and from a theoretical point of view, all the preceding work applies immediately to regular families of hexahedra with straight edges but non planar faces, with the possible exception of the inf-sup condition. To establish the inf-sup condition in this case, we proceed as in 10, namely use Fortin’s lemma. This amounts to constructing a suitable approximation operator Ph∈LH01(Ω)d;Vh,k,0EG, satisfying for all v∈H01(Ω)d,

In 10, such an operator is constructed by correcting a Scott–Zhang approximation operator, say Π_h ∈ LH01(Ω)d;MkTh∩H01(Ω)d by means of piecewise constants,

where c(v) is a constant vector in each cell. It is shown in 10 that a sufficient condition for the equality in (C.1) is that c(v) satisfy on all faces γ_i of E with interior normal , n_γ

In the case when γ_i is a plane face, (C.2) reads

The regularity of the element E implies that this system is solvable provided the maximum number of such faces per element is no more than three in 3D or two in 2D and none are opposite faces, in other words, these faces are part of a conical angle. The situation is the same when γ_i is a curved face and the mesh is regular, because the normal vector to a face varies little. In this case (C.3) is replaced by

Again, we assume that an element has no more than d faces with interior normals and none are opposite faces. Take the 3D case of exactly three such faces, say γ₁, γ₂, γ₃. When written in matrix form, (C.4) reads

where b is the vector with ith component

and the ith line of the matrix N is

Of course, the fact that n_γi are unit normals imply that N is uniformly bounded. Moreover, the regularity assumption on the family of triangulations implies that N is nonsingular and its inverse is uniformly bounded. The second part of (C.1) follows from this.