• Views 1,384
• Citations 0
• ePub 56
• PDF 687
`Advances in Numerical AnalysisVolume 2013 (2013), Article ID 189045, 9 pageshttp://dx.doi.org/10.1155/2013/189045`
Research Article

## Mixed Finite Element Methods for the Poisson Equation Using Biorthogonal and Quasi-Biorthogonal Systems

School of Mathematical and Physical Sciences, University of Newcastle, Callaghan, NSW 2308, Australia

Received 10 October 2012; Revised 20 February 2013; Accepted 25 February 2013

Academic Editor: Norbert Heuer

Copyright © 2013 Bishnu P. Lamichhane. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

We introduce two three-field mixed formulations for the Poisson equation and propose finite element methods for their approximation. Both mixed formulations are obtained by introducing a weak equation for the gradient of the solution by means of a Lagrange multiplier space. Two efficient numerical schemes are proposed based on using a pair of bases for the gradient of the solution and the Lagrange multiplier space forming biorthogonal and quasi-biorthogonal systems, respectively. We also establish an optimal a priori error estimate for both finite element approximations.

#### 1. Introduction

In many practical situations, it is important to compute dual variables of partial differential equations more accurately. For example, the gradient of the solution is the dual variable in case of the Poisson equation, whereas the stress or pressure variable is the dual variable in case of elasticity equation. Working with the standard finite element approach these variables should be obtained a posteriori by differentiation, which will result in a loss of accuracy. In these situations, a mixed method is often preferred as these variables can be directly computed using a mixed method.

In this paper, we introduce two mixed finite element methods for the Poisson equation using biorthogonal or quasi-biorthogonal systems. Both formulations are obtained by introducing the gradient of the solution of Poisson equation as a new unknown and writing an additional variational equation in terms of a Lagrange multiplier. This gives rise to two additional vector unknowns: the gradient of the solution and the Lagrange multiplier. In order to obtain an efficient numerical scheme, we carefully choose a pair of bases for the space of the gradient of the solution and the Lagrange multiplier space in the discrete setting. Choosing the pair of bases forming a biorthogonal or quasi-biorthogonal system for these two spaces, we can eliminate the degrees of freedom associated with the gradient of the solution and the Lagrange multiplier and arrive at a positive definite formulation. The positive definite formulation involves only the degrees of freedom associated with the solution of the Poisson equation. Hence a reduced system is obtained, which is easy to solve. The first formulation is discretized by using a quasi-biorthogonal system, whereas the second one, which is a stabilized version of the first one, is discretized using a biorthogonal system.

There are many mixed finite element methods for the Poisson equation [18]. However, all of them are based on the two-field formulation of the Poisson equation and hence are not amenable to the application of the biorthogonal and quasi-biorthogonal systems. We need a three-field formulation to apply the biorthogonal and quasi-biorthogonal systems which leads to a symmetric formulation (see [9] for a three-field formulation in linear elasticity). The use of biorthogonal and quasi-biorthogonal systems allows an easy static condensation of the auxiliary variables (gradient of the solution and Lagrange multiplier) leading to a reduced linear system. These variables can be recovered just by inverting a diagonal matrix. Therefore, in this paper we present three-field formulations of the Poisson equation to apply the biorthogonal and quasi-biorthogonal systems.

The structure of the rest of the paper is as follows. In Section 2, we introduce our three-field formulations and show their well-posedness. We propose finite element methods for both formulations and prove a priori error estimates in Section 3. Finally, a short conclusion is drawn in Section 4.

#### 2. Two Three-Field Formulations of the Poisson Equation

In this section we introduce two three-field mixed formulations of the Poisson problem. Let , , be a bounded convex domain with polygonal or polyhedral boundary with the outward pointing normal on . There are many mixed formulations of the Poisson equation with Dirichlet boundary condition Neumann boundary condition or a mixture of these two (see [26, 8, 10]).

