A new and simple implementation of the element-local L²-projected continuous finite element method

Abstract

In this paper we present a new and simple method, modified from the original element-local L²-projected continuous finite element method, to resolve some static electromagnetic problems in two dimension. This numerical method is established on the least-square process to minimize the total error’s energy, where the residual due to an additional element-local projection is considered. Many numerical experiments are presented for the electrostatic problem, the Maxwell’s equations and the Maxwell’s eigenvalue problem, and the numerical performances are excellent with the optimal convergence orders even for the non-H¹ solution.

Introduction

In this paper we would like to continue to carry out some deep studies on the continuous finite element (C⁰-FE) method for the problems of static electromagnetic fields, including the electrostatic problem, Maxwell’s equations, and the Maxwell’s eigenvalue problems. These models are famous and important, in both theory and numerics.

As a typical problem, let us begin this topic with the electrostatic problem; the others are similar. The electrostatic problem in an isotropic homogeneous medium is described as the div-curl system

$$\mathrm{curl}\mathit{u}=f\text{,}\mathrm{div}\mathit{u}=g\text{,}\mathrm{in}\mathrm{\Omega }\text{;}$$

with a suitable boundary condition, for example, $\mathit{u}\times \mathit{n}=0$ on $\mathrm{\Gamma }=\partial \mathrm{\Omega }$ will be considered in this paper. Here $\mathrm{\Omega }$ is a two-dimensional bounded connected domain with Lipschitz continuous boundary, and f and g are given data in $L^2(\mathrm{\Omega })$.

Till now, many numerical methods have been given and studied for (1) in lots of literatures. If the exact solution of (1) belongs to ${(H^1(\mathrm{\Omega }))}^2$, it is natural to adopt the C⁰-FE method. The existing works [1], [2], [3] have pointed out that this method performs very well for the smooth H¹-solution which means $\mathit{u}\in {(H^r(\mathrm{\Omega }))}^2$ with $r⩾1$ in this paper. However, the exact solution of (1) may have lower regularity [4], i.e., it merely belongs to ${(H^r(\mathrm{\Omega }))}^2$ with a certain constant $r\in (0\text{,}1)$, and named as the non-H¹ solution in this paper. In this situation, the C⁰-FE approximation may give a wrong or very bad convergence; please refer to the papers [5], [6], [7], [8]. In order to overcome this obstacle and obtain the satisfying approximation solution to the non-H¹ solution, some nonstandard finite element approximation techniques using C⁰-FE space are presented in last decade, for example, the singular field method [9], the weighted formulation [10], the negative method [8], [11], and the first order system LL^∗ method [7].

Recently, the element-local ${\mathrm{L}}^2$-projected stabilized C⁰-FE method are presented in [6], [12], which is referred to as the original EL²P method in this paper. Theory analysis and numerical experiments have shown this method works very well to solve the non-H¹ solution for the problems of static electromagnetic fields. In this paper we would like to focus our eyesight on this original EL²P method, and make a deep study on the implementation of this method. We want to establish a new and simple method based on the original EL²P method, which is easy to understand and implement.

To state our idea clearly, let us review here the implementation of the original EL²P method in [6] to resolve problem (1). Let $U_h$ be the C⁰-FE space satisfying the boundary condition. The approximation solution ${\mathit{u}}_h\in U_h$ are sought by the variational form for all ${\mathit{v}}_h\in U_h$,

$${\mathcal{M}}_h({\mathit{u}}_h\text{,}{\mathit{v}}_h)+{\mathcal{S}}_h({\mathit{u}}_h\text{,}{\mathit{v}}_h)=Q_h({\mathit{v}}_h)+Z_h(f\text{,}g\text{;}{\mathit{v}}_h)\text{.}$$

Here, two element-local L²-projections, denoted by $R_h^{\mathrm{c}}$ and $R_h^{\mathrm{d}}$ respectively, are embedded in the standard variational form of this problem, in order to overcome the difficulty due to the large kernel space about the div and curl operators. This treatment leads to

$${\mathcal{M}}_h({\mathit{u}}_h\text{,}{\mathit{v}}_h)={(R_h^{\mathrm{c}}\mathrm{curl}{\mathit{u}}_{\mathit{h}}\text{,}{R}_h^{\mathrm{c}}\mathrm{curl}{\mathit{v}}_{\mathit{h}})}_0+{(R_h^{\mathrm{d}}\mathrm{div}{\mathit{u}}_h\text{,}{R}_h^{\mathrm{d}}\mathrm{div}{\mathit{v}}_h)}_0\text{,}$$

