A simple nonconforming tetrahedral element for the Stokes equations

In this paper we apply a nonconforming rotated bilinear tetrahedral element to the Stokes problem in $\mathbb{R}^3$. We show that the element is stable in combination with a piecewise linear, continuous, approximation of the pressure. This gives an approximation similar to the well known continuous $P^2-P^1$ Taylor$-$Hood element, but with fewer degrees of freedom. The element is a stable non-conforming low order element which fulfils Korn's inequality, leading to stability also in the case where the Stokes equations are written on stress form for use in the case of free surface flow.


INTRODUCTION
The nonconforming rotated Q 1 tetrahedron is derived from the nonconforming hexahedral element proposed by Rannacher and Turek [6] and was applied to linear and nonlinear elasticity problems in [3,4,5].It has properties similar to the hexahedral element, with improved bending behaviour in elasticity, compared to the P 1 tetrahedron, and it allows for diagonal mass lumping for explicit time-stepping in dynamic problems.
In this paper, we investigate its properties as a Stokes element, and show that in combination with linear, continuous pressures, it is inf-sup stable.In a sense it is thus a reduced Taylor-Hood [7] element with fewer degrees of freedom.It is one of the lowest order elements that is stable for Stokes, and, unlike some other low order non-conforming elements [2,6], it fulfills Korn's inequality and can thus handle the strain form of Stokes, cf.[3].
An outline of the paper is as follows: in Section 2 we recall the rotated Q 1 element; in Section 3 we apply it to the Stokes equations and prove stability and convergence; in Section 4 we give some numerical examples to show the properties of the approximation.

THE ROTATED Q 1 APPROXIMATION
In order to define a low order approximation which is stable for the Stokes problem, Rannacher and Turek [6] constructed a hexahedral element with nodes on the faces.Two different kinds of continuity can now be imposed: point-wise continuity and average continuity and we will in this paper consider point-wise continuity.The Rannacher-Turek element may be viewed as an extension of the classical Crouzeix-Raviart nonconforming tetrahedral element [2], for which point-wise and average continuity is identical.
The difficulty in the construction of the Rannacher-Turek element compared to the Crouzeix-Raviart element is that we have six degrees of freedom while the dimension of the space of linears is four, trilinears eight, and quadratics ten.Therefore we shall start with linears and add two suitable quadratic functions in such a way that the nodal mapping is invertible.To that end let {e i } 3 i=1 be an orthonormal coordinate system in R 3 and let K = [−1, 1] 3 be the reference cube.The midpoints of the faces of K takes the form {x i } 6 i=1 = {±e i } 3 i=1 , and therefore has one coordinate equal to ±1 and the other two equal to 0. See Fig. 1 for the enumeration of the midpoints.We now seek a space V( K) of shape functions on K such that P 1 ( K) ⊂ V( K) ⊂ P 2 ( K) and that has the values at the six midpoints {x i } 6 i=1 of the faces as degrees of freedom.The linear functions P 1 ( K) = span(1, x 1 , x 2 , x 3 ) is a four dimensional vector space and therefore we need to add two quadratic polynomials to obtain a six dimensional space.To find these quadratic polynomials let N : i=1 ∈ R 6 be the nodal mapping and note that N(x i x j ) = 0, for i j, and therefore we are restricted to adding polynomials in span(x 2  1 , x 2 2 , x 2 3 ).Now N(x 2 1 + x 2 2 + x 2 3 ) = N(1) and it is therefore natural to consider the two dimensional space where the two (non unique) basis functions on the right hand side is easily chosen by observing that and verify by explicit calculation that the coordinate mapping N : where The terminology rotated Q 1 elements is motivated by the fact that = ξ 1 ξ 2 where ξ 1 = x 1 − x 2 and ξ 1 + ξ 2 are the degrees of freedom in a coordinate system rotated π/2 around the e 3 axis.
In [3], we made the observation that there is a reference tetrahedron T inscribed in the reference hexahedron K, with edges that are diagonals of the faces of K, and thus the midpoints of the faces of K are precisely the midpoints of the edges of T , see Fig. 1.We note that for the reference element all edges have the same length and the centre of gravity is the origin.Using an affine map F : T → T , we map mid side nodes in the reference configuration to mid edge nodes in the physical configuration.
Let T h := {T } be a conforming, shape regular tetrahedrization of Ω ⊂ R 3 with mesh parameter h ∈ (0, h 0 ].We also let F h be the set of faces and E h the set of edges.We make the standard assumption for the Taylor-Hood approximation [1], that every T ∈ T h has at least three internal edges (6) We define the non-conforming finite element space with midpoint continuity for all interior edges.We note that dim(V h ) = |E h |, the number of edges in T h .