We start with the following minimization problem for the Poisson problem: with where, for simplicity, we have assumed the Dirichlet boundary condition, , and Let and , and for two vector-valued functions   and  , the Sobolev inner product on the Sobolev space is defined as where with , for , and the norm is induced from this inner product, where . Note that is the standard inner product on the Sobolev space .

Our new formulation is obtained by introducing an auxiliary variable   such that the minimization problem (4) is rewritten as the following constrained minimization problem: where the natural norm for the product space is for . We write a weak variational equation for equation   in terms of the Lagrange multiplier space to obtain the saddle-point problem of the minimization problem (8). The saddle-point formulation is to find such that where A similar three-field formulation—called the Hu-Washizu formulation—is very popular in designing finite element methods to alleviate locking effect in linear elasticity [1113].

In order to show that the saddle-point problem (9) has a unique solution, we want to apply a standard saddle-point theory [3, 4, 10]. To this end, we need to show the following  three conditions of well-posedness. The linear form , the bilinear forms and are continuous on the spaces on which they are defined. The bilinear form is coercive on the kernel space defined as The bilinear form satisfies the inf-sup condition: for a constant .

The linear form   and the two bilinear forms and are continuous by the Cauchy-Schwarz inequality. The coercivity of the bilinear form on the kernel space follows from Poincaré inequality: It remains to show that the bilinear form satisfies the inf-sup condition.

Lemma 1. The bilinear form satisfies the inf-sup condition; that is,

Proof. Using in the expression on the left above, we have

To summarize we have proved the following theorem.

Theorem 2. The saddle-point problem (9) has a unique solution and

The main difficulty in this formulation is that the bilinear form is not coercive on the whole space . It is only coercive on the subspace of . One has to consider this difficulty in developing a finite element method for this problem. One approach to make the bilinear form coercive on the whole space is to stabilize it. There are many approaches to stabilize this formulation (see [2, 5, 6]). Here we modify the bilinear form to get a consistent stabilization as proposed in [14] for the Mindlin-Reissner plate so that it is now coercive on the whole space . This is obtained by adding the term to the bilinear form . Thus our problem is to find such that where

Theorem 3. The saddle-point problems (9) and (18) yield the same solution.

Proof. By construction the solution to (18) is also the solution to (9). Moreover, as the second equation of (18) yields , we can substitute this into the first equation and set to get . This proves that the solution to (18) is also the solution to (9).

Here the bilinear form is coercive on the whole space from the triangle and Poincaré inequalities.

Lemma 4. The bilinear form is coercive on . That is,

Proof. We start with the triangle inequality From Poincaré inequality, there exists a constant such that and hence we have

Moreover, the other conditions of well-posedness for this new saddle-point formulation (18) can also be shown exactly as for (9).

In the following, we present finite element approximations for both proposed three-field formulations. In the first step, we propose a finite element method for the nonstabilized formulation (9), which is based on quasi-biorthogonal systems. In the second step, we present a finite element method for the stabilized formulation (18), which is based on biorthogonal systems. The first method works only for simplicial elements, whereas the second method works for meshes of parallelotopes and simplices.

#### 3. Finite Element Approximation and a Priori Error Estimate

Let be a quasi-uniform partition of the domain in simplices or parallelotopes having the mesh size . Let be a reference simplex or -cube, where the reference simplex is defined as and a -cube . Note that for a simplex is a triangle, and for it is a tetrahedron. Similarly, a -cube is a square, whereas, -cube is a cube.

The finite element space is defined by affine maps from a reference element to a physical element . Let be the space of bilinear/trilinear polynomials in if is a reference square/cube or the space of linear polynomials in if is a reference simplex. Then the finite element space based on the mesh is defined as the space of continuous functions whose restrictions to an element are obtained by maps of bilinear, trilinear, or linear functions from the reference element; that is, (See [3, 10, 15]). We start with some abstract assumptions for two discrete spaces and to discretize the gradient of the solution and the Lagrange multiplier spaces, where and both are piecewise polynomial spaces with respect to the mesh and are both subsets of . Explicit construction of local basis functions for these spaces for the mixed finite element method (9) is given in Section 3.1 and for the mixed finite element method (18) is given in Section 3.2.

