A mixed finite element method on polytopal mesh
AA MIXED FINITE ELEMENT METHOD ON POLYTOPAL MESH
YANPING LIN ∗ , XIU YE † , AND
SHANGYOU ZHANG ‡ Abstract.
In this paper, we introduce new stable mixed finite elements of any order on polytopalmesh for solving second order elliptic problem. We establish optimal order error estimates for velocityand super convergence for pressure. Numerical experiments are conducted for our mixed elements ofdifferent orders on 2D and 3D spaces that confirm the theory.
Key words.
Mixed finite element methods, second order elliptic problem
AMS subject classifications.
Primary, 65N15, 65N30; Secondary, 35B45, 35J50
1. Introduction.
The considered model problem seeks a flux function q = q ( x )and a scalar function u = u ( x ) defined in an open bounded polygonal or polyhedraldomain Ω ⊂ R d ( d = 2 ,
3) satisfying a q + ∇ u = 0 in Ω , (1.1) ∇ · q = f in Ω , (1.2) u = − g on ∂ Ω , (1.3)where a is a symmetric, uniformly positive definite matrix on the domain Ω. A weakformulation for (1.1)-(1.3) seeks q ∈ H (div , Ω) and u ∈ L (Ω) such that( α q , v ) − ( ∇ · v , u ) = (cid:104) g v · n (cid:105) ∂ Ω ∀ v ∈ H (div , Ω) , (1.4) ( ∇ · q , w ) = ( f, w ) ∀ w ∈ L (Ω) . (1.5)Here L (Ω) is the standard space of square integrable functions on Ω, ∇ · v is thedivergence of vector-valued functions v on Ω, H (div , Ω) is the Sobolev space consistingof vector-valued functions v such that v ∈ [ L (Ω)] d and ∇ · v ∈ L (Ω), ( · , · ) standsfor the L -inner product in L (Ω), and (cid:104)· , ·(cid:105) ∂ Ω is the inner product in L ( ∂ Ω).Finite element methods based on the weak formulation (1.4)-(1.5) and finite di-mensional subspaces of H (div , Ω) × L (Ω) with piecewise polynomials are known as ∗ Department of Applied Mathematics, Hong Kong Polytechnic University, Hung Hom, HongKong, China ([email protected]). † Department of Mathematics, University of Arkansas at Little Rock, Little Rock, AR 72204([email protected]). This research was supported in part by National Science Foundation Grant DMS-1620016. ‡ Department of Mathematical Sciences, University of Delaware, Newark, DE 19716([email protected]). 1 a r X i v : . [ m a t h . NA ] S e p mixed finite element methods (MFEM). The mixed finite element methods have beenintensively studied [1, 2, 3, 4, 5, 7, 8] and many stable mixed finite elements have beendeveloped such as Raviart-Thomas (RT) and Brezzi-Douglas-Marini (BDM) elements.However most of the existing mixed elements are defined on triangle/rectangle in twodimensional space and tetrahedron/cuboid in three dimensional space.Construction of stable mixed finite elements on general polytopal mesh can bevery challenging. Recently, a lowest order mixed element on polytopal mesh wasintroduced in [6] by using rational Wachspress coordinates. The goal of this paper is toconstruct stable mixed elements of any order on polytopal mesh. Optimal convergencerate for velocity and superconvergence for pressure are obtained. Extensive numericalexamples are tested for the new mixed finite elements of different degrees in two andthree dimensions.
2. Construction of a H (div , Ω) Element.
Let T h be a partition of the domainΩ consisting of polygons in two dimension or polyhedra in three dimension satisfyinga set of conditions specified in [9]. Denote by E h the set of all edges/faces in T h , andlet E h = E h \ ∂ Ω be the set of all interior edges/faces. For simplicity, we will use termedge for edge/face without confusion. Let P k ( K ) consist all the polynomials degreeless or equal to k defined on T .The space H (div; Ω) is defined as the set of vector-valued functions on Ω which,together with their divergence, are square integrable; i.e., H (div; Ω) = (cid:8) v ∈ [ L (Ω)] d : ∇ · v ∈ L (Ω) (cid:9) . For any T ∈ T h , we divide it in to a set of disjoint triangles/tetrahedra T i with T = ∪ T i . We define Λ h ( T ) asΛ k ( T ) = { v ∈ H (div; T ) : v | T i ∈ RT k ( T i ) , ∇ · v ∈ P k ( T ) } , (2.1)where RT k ( T i ) = { [ P k ( T i )] d ⊕ x (cid:80) | α | = k a α x α } is the usual Raviart-Thomas elementof order k .Associated with the given mesh, we introduce two finite element spaces(2.2) V h = { v ∈ H (div; Ω) : v | T ∈ Λ k ( T ) , T ∈ T h } , and(2.3) W h = { w ∈ L (Ω) : w | T ∈ P k ( T ) , T ∈ T h } . Lemma 2.1.
For the projection Π h defined in (2.7) below and for τ ∈ H (div; Ω) and v ∈ P k ( T ) , we have ( ∇ · τ, v ) T = ( ∇ · Π h τ, v ) T , (2.4) (cid:107) Π h τ − τ (cid:107) ≤ Ch k +1 | τ | k +1 . (2.5) Proof . We assume no additional inner vertex/edges is introduced in subdividing apolygon/polyhedron T in to n triangles/tetrahedrons { T i } . That is, we have precisely n − T into n parts. We limit the proof to3D. We need only omit the fourth equation in (2.7) to get a 2D proof.On n tetrahedrons, a function of Λ k can be expressed as v h | T i = (cid:88) i + j + l ≤ k a ,ijl a ,ijl a ,ijl x i y j z l + (cid:88) i + j + l = k xyz a ,ijl x i y j z l , i = 1 , ...n. (2.6) v h | T i is determined by n ( k + 1)( k + 2)( k + 3)2 + n ( k + 1)( k + 2)2 = n ( k + 1)( k + 2)( k + 4)2coefficients. For any v ∈ H (div; T ), Π h v ∈ Λ k ( T ) is defined by (cid:90) F ij ⊂ ∂T (Π h v − v ) · n ij p k dS = 0 ∀ p k ∈ P k ( F ij ) , (cid:90) T (Π h v − v ) · n p k − d x = 0 ∀ p k − ∈ P k − ( T ) , (cid:90) T i (Π h v − v ) · n p k − d x = 0 ∀ p k − ∈ P k − ( T i ) , i = 1 , ...n, (cid:90) T i (Π h v − v ) · n p k − d x = 0 ∀ p k − ∈ P k − ( T i ) , i = 1 , ...n, (cid:90) F ij ⊂ T [Π h v ] · n ij p k dS = 0 ∀ p k ∈ P k ( F ij ) , (cid:90) T ∇ · (Π h v | T i − Π h v | T ) p k d x = 0 ∀ p k ∈ P k ( T ) , i = 2 , ..., n, (2.7)where F ij is the j -th face triangle of T i with a fixed normal vector n ij , n is aunit vector not parallel to any internal face normal n ij , ( n , n , n ) forms a right-hand orthonormal system, [ · ] denotes the jump on a face triangle, and Π h v | T i isunderstood as a polynomial vector which can be used on another tetrahedron T .The linear system (2.7) of equations has the following number of equations,(2 n + 2) ( k + 1)( k + 2)2 + (2 n + 1) k ( k + 1)( k + 2)6+ ( n −
1) ( k + 1)( k + 2)2 + ( n −
1) ( k + 1)( k + 2)( k + 3)6= n ( k + 1)( k + 2)( k + 4)2 , which is exactly the number of coefficients for a v h function in (2.6). Thus we havea square linear system. The system has a unique solution if and only if the kernel is { } .Let v = 0 in (2.7). Though Π h v is a P k +1 polynomial, Π h v · n ij is a P k polynomialwhen restricted on F ij . This can be seen by the normal format of plane equation fortriangle F ij . By the first equation of (2.7), Π h v · n ij = 0 on F ij . By the sixthequation of (2.7), ∇ · Π h v is a one-piece polynomial on the whole T . Because ∇ · Π h v is continuous on inner interface triangles and is a P k ( F ij ) polynomial on the outerface triangles, by the first five equations in (2.7), we have (cid:90) T ( ∇ · Π h v ) d x = n (cid:88) i =1 (cid:16) (cid:90) T i − Π h v · ∇ ( ∇ · Π h v ) d x + (cid:90) ∂T i Π h v · n ( ∇ · Π h v ) dS (cid:17) = n (cid:88) i =1 3 (cid:88) j =1 (cid:90) T i − (Π h v · n j )( n j · ∇ ( ∇ · Π h v )) d x = 0 . That is, ∇ · Π h v = 0 on T. (2.8)Starting from a corner tetrahedron T , we have its three face triangles, F , F and F , on the boundary of T . The forth face triangle F of T is shared by T . Bythe selection of n , the normal vector n = c n + c n + c n of F has a non zero c (cid:54) = 0. a 2D polynomial p k ∈ P k ( F ) can be expressed as p k ( x , x ), where we use( x , x , x ) as the coordinate variables under the system ( n , n , n ). Viewing thispolynomial as a 3D polynomial, i.e. extending it constantly in x -direction, we have p k ( x , x , x ) = p k ( x , x ) , ( x , x , x ) ∈ T . By (2.8) and the third and fourth equations of (2.7), it follows that0 = (cid:90) T ( ∇ · Π h v ) p k d x = − (cid:90) T (cid:16) (Π h v · n ) ∂ x p k + (Π h v · n ) ∂ x p k + (Π h v · n ) ∂ x p k (cid:17) d x + (cid:90) F (Π h v ) · n p k dS = − (cid:90) T (Π h v · n ) · d x + 0 + 0 + (cid:90) F (Π h v ) · n p k dS = (cid:90) F (Π h v ) · n p k dS ∀ p k ∈ P k ( F ) . (2.9)Next, for any p k − ∈ P k − ( T ), we let p k ∈ P k ( T ) be one of its anti- x -derivative,i.e., ∂ x p k = p k − . Thus, by (2.8), the third and fourth equations of (2.7) and (2.9),we get 0 = (cid:90) T ∇ · Π h v p k d x = − (cid:90) T (cid:16) (Π h v · n ) ∂ x p k + 0 + 0 (cid:17) d x + (cid:90) F (Π h v ) · n p k dS = − (cid:90) T (Π h v · n ) p k − d x ∀ p k − ∈ P k − ( T ) . (2.10)Continuing work on T , by ∇· Π h v = 0, all a ,ijl = 0 in (2.6), since the divergenceof each such term is non-zero and independent of the divergence of other terms.Thus Π h v | T i is in [ P k ( T i )] d , instead of RT k ( T i ). It can be linearly expanded by thethree projections on three linearly independent directions. In particular, on a cornertetrahedron T we have three outer triangles F j on ∂T . On T ,Π h v = A Π h v · n Π h v · n Π h v · n = A p p p , where p , p and p are scalar P k polynomials, and A is a 3 × p vanishes on F and p = λ q k − on T , where λ is a barycentric coordinate of T (which is a linear function assuming 0 on F ), and q k − is a P k − ( T ) polynomial. Let p k ∈ P k ( T ) be an anti- x -derivative of( n ) q k − , i.e., ( ∇ p k ) = ( n ) q k − . Note that ( ∇ p k ) and ( ∇ p k ) can be anything(of y and z functions) which result in zero integrals below. By (2.10) and the thirdand the fourth equations of (2.7), since ∇ · Π h v = 0, we get (cid:90) T λ q k − d x = (cid:90) T Π h v · ( n q k − ) d x = 0 . Since λ > T , we conclude with q k − = 0 and p = 0. Repeating the analysiswe get p = p = 0 and Π h v = 0 on T .Adding the equations (2.9) and (2.10) to (2.7), T would be a new corner tetrahe-dron with three no-flux boundary triangles. Repeating the estimates on T , it wouldlead Π h v = 0 on T . Sequentially, we obtain Π h v = 0 on all T i , i.e., on the whole T .For a τ ∈ H (div; Ω) and a v ∈ P k ( T ), we have, by (2.7), (2.9) and (2.10),( ∇ · ( τ − Π h τ ) , v ) T = n (cid:88) i =1 (cid:16) (cid:90) T i ( τ − Π h τ ) · ∇ vd x + (cid:90) ∂T i ( τ − Π h τ ) · n vdS (cid:17) = n (cid:88) i =1 (cid:90) ∂T ( τ − Π h τ ) · n vdS (cid:17) = 0 . That is, (2.4) holds.Since [ P k ( T )] ⊂ Λ k and Π h is uni-solvent, Π h v = v for all v ∈ [ P k ( T )] . Onone size 1 T , by the finite dimensional norm-equivalence and the shape-regularityassumption on sub-triangles, the interpolation is stable in L ( T ), i.e., (cid:107) Π h τ (cid:107) T ≤ C (cid:107) τ (cid:107) T . (2.11)After a scaling, the constant C in (2.11) remains same. It follows that (cid:107) Π h τ − τ (cid:107) ≤ C (cid:88) T ∈T h ( (cid:107) Π h ( τ − p k,T ) (cid:107) T + (cid:107) p k,T − τ (cid:107) T ) ≤ C (cid:88) T ∈T h ( C (cid:107) τ − p k,T (cid:107) T + (cid:107) p k,T − τ (cid:107) T ) ≤ C (cid:88) T ∈T h h k +2 | τ | k +1 ,T = Ch k +2 | τ | k +1 , where p k,T is a k -th Taylor polynomial of τ on T .
3. Mixed Finite Element Method.
In this section, we develop a mixed finiteelement method on polytopal mesh by employing our new mixed elements and obtainoptimal order error estimates for the method. First let V = H (div; Ω) and W = L (Ω). Algorithm 1.
A mixed finite element method for the problem (1.4)-(1.5) seeks ( q h , u h ) ∈ V h × W h satisfying ( a q h , v ) − ( ∇ · v , u h ) = (cid:104) g, v · n (cid:105) ∂ Ω ∀ v ∈ V h , (3.1) ( ∇ · q h , w ) = ( f, w ) ∀ w ∈ W h . (3.2)We introduce a norm (cid:107) v (cid:107) V for any v ∈ V as follows: (cid:107) v (cid:107) V = (cid:107) v (cid:107) + (cid:107)∇ · v (cid:107) . (3.3) Lemma 3.1.
There exists a positive constant β independent of h such that for all ρ ∈ W h , (3.4) sup v ∈ V h ( ∇ · v , ρ ) (cid:107) v (cid:107) V ≥ β (cid:107) ρ (cid:107) . Proof . For any given ρ ∈ W h ⊂ L (Ω), it is known [3] that there exists a function˜ v ∈ V such that(3.5) ( ∇ · ˜ v , ρ ) (cid:107) ˜ v (cid:107) V ≥ C (cid:107) ρ (cid:107) , where C > h . By setting v = Π h ˜ v ∈ V h and using(2.5), we have(3.6) (cid:107) v (cid:107) V = (cid:107) Π h ˜ v (cid:107) V ≤ C (cid:107) ˜ v (cid:107) V . Using (2.4), (3.6) and (3.5, we have | ( ∇ · v , ρ ) |(cid:107) v (cid:107) V = | ( ∇ · Π h ˜ v , ρ ) |(cid:107) v (cid:107) V ≥ | ( ∇ · ˜ v , ρ ) | C (cid:107) ˜ v (cid:107) V ≥ β (cid:107) ρ (cid:107) , for a positive constant β . This completes the proof of the lemma. Theorem 3.2.
Let ( q h , u h ) ∈ V h × W h be the mixed finite element solution of(3.1)-(3.2). Then, there exists a constant C such that (3.7) (cid:107) q − q h (cid:107) V + (cid:107) u − u h (cid:107) ≤ Ch k +1 ( | q | k +1 + | u | k +1 ) . Proof . Let e h = Π h q − q h and (cid:15) h = Q h u − u h , where Q h is the element-wisedefined L projection onto P k ( T ) on each element T . The differences of (1.4)-(1.5)and (3.1)-(3.2) imply( a ( q − q h ) , v ) − ( ∇ · v , u − u h ) = 0 ∀ v ∈ V h , (3.8) ( ∇ · ( q − q h ) , w ) = 0 ∀ w ∈ W h . (3.9)By adding ( a Π h q h , v ) to the both sides of (3.8) and using the definition of Q h , (3.8)becomes ( a e h , v ) − ( ∇ · v , (cid:15) h ) = ( a (Π h q − q ) , v ) . (3.10)It follows from (2.4) and (3.9) that for w ∈ W h ( ∇ · e h , w ) = ( ∇ · (Π h q − q h ) , w ) = ( ∇ · ( q − q h ) , w ) = 0 . (3.11)Combining (3.10)-(3.11), we have for all ( v , w ) ∈ V h × W h ( a e h , v ) − ( ∇ · v , (cid:15) h ) = ( a (Π h q − q ) , v ) , (3.12) ( ∇ · e h , w ) = 0 . (3.13)Letting v = e h in (3.12) and using (3.11), we have( a e h , e h ) = ( a (Π h q − q ) , e h ) , which gives (cid:107) Π h q − q h (cid:107) V ≤ Ch k +1 | q | k +1 . (3.14)It follows from (3.12) and (3.14) that for all v ∈ V h ( ∇ · v , (cid:15) h ) ≤ | ( a e h , v ) | + | ( a (Π h q − q ) , v ) | ≤ Ch k +1 (cid:107) q (cid:107) k +1 (cid:107) v (cid:107) V . (3.15)The inf-sup condition (3.4) and the estimate (3.15) yield (cid:107) Q h u − u h (cid:107) ≤ Ch k +1 (cid:107) q (cid:107) k +1 . (3.16)It follows from (3.14) and (3.16)(3.17) (cid:107) Π h q − q h (cid:107) V + (cid:107) Q h u − u h (cid:107) ≤ Ch k +1 | q | k +1 . The error bound (3.7) follows from the triangle inequality and (3.17) and we haveproved the theorem.To obtain superconvergence for u h , we consider the dual system: seek ( ψ , θ ) ∈ H (div; Ω) × L (Ω) such that( a ψ , v ) − ( ∇ · v , θ ) = 0 ∀ v ∈ H (div; Ω) , (3.18) ( ∇ · ψ , w ) = ( Q h u − u h , w ) ∀ w ∈ L (Ω) . (3.19)Assume that the following regularity holds(3.20) (cid:107) ψ (cid:107) + (cid:107) θ (cid:107) ≤ C (cid:107) Q h u − u h (cid:107) . Theorem 3.3.
Let ( q h , u h ) ∈ V h × W h be the mixed finite element solution of(3.1)-(3.2). Assume that (3.20) holds true. Then, there exists a constant C such that (3.21) (cid:107) Q h u − u h (cid:107) ≤ Ch k +2 ( | q | k +1 + | u | k +1 ) . Proof . Letting w = Q h u − u h in (3.19) and using (2.4), (3.8), (3.18), (3.9), (3.7)and (3.20), we have (cid:107) Q h u − u h (cid:107) = ( ∇ · ψ , Q h u − u h )= ( ∇ · Π h ψ , Q h u − u h )= (Π h ψ , a ( q − q h ))= (Π h ψ − ψ , a ( q − q h )) + ( ψ , a ( q − q h ))= (Π h ψ − ψ , a ( q − q h )) + ( ∇ · ( q − q h ) , θ )= (Π h ψ − ψ , a ( q − q h )) + ( ∇ · ( q − q h ) , θ − Q h θ ) ≤ Ch k +2 (cid:107) q (cid:107) k +1 (cid:107) Q h u − u h (cid:107) , which implies (3.21) and we have proved the theorem.
4. Numerical Example.
We solve problem (1.1)–(1.3) on the unit square do-main with the exact solution q = (cid:18) π sin( πy ) cos( πx ) π sin( πx ) cos( πy ) (cid:19) , u = sin( πx ) sin( πy ) . (4.1)We first use quadrilateral grids. To avoid asymptotic parallelograms under nestedrefinements, we use fixed types of quadrilaterals in our multi-level grids, shown inFigure 4.1. We list the computational results in Table 4.1. As proved, we have oneorder of super-convergence for both u h and q h .0 Fig. 4.1 . The first three levels of grids, for Table 4.1.
Table 4.1
Error profiles and convergence rates on grids shown in Figure 4.1 for (4.1) . level (cid:107) Q h u − u h (cid:107) rate (cid:107) Π h q − q h (cid:107) V rateby the Λ - P mixed element6 0.1464E-03 2.00 0.5185E-01 1.007 0.3660E-04 2.00 0.2593E-01 1.008 0.9151E-05 2.00 0.1296E-01 1.00by the Λ - P mixed element6 0.1072E-05 3.00 0.4103E-03 2.007 0.1340E-06 3.00 0.1025E-03 2.008 0.1674E-07 3.00 0.2563E-04 2.00by the Λ - P mixed element5 0.5704E-06 4.00 0.1878E-03 3.006 0.3567E-07 4.00 0.2349E-04 3.007 0.2231E-08 4.00 0.2937E-05 3.00by the Λ - P mixed element3 0.1403E-05 5.84 0.1559E-03 4.884 0.2765E-07 5.66 0.5969E-05 4.715 0.6808E-09 5.34 0.2837E-06 4.39Next we solve the same problem (4.1) on a type of grids with quadrilaterals andhexagons, shown in Figure 4.2. We list the result of computation in Table 4.2 wherewe obtain one order of superconvergence in all cases.We solve 3D problem (1.1)–(1.3) on the unit cube domain Ω = (0 , with theexact solution q = (1 − x ) ( y − y )( z − z )2 ( x − x ) (1 − y )( z − z )2 ( x − x ) ( y − y )(1 − z ) ,u = 2 ( x − x ) ( y − y )( z − z ) . (4.2)Here we use a uniform wedge-type (polyhedron with 2 triangle faces and 3 rectanglefaces) grids, shown in Figure 4.3. Here each wedge is subdivided in to three tetra-hedrons with three rectangular faces being cut, when defining piecewise RT k element1 Fig. 4.2 . The first three levels of quadrilateral-hexagon grids, for Table 4.2.
Table 4.2
Error profiles and convergence rates on grids shown in Figure 4.2 for (4.1) . level (cid:107) Q h u − u h (cid:107) rate (cid:107) Π h q − q h (cid:107) V rateby the Λ - P mixed element6 0.1523E-03 2.00 0.5282E-01 1.007 0.3808E-04 2.00 0.2641E-01 1.008 0.9520E-05 2.00 0.1321E-01 1.00by the Λ - P mixed element6 0.1015E-05 3.00 0.3958E-03 2.007 0.1269E-06 3.00 0.9893E-04 2.008 0.1586E-07 3.00 0.2473E-04 2.00by the Λ - P mixed element5 0.8069E-07 4.01 0.2283E-04 3.056 0.5038E-08 4.00 0.2830E-05 3.017 0.3149E-09 4.00 0.3530E-06 3.00by the Λ - P mixed element3 0.9106E-06 5.72 0.1176E-03 4.834 0.2080E-07 5.45 0.4844E-05 4.605 0.5735E-09 5.18 0.2481E-06 4.29 (cid:0)(cid:0) (cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64) (cid:0)(cid:0) (cid:0)(cid:0) (cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64) (cid:64)(cid:64)(cid:64)(cid:64) (cid:64)(cid:64)(cid:64)(cid:64) (cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:0)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64) (cid:64)(cid:64)(cid:64)(cid:64) (cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64) (cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64) (cid:64)(cid:64)(cid:64)(cid:64)(cid:64)(cid:64) Fig. 4.3 . The first three levels of wedge grids used in Table 4.3. Λ k . The results are listed in Table 4.3, confirming the theory. REFERENCES Table 4.3
Error profiles and convergence rates on grids shown in Figure 4.3 for (4.2) . level (cid:107) Q h u − u h (cid:107) rate (cid:107) Π h q − q h (cid:107) V rateby the 3D Λ - P mixed element5 0.0044197 2.0 0.5802877 1.06 0.0011145 2.0 0.2909994 1.07 0.0002793 2.0 0.1456072 1.0by the 3D Λ - P mixed element4 0.0049106 2.9 0.5234688 2.05 0.0006228 3.0 0.1317047 2.06 0.0000782 3.0 0.0329830 2.0by the 3D Λ - P mixed element4 0.0004943 4.0 0.1005523 3.05 0.0000310 4.0 0.0126031 3.06 0.0000019 4.0 0.0015765 3.0by the 3D Λ - P mixed element3 0.0006668 5.0 0.1986805 3.94 0.0000207 5.0 0.0125641 4.05 0.0000006 5.0 0.0007875 4.0 [1] D. Arnold and F. Brezzi, Mixed and nonconforming finite element methods: implementation,postprocessing and error estimates, RAIRO Model. Math. Anal. Numer., 19 (1985), 7-32.[2] I. Babuska, The finite element method with Lagrange multipliers, Numer. Math., 20 (1973),179-192.[3] F. Brezzi and M. Fortin, Mixed and Hybrid Finite Elements, Springer-Verlag, New York, 1991.[4] F. Brezzi, J. Douglas, R. Duran and M. Fortin, Mixed finite elements for second order ellipticproblems in three variables, Numer. Math., 51 (1987), 237-250.[5] F. Brezzi, J. Douglas, and L.D. Marini, Two families of mixed finite elements for second orderelliptic problems, Numer. Math., 47 (1985), 217-235.[6] W. Chen and Y. Wang, Minimal Degree H ( curl ) and H ( divdiv