$$Q_h({\mathit{v}}_h)={(f\text{,}{R}_h^{\mathrm{c}}\mathrm{curl}{\mathit{v}}_h)}_0+{(g\text{,}{R}_h^{\mathrm{d}}\mathrm{div}{\mathit{v}}_h)}_0\text{,}$$

where ${(·\text{,}·)}_0$ is the usual L²-inner product in $\mathrm{\Omega }$. Furthermore, in order to remedy the instability resulted from the above two element-local L²-projections, the authors introduced two mesh-dependent stabilization terms ${\mathcal{S}}_h({\mathit{u}}_h\text{,}{\mathit{v}}_h)$ and $Z_h(f\text{,}g\text{;}{\mathit{v}}_h)$. For more information, please see [6].

The success of the original EL²P method comes from the element-local L²- projections, as well as the stabilization mechanism. However, we think that the latter is the Achilles’ heel of this method, since the numerical treatments involved in ${\mathcal{S}}_h({\mathit{u}}_h\text{,}{\mathit{v}}_h)$ and $Z_h(f\text{,}g\text{;}{\mathit{v}}_h)$ are a little complex. In the original EL²P method, the authors introduced additional auxiliary shape functions, precisely the bubble functions associated with each face and/or each element to the finite element subspace $U_h$. Thus we have to take more efforts to define those complex bubbles and take more time cost to derive the element stiff matrices in the finite element subroutines.

In this paper, we would like to present a new and simple treatment to achieve enough stability. This new development comes from the techniques on the element-local L²-projections, and the stabilization mechanism is established completely in the simple framework of least-squares method. Based on these ideas, we would like to seek the approximation solution ${\mathit{u}}_h\in U_h$, such that for any test function ${\mathit{v}}_h\in U_h$, there holds the variational form

$${\mathcal{M}}_h({\mathit{u}}_h\text{,}{\mathit{v}}_h)+{({\mathit{u}}_h-R_h^{\mathrm{s}}{\mathit{u}}_h\text{,}{\mathit{v}}_h-R_h^{\mathrm{s}}{\mathit{v}}_h)}_0=Q_h({\mathit{v}}_h)\text{,}$$

where a new element-local L²-projection $R_h^{\mathrm{s}}$ is introduced here, in additional to two element-local L²-projections, $R_h^{\mathrm{c}}$ and $R_h^{\mathrm{d}}$, same as those in the original EL²P method. The detailed implementation and discussion will be given in Section 2. For convenience, in this paper we would like to refer to the above method as the new EL²P method.

The stability mechanism based on the least-square method is the highlights of this new EL²P method. On comparison with the original method, the stability mechanism now is reduced to a compact form, which depends only on the element-local L²-projection $R_h^{\mathrm{s}}$ and is independent of the given data f and g. As a result, the new method is simpler and easier to be implemented, and much time cost can be saved. Also, this simple implementation in the new EL²P method can be applied directly and easily to the numerical simulation on the Maxwell’s equations and the Maxwell’s eigenvalue problems. Further, this simple treatment in the new EL²P method still provides excellent performances for both the H¹-solution and the non-H¹-solution. In many numerical tests, the numerical result is almost same as or better than that given by the original EL²P method. For example, when this new method is applied to some benchmarks for the Maxwell’s equations and the Maxwell’s eigenvalues problem, the relative errors seem more superior to that given in [12]. Finally, it is worthy to mention that this new method also has the coercivity property and the optimal error estimate as the original EL²P method, in numerical experiments, although in this paper we only prove this issue for a special case.

The remaining part of this paper is organized as follows. Section 2 is the main body of this paper, where the new EL²P method is established for the electrostatics problem. Two theoretical conclusions about the coercivity property and the optimal error estimate for the non-H¹ solution are then presented in a special case. The detailed proofs are derived in the appendix. Finally we present some numerical tests to verify the optimal convergence order and show the nice saving in time cost. As an extension and application of this new method, we would like in Section 3 to apply this new and simple implementation to resolve the Maxwell’s equations and its eigenvalue problems. The remark conclusion and future work are given in Section 4.