We assume that and satisfy the following assumptions (see also [16, 17]).

Assumption 1. (i)  .
(ii) There is a constant independent of the mesh such that
(iii) The space has the approximation property:
(iv) The space has the approximation property:
Now we define biorthogonality and quasi-biorthogonality (see [16, 18, 19]).

Definition 5. The pair of bases of and of is called  biorthogonal if the resulting Gram matrix is diagonal. The pair of bases of and of is called  quasi-biorthogonal and the resulting Gram matrix    is called  quasidiagonal if    is of the form where and are diagonal matrices and is a sparse rectangular matrix.

We recall that a Gram matrix of two sets of basis functions and is the matrix whose th entry is It is worth noting that the quasi-diagonal matrix is inverted by inverting two diagonal matrices and scaling the rectangular matrix.

##### 3.1. Finite Element Approximation for the Nonstabilized Formulation (9)

We now turn our attention to a finite element approximation for the nonstabilized formulation (9). This approximation works only for simplicial meshes and is based on a quasi-biorthogonal system as in [19]. Let be the space of bubble functions, where is a constant and is the set of barycentric coordinates on . Let . Explicit construction of basis functions of and on a reference element is provided below.

The finite element problem is to find such that We note that the standard finite element space is enriched with elementwise defined bubble functions to obtain the space , which is done to obtain the coercivity of the bilinear form on the kernel space defined as (see Lemma 7) To simplify the numerical analysis of the finite element problem, we introduce a quasiprojection operator: , which is defined as   This type of operator was introduced in [9] to obtain the finite element interpolation of nonsmooth functions satisfying boundary conditions and is used in [20] in the context of mortar finite elements. The definition of allows us to write the weak gradient as where operator is applied to the vector component-wise. We see that is well defined due to Assumption 1(i) and (ii). Furthermore, the restriction of to is the identity. Hence is a projection onto the space . We note that is not the orthogonal projection onto but an oblique projection onto see [21, 22]. The stability and approximation properties of in and -norm can be shown as in [16, 23]. In the following, we will use a generic constant , which will take different values at different places but will be always independent of the mesh size .

Lemma 6. Under Assumption 1(i)–(iii) and for and

We need to construct another finite element space that satisfies Assumption 1(i)–(iii) and also leads to an efficient numerical scheme. The main goal is to choose the basis functions for and so that the matrix associated with the bilinear form is easy to invert. To this end, we consider the algebraic form of the weak equation for the gradient of the solution, which is given by where and  are the solution vectors for and   and is the Gram matrix between the bases of and . Hence if the matrix is diagonal or quasi-diagonal, we can easily eliminate the degrees of freedom corresponding to   and  . This then leads to a formulation involving only one unknown . Statically condensing out variables and we arrive at the variational formulation of the reduced system given by (68) (see Section 3.3). Note that the matrix is the Gram matrix between the bases of and .

We cannot use a biorthogonal system in this case as the stability will be lost (see the Appendix). This motivates our interest in a quasi-biorthogonal system. We now show the construction of the local basis functions for and so that they form a quasi-biorthogonal system. The construction is shown on the reference triangle in the two-dimensional case and on the reference tetrahedron in the three-dimensional case.

We now start with the reference triangle in the two-dimensional case. Let defined as be the local basis functions of on the reference triangle associated with three vertices , and , and one barycenter . Note that is the standard bubble function used, for example, in the minielement for the Stokes equations [24].

Then if we define the local basis functions of as the two sets of basis functions and form a quasi-biorthogonal system on the reference triangle. These four basis functions of are associated with three vertices , and , and the barycenter of the reference triangle.

Similarly, on the reference tetrahedron, the local basis functions of defined as and the local basis functions of defined as are quasi-biorthogonal. These five basis functions of or of are associated with four vertices , , , and and the barycenter of the reference tetrahedron. Using the finite element basis functions for , the finite element basis functions for are constructed in such a way that the Gram matrix on the reference element is given by for the two-dimensional case and for the three-dimensional case. After ordering the degrees of freedom in and so that the degrees of freedom associated with the barycenters of elements come after the degrees of freedom associated with the vertices of elements, the global Gram matrix is quasi-diagonal.