Problem Formulation and Finite Element Approximation
We consider the Stokes equations in a domain Ω in R 3 : find the velocity u = [u i ] 3 i=1 and the pressure p such that Let us define the spaces and where (•, •) Ω is the standard L 2 scalar product.Then we have the weak form of (8 where the forms are given by To define the finite element method, we introduce the non-conforming finite element space constructed from the space V h in (7) by defining and the space of continuous piecewise linear polynomials The finite element method is to find (u Here the form is defined by with a h (v, w) := and ∇v = v ⊗ ∇ is the tensor with elements (v ⊗ ∇) i j = ∂ j v i .
Remark 3.1.We will also consider the alternative form This form is preferable in the presence of natural boundary conditions, and requires less numerical computations, but we will see that the proof of the inf-sup condition is more complicated due to the presence of the trace term on ∂T in the right hand side of (18).Throughout the paper we will focus our presentation on the form b and we will for each result add a remark on the modifications necessary to obtain the corresponding result for b.Finally, we will present an inf-sup result for b in Appendix A.
Remark 3.2.Unlike some nonconforming approximations, the rotated Q 1 approximation fulfills Korn's inequality [3], which means that we may also use the strain form of Stokes: find the velocity u and the pressure p such that where is the strain tensor and ∇• denotes matrix divergence.This is of interest in free surface flows where we need zero stress as a natural boundary condition, cf.Section 4.2.

Norms and Continuity of the Forms
Define the norms and We let a b denote a ≤ Cb with C a positive constant independent of the mesh parameter.Then we have the following continuities of the forms b and A h .

Lemma 3.1.
There are constants such that for all functions in for an interior face shared by elements T 1 and T 2 and [n • v] = n • v for a face at the boundary belonging to element T .Noting that for a face F we then have [n • v(x E )] = 0 for the midpoint x E of the edge E. Therefore ([n • v], 1) F = 0, since the quadrature formula based on the midpoints of the edges on a triangle is exact for quadratic polynomials.We therefore have ( , where P 0,F is the L 2 projection on constants on the face F. With these preparations at hand we obtain the following bound which proves (23).Finally, (24) follows directly from (23) and the Cauchy-Schwarz inequality.
Remark 3.3.For b defined in (18) we directly have since there are no trace terms on the boundary of the elements.

Interpolation
We shall now define interpolants for the finite element space.Starting with the pressure space we let π h,p : L 2 (Ω) → Q h be a Clement interpolant.We then have the standard estimate To construct an interpolant for the velocity space we use component-wise Scott-Zhang interpolation to satisfy the Dirichlet boundary conditions (in the nodes), Here we also have the interpolation estimate From here on we simplify the notation and write π h,p = π h and π h,u = π h and interpret the operator in the correct way depending on in which space the argument reside.Combining the estimates we get

Stability Analysis
We first recall the following standard result from [1].
Then there is a constant such that We shall now prove that (34) holds for the nonconforming space W h × Q h using an approach called Verfürth's trick [8], which proceeds in two steps.

Lemma 3.2 (Step 1).
There are constants c 1 and c 2 such that We shall now replace v ∈ W by the interpolant π h v and estimate the remainder term as follows Here we used partial integration together with the continuity of q and the boundary condition v = 0 on ∂Ω, the interpolation estimate (32) and the boundedness ∇π h v Ω ∇v Ω of the interpolation operator.
Lemma 3.3 (Step 2).Under the assumption (6) there exists c 3 > 0 such that PROOF.Using partial integration we have the identity Observing that ∇p is element-wise constant we may apply the quadrature formula which is exact for w ∈ V h , to obtain where {ϕ E } E∈E h is the global basis in V h , E h (T ) is the set of edges belong to element T , and Since q is continuous it follows that the tangent derivative t E • ∇q along each edge E, with unit tangent vector t E , is also continuous and thus taking where E I h ⊂ E h is the set of interior edges, we get since ∇q is element-wise constant and there are, by assumption, three linearly independent tangent vectors in the set {t E } E∈E h (T ) for each element T ∈ T h .Finally, noting that where we used an inverse bound to conclude that ∇ϕ E h −2 h 3 , and summing over T gives Combining the estimates we get the desired result since Remark 3.4.For the alternative form b we get the more complicated expression where the trace term on the boundary of T does not vanish.We will however show in Appendix A that with the same choice of v * the trace term can be shown to be dominated by the bulk term on each element.The proof is based on mapping to the reference element and explicit computation of the two integrals.Thus Lemma 3.3 also holds for the form b.

Lemma 3.4.
There is a constant such that the inf-sup condition (34) holds.

Error Analysis Theorem 3.2.
There is a constant such that PROOF.We first split the error in an interpolation error part and a discrete part The first term can be estimated using the interpolation error estimate (33).To estimate the discrete part of the error we employ the inf-sup condition (35) to obtain Here we used continuity (24) of the form A h for the first term which can now be estimated using the interpolation error estimate (33).The second term accounts for the consistency error and using partial integration we find that where for a face at the boundary.Then using the fact that [v] = 0 in the midpoints of the edges and that midpoint quadrature is exact for quadratic polynomials on a triangle it follows that (1, ([v ⊗ n]) i j ) F = 0 and therefore we may subtract the L 2 -projection on constant functions on the faces and then estimate the contributions using the following standard bounds where ∇ F = (I − n F ⊗ n)∇, with n F a unit normal to the F, is the tangential gradient to the face F. Combing the bounds gives the desired estimate.
Remark 3.5.For the form b we get a consistency error of the form which, using the same approach as in estimates (63-65), can be estimated by

Convergence
We consider a problem in the ball with radius 1 and with center at the origin.A fabricated solution is given by u with f = 0.The exact solution is used as Dirichlet data and zero mean pressure is imposed by a Lagrange multiplier.We compare the convergence of the pressure inconsistent method, using b(u, q), to the pressure consistent method, using b(u, q), in Fig. 3.The convergence is shown in L 2 for the pressure and the velocity and in broken H 1 semi-norm for the velocity.We note that the methods converge at the same rate, albeit with a slightly larger error constant for the inconsistent method.The observed rates from Fig. 3 are As we do not have full P 2 approximation for the velocity, the observed rate of convergence of pressure is better than expected.In Fig. 4 we compare the pressure solutions to the interpolated pressure on the boundary, shown on the finest mesh in the sequence used for the convergence study.

Laplacian Form vs. Strain Form
In this example we show that the strain formulation of Stokes equations poses no problem as we have a Korn inequality for our approximation.We consider a Poiseuille type problem in the domain Ω = (0, 3) × (0, 1) × (0, 1/10) with boundary conditions u = (x 2 (1 − x 2 ), 0, 0) at x 1 = 0, at x 2 = 0, and at x 2 = 1, and with u 3 = 0 at x 3 = 0 and at x 3 = 1/10.In Fig. 5 we show the velocity field in the (x 1 , x 2 )-plane and we note that the strain formulation gives a stress free condition at the outflow.In Fig. 6 we show the corresponding pressure.The computations were made using the b-form of the side condition.
< l a t e x i t s h a 1 _ b a s e 6 4 = " Z d b 7 9 L + z k n i g 4 R 9 a x I V C o o B H 2 p Q = " > A A A B 6 H i c b V D L T g J B E O z F F + I L 9 e h l I j H x R H Y N Q Y 9 E L x 4 h k U c C G z I 7 9 M L I 7 O x m Z t a E E L 7 A i w e N 8 e o n e f N v H G A P C l b S S a W q O 9 1 d Q S K 4 N q 7 7 7 e Q 2 N r e 2 d / K 7 h b 3 9 g 8 O j 4 v F J S 8 e p Y t h k s X e n E f n x X l 3 P p a t O S e b O Y U / c D 5 / A I U f j M M = < / l a t e x i t > 6 < l a t e x i t s h a 1 _ b a s e 6 4 = " p u k K C 4 / a / 4 V u m d v 9 C 1 r D g x p F P r o = " > A A A B 6 H i c b V D L T g J B E O z F F + I L 9 e h l I j H x R H Y N P o 5 E L x 4 h k U c C G z I 7 9 M L I 7 O x m Z t a E E L 7 A i w e N 8 e o n e f N v H G A P C l b S S a W q O 9 1 d Q S K 4 N q 7 7 7 e T W 1 j c 2 t / L b h Z 3 d v f 2 D 4 u F R U 8 e p Y t h g s Y h V O 6 A a B Z f Y M N w I b C c K a R Q I b A W j u 5 n f e k K l e S w f z D h B P 6 I D y U P O q L F S / b J X L L l l d w 6 y S r y M l C B D r V f 8 6 v Z j l k Y o D R N U 6 4 7 n J s a f U G U 4 E z g t d F O N C W U j O s C O p Z J G q P 3 J / N A p O b N K n 4 S x s i U N m a u / J y Y 0 0 n o c B b Y z o m a o l 7 2 Z + J / X S U 1 4 4 0 + 4 T F K D k i 0 W h a k g J i a z r 0 m f K 2 R G j C 2 h T H F 7 K 2 F D q i g z N p u C D c F b f n m V N C / K 3 l W 5 U q + U q r d Z H H k 4 g V M 4 B w + u o Q r 3 U I M G M E B 4 h l d 4 c x 6 d F + f d + V i 0 5 p x s 5 h j + w P n 8 A Y O b j M I = < / l a t e x i t > 5 < l a t e x i t s h a 1 _ b a s e 6 4 = " y + W D o q t t s z d e x E L Z g 6 E J P D Y o r d 8 7 3 U B D l V h j O B M 7 e f a U w p m 9 A R 9 i y V N E Y d 5 I t D Z + T c K k M S J c q W N G S h / p 7 I a a z 1 N A 5 t Z 0 z N W K 9 6 c / E / r 5 e Z 6 C b I u U w z g 5 I t F 0 W Z I C Y h 8 6 / J k C t k R k w t o U x x e y t h Y 6 o o M z Y b 1 4 b g r 7 6 8 T t q X N f + q V q 8 2 b o s w y n A K Z 3 A B P l x D A + 6 h C S 1 g g P A C b / D u P D m v z s e y s e Q U E y f w B 8 7 n D x e d i 5 c = < / l a t e x i t > 4 < l a t e x i t s h a 1 _ b a s e 6 4 = " P v k q X e n E f n x X l 3 P p a t O S e b O Y U / c D 5 / A H 8 P j L 8 = < / l a t e x i t > 2 < l a t e x i t s h a 1 _ b a s e 6 4 = " E f b C a c I R g j F v 3 s x P W p z 9 9 K W 1 5 4 o = " > A A A B 6 X i c b V B N S 8 N A E J 3 U r 1 q / q h 6 9 L B b B U 0 m k q M e i F 4 9 V 7 A e 0 o W y 2 k 3 b p Z h N 2 N 0 I p / Q d e P C j i 1 X / k z X / j p s 1 B W x 8 M P N 6 b Y W Z e k A i u j e t + O 4 W 1 9 Y 3 N r e J 2 a W d 3 b / X e n L H z 4 r w 7 H 4 v W g p P P H M M f O J 8 / s c i M 0 g = = < / l a t e x i t > 1 < l a t e x i t s h a 1 _ b a s e 6 4 = " m s s F Z 8 h 6 g J k /D r c R + g C T i y W x V Z I = " > A A A B + n i c b V C 7 T s M w F L 3 h W c o r h Z H F o k J i q h K o g L G C h b F I 9 C G 1 U e Q 4 b m v V c S L b A V W h n 8 L C A E K s f A k b f 4 P T Z o C W I 1 k + O u d e + f g E C W d K O 8 6 3 t b K 6 t r 6 x W d o q b + / s 7 u 3 b l Y O 2 i l N J a I v E P J b d A C v K m a A t z T S n 3 U R S H A W c d o L x T e 5 3 H q h U L B b 3 e p J Q L 8 J D w Q a M Y G 0 k 3 6 7 0 g 5 i H a h K Z K 6 N T / 7 z s 2 1 W n 5 s y A l o l b k C o U a P r 2 V z + M S R p R o Q n H S v V c J 9 F e h q V m h N N p u Z 8 q m m A y x k P a M 1 T g i C o v m 0 W f o h O j h G g Q S 3 O E R j P 1 9 0 a G I 5 W n M 5 M R 1 i O 1 6 O X i f 1 4 v 1 Y M r L 2 M i S T U V Z P 7 Q I O V I x y j v A Y V M U q L 5 x B B M J D N Z E R l h i Y k2 b e U l u I t f X i b t s 5 p 7 U a v f 1 a u N 6 6 K O E h z B M Z y C C 5 f Q g F t o Q g s I P M I z v M K b 9 W S 9 W O / W x 3 x 0 x S p 2 D u E P r M 8 f H 5 2 T 6 g = = < / l a t e x i t > e 3 < l a t e x i t s h a 1 _ b a s e 6 4 = " e w N f V J R D Y n S y H 6 P 6 H j L y W k t 7 U + I = " > A A A B + n i c b V C 9 T s M w G P x S / k r 5 S 2 F k s

2 <
t 4 G p 8 0 A L S d Z P t 1 9 n 3 y + I O F M a c f 5 t k p r 6 x u b W + X t y s 7 u 3 v 6 B X T 3 s qj i V h H Z I z G N 5 H 2 B F O R O 0 o 5 n m 9 D 6 R F E c B p 7 1 g c p 3 7 v Q c q F Y v F n Z 4 m 1 I v w S L A h I 1 g b y b e r g y D m o Z p G 5 s r o z G 9 U f L v m 1 J 0 5 0 C p x C 1 K D A m 3 f / h q E M U k j K j T h W K m + 6 y T a y 7 D U j H A 6 q w x S R R N M J n h E + 4 Y K H F H l Z f P o M 3 R q l B A N Y 2 m O 0 G i u / t 7 I c K T y d G Y y w n q s l r 1 c / M / r p 3 p 4 6 W V M J K m m g i w e G q Y c 6 R j l P a C Q S U o 0 n x q C i W Q m K y J j L D H R p q 2 8B H f 5 y 6 u k 2 6 i 7 5 / X m b b P W u i r q K M M x n M A Z u H A B L b i B N n S A w C M 8 w y u 8 W U / W i / V u f S x G S 1 a x c w R / Y H 3 + A B 4 Y k + k = < / l a t e x i t > e l a t e x i t s h a 1 _ b a s e 6 4 = " A a S E e f E q U J L 2 H t G E a J i 7 q 3 z G B r U = " > A A A B + n i c b V C 7 T s M w F L 0 p r 1 J e K Y w s F h U S U 5 W g C h g r W B i L R B 9 S G 0 W O 4 7 Z W H S e y H V A V + i k s D C D E y p e w 8 T c 4 b Q Z o O Z L l o 3 P u l Y 9 P k H C m t O N 8 W 6 W 1 9 Y 3 N r f J 2 Z W d 3 b / / A r h 5 2 2 s p L c J e / v E o 6 5 3 X 3 o t 6 4 a 9 S a 1 0 U d Z T i G E z g D F y 6 h C b f Q g j Y Q e I R n e I U 3 6 8 l 6 s d 6 t j 8 V o y S p 2 j u A P r M 8 f H J O T 6 A = = < / l a t e x i t > e 1

Figure 1 :
Figure 1: The reference element and the enumeration of the six degrees of freedom.
q 4 7 a o o w w n c A r n 4 M E 1 N O A B m t A C D G N 4 h l d 4 s 3 L r x X q 3 P h a j J a v Y O Y Y / s D 5 / A L P j k 7 Y = < / l a t e x i t >

Figure 3 :
Figure 3: Convergence of velocity and pressure for the pressure inconsistent b-form and the consistent b-form.

Figure 5 :
Figure 5: Velocity for the Laplacian form (left) and strain form (right).

Figure 6 :
Figure 6: Pressure for the Laplacian form (left) and strain form (right).