Thus in the two-dimensional case, the local basis functions , and are associated with the vertices of the reference triangle, and the function is associated with the barycenter of the triangle and defined elementwise. That means is supported only on the reference element. The situation is similar in the three-dimensional case. Hence using a proper ordering, the support of a global basis function of is the same as that of the global basis function of for . However, the global basis functions of are not continuous on interelement boundaries.

Now we show that satisfies Assumption 1(i)–(iii). As the first assumption is satisfied by construction, we consider the second assumption. Let and set . By using the quasibiorthogonality relation (1) and the quasiuniformity assumption, we get where denotes the meshsize at th vertex. Taking into account the fact that , we find that Assumption 1(ii) is satisfied. Since the sum of the local basis functions of is one, Assumption 1(iii) can be proved as in [23, 25]. Assumption 1(iv) follows by standard arguments.

There are local basis functions of , where of them are associated with vertices and one is associated with the element barycenter. Then the local basis functions of are also divided into two parts: basis functions associated with the vertices and basis function associated with the element. Let be the local basis function associated with the element barycenter of for . Let be the space of piecewise constant functions with respect to , and is defined as Then is a projection operator onto . Moreover, and this implies that This yields the following coercivity result.

Lemma 7. The bilinear form is coercive on the kernel space . In fact,

Proof. Since on the kernel space , we have The proof follows from the result

Remark 8. Using a quadrilateral or hexahedral meshes, we cannot have , and hence the coercivity fails. One needs to modify the finite element method for these meshes.

The inf-sup condition follows exactly as in the continuous setting using Assumption 1(ii).

Lemma 9. There exists a constant such that

Proof. The proof follows exactly as in the continuous setting. We use in the expression to get Assumption 1(ii) then yields the final result.

Hence we have the following stability and approximation result.

Theorem 10. The discrete saddle-point problem (32) has a unique solution , and Moreover, the following estimate holds for the solution of (32)

Owing to Assumption 1(iii)-(iv) the discrete solution converges optimally to the continuous solution.

##### 3.2. Finite Element Approximation for the Stabilized Formulation (18)

After stabilizing formulation (9), the bilinear form (18) is coercive on the whole product space . Therefore, we can use the standard finite element space to discretize the gradient of the solution for the saddle-point problem (18) so that . Thus our discrete saddle-point formulation of (18) is to find so that We can now verify the three conditions of well-posedness for the discrete saddle point problem. The continuity of , , and follows as in the continuous setting, and coercivity is immediate from Lemma 4 as . Explicit construction of basis functions of on a reference element will be given below.

Assumption 1(ii) guarantees that the bilinear form satisfies the inf-sup condition as in the case of the previous finite element method.

We have thus proved the following theorem from the standard theory of saddle-point problem  (see [4]).

Theorem 11. The discrete saddle-point problem (57) has a unique solution and Moreover, the following estimate holds for the solution of (9)

The optimal convergence of the discrete solution is then guaranteed by the standard approximation property of and and Assumption 1(iii).

In the following, we give these basis functions for linear simplicial finite elements in two and three dimensions [18]. The parallelotope case is handled by using the tensor product construction of the one-dimensional basis functions presented in [20, 26]. These basis functions are very useful in the mortar finite element method and extensively studied in one- and two-dimensional cases [16, 23, 26].

Note that we have imposed the condition to get that our Gram matrix is a square matrix. Therefore, the local basis functions for linear/bilinear/trilinear finite spaces are associated with the vertices of triangles, tetrahedra, quadrilaterals, and hexahedra. For example, for the reference triangle , we have where the basis functions , and are associated with three vertices , , and of the reference triangle. For the reference tetrahedron , we have where the basis functions , , , and , associated with four vertices , , , and of the reference tetrahedron (see also [17, 18]). Note that the basis functions are constructed in such a way that the biorthogonality condition is satisfied on the reference element.

The global basis functions for the test space are constructed by gluing the local basis functions together following exactly the same procedure of constructing global finite element basis functions from the local ones. These global basis functions then satisfy the condition of biorthogonality with global finite element basis functions of . These basis functions also satisfy Assumption 1(i)–(iv) (see [17, 18]).

##### 3.3. Static Condensation

Here we want to eliminate   and Lagrange multiplier   from the saddle-point problem (57). We make use of the operator defined previously. Using the biorthogonality relation between the basis functions of and , the action of operator on a function can be written as which tells that the operator is local in the sense to be given below (see also [27]).

Using the property of operator , we can eliminate the degrees of freedom corresponding to   and   so that our problem is to find such that where for the nonstabilized formulation (9) and for the stabilized formulation (18).

Let the bilinear form be defined as for the non-stabilized formulation (9) and for the stabilized formulation (18) with   and  . Due to the biorthogonality or quasi-biorthogonality condition, the action of is computed by inverting a diagonal matrix.

Since the bilinear form is symmetric, the minimization problem (63) is equivalent to the variational problem of finding such that [3, 15]

The coercivity of and the stability result of Lemma 9 immediately yield the coercivity of the bilinear form on the space with the norm for . That is,

#### 4. Conclusion

We have considered two three-field formulations of Poisson problem and shown the well-posedness of them. We have proposed efficient finite element schemes for both of these formulations: one for the stabilized three-field formulation based on biorthogonal systems and the other for the non-stabilized three-field formulation based on quasi-biorthogonal systems. Optimal a priori error estimates are proved for both approaches, and a reduced discrete problem is presented.

#### Appendix

If we want to enrich the finite element space with the bubble function to attain stability, it is not possible to use a biorthogonal system. We consider the local basis functions of for the reference triangle . Let be defined as in (39). Then using the biorthogonal relation for the four local basis functions of with the four basis functions of for the triangle and on , should satisfy since it is associated with the barycenter of the triangle. This leads to instability in the formulation as we suggest that the function should have nonzero integral over (see (46)). This explains why we cannot use a biorthogonal system. We have a similar issue in the three-dimensional case as , resulting in if we impose biorthogonality.

#### Acknowledgments

Support from the new staff Grant of the University of Newcastle is gratefully acknowledged. The author is grateful to the anonymous referees for their valuable suggestions to improve the quality of the earlier version of this work.

#### References

1. Y. Achdou, C. Bernardi, and F. Coquel, “A priori and a posteriori analysis of finite volume discretizations of Darcy's equations,” Numerische Mathematik, vol. 96, no. 1, pp. 17–42, 2003.
2. P. B. Bochev and C. R. Dohrmann, “A computational study of stabilized, low-order ${C}^{0}$ finite element approximations of Darcy equations,” Computational Mechanics, vol. 38, no. 4-5, pp. 310–333, 2006.
3. D. Braess, Finite Elements. Theory, Fast Solver, and Applications in Solid Mechanics, Cambridge University Press, Cambridge, UK, 2nd edition, 2001.
4. F. Brezzi and M. Fortin, Mixed and Hybrid Finite Element Methods, Springer, New York, NY, USA, 1991.
5. F. Brezzi, T. J. R. Hughes, L. D. Marini, and A. Masud, “Mixed discontinuous Galerkin methods for Darcy flow,” Journal of Scientific Computing, vol. 22-23, pp. 119–145, 2005.
6. K. A. Mardal, X.-C. Tai, and R. Winther, “A robust finite element method for Darcy-Stokes flow,” SIAM Journal on Numerical Analysis, vol. 40, no. 5, pp. 1605–1631, 2002.
7. A. Masud and T. J. R. Hughes, “A stabilized mixed finite element method for Darcy flow,” Computer Methods in Applied Mechanics and Engineering, vol. 191, no. 39-40, pp. 4341–4370, 2002.
8. J. M. Urquiza, D. N'Dri, A. Garon, and M. C. Delfour, “A numerical study of primal mixed finite element approximations of Darcy equations,” Communications in Numerical Methods in Engineering with Biomedical Applications, vol. 22, no. 8, pp. 901–915, 2006.
9. B. P. Lamichhane and E. P. Stephan, “A symmetric mixed finite element method for nearly incompressible elasticity based on biorthogonal systems,” Numerical Methods for Partial Differential Equations, vol. 28, no. 4, pp. 1336–1353, 2012.
10. S. C. Brenner and L. R. Scott, The Mathematical Theory of Finite Element Methods, Springer, New York, NY, USA, 1994.
11. B. P. Lamichhane, B. D. Reddy, and B. I. Wohlmuth, “Convergence in the incompressible limit of finite element approximations based on the Hu-Washizu formulation,” Numerische Mathematik, vol. 104, no. 2, pp. 151–175, 2006.
12. J. C. Simo and M. S. Rifai, “A class of mixed assumed strain methods and the method of incompatible modes,” International Journal for Numerical Methods in Engineering, vol. 29, no. 8, pp. 1595–1638, 1990.
13. K. Washizu, Variational Methods in Elasticity and Plasticity, Pergamon Press, 3rd edition, 1982.
14. D. N. Arnold and F. Brezzi, “Some new elements for the Reissner-Mindlin plate model,” in Boundary Value Problems for Partial Differerntial Equations and Applications, pp. 287–292, Masson, Paris, France, 1993.
15. P. G. Ciarlet, The Finite Element Method for Elliptic Problems, North-Holland, Amsterdam, The Netherlands, 1978.
16. B. P. Lamichhane, Higher order mortar finite elements with dual lagrange multiplier spaces and applications [Ph.D. thesis], University of Stuttgart, 2006.
17. B. P. Lamichhane, “A stabilized mixed finite element method for the biharmonic equation based on biorthogonal systems,” Journal of Computational and Applied Mathematics, vol. 235, no. 17, pp. 5188–5197, 2011.
18. B. P. Lamichhane, “A mixed finite element method for the biharmonic problem using biorthogonal or quasi-biorthogonal systems,” Journal of Scientific Computing, vol. 46, no. 3, pp. 379–396, 2011.
19. B. P. Lamichhane, S. Roberts, and L. Stals, “A mixed finite element discretisation of thin-plate splines,” ANZIAM Journal, vol. 52, pp. C518–C534, 2011, Proceedings of the 15th Biennial Computational Techniques and Applications Conference.
20. B. P. Lamichhane and B. I. Wohlmuth, “Biorthogonal bases with local support and approximation properties,” Mathematics of Computation, vol. 76, no. 257, pp. 233–249, 2007.
21. A. Galántai, Projectors and Projection Methods, vol. 6, Kluwer Academic Publishers, Dordrecht, The Netherlands, 2004.
22. D. B. Szyld, “The many proofs of an identity on the norm of oblique projections,” Numerical Algorithms, vol. 42, no. 3-4, pp. 309–323, 2006.
23. C. Kim, R. D. Lazarov, J. E. Pasciak, and P. S. Vassilevski, “Multiplier spaces for the mortar finite element method in three dimensions,” SIAM Journal on Numerical Analysis, vol. 39, no. 2, pp. 519–538, 2001.
24. D. N. Arnold, F. Brezzi, and M. Fortin, “A stable finite element for the Stokes equations,” Calcolo, vol. 21, no. 4, pp. 337–344, 1984.
25. B. P. Lamichhane, “A mixed finite element method for non-linear and nearly incompressible elasticity based on biorthogonal systems,” International Journal for Numerical Methods in Engineering, vol. 79, no. 7, pp. 870–886, 2009.
26. B. I. Wohlmuth, Discretization Methods and Iterative Solvers Based on Domain Decomposition, vol. 17 of Lecture Notes in Computational Science and Engineering, Springer, Heidelberg, Germany, 2001.
27. M. Ainsworth and J. T. Oden, A Posteriori Error Estimation in Finite Element Analysis, Wiley-Interscience, New York, NY, USA, 2000.