A note on the optimal degree of the weak gradient of the stabilizer free weak Galerkin finite element method

TL;DR

This paper investigates the optimal polynomial degree for the stabilizer free weak Galerkin method to balance convergence rate and computational efficiency, providing rigorous proof for the optimal choice.

Contribution

It determines the optimal polynomial degree j0 for the stabilizer free weak Galerkin method with mathematical validation, enhancing its practical implementation.

Findings

Optimal degree j0 ensures convergence without numerical locking

Mathematical proof of the optimal degree j0

Improved efficiency of the weak Galerkin method

Abstract

Recently, a new stabilizer free weak Galerkin method (SFWG) is proposed, which is easier to implement. The idea is to raise the degree of polynomials j for computing weak gradient. It is shown that if j>=j0 for some j0, then SFWG achieves the optimal rate of convergence. However, large j will cause some numerical difficulties. To improve the efficiency of SFWG and avoid numerical locking, in this note, we provide the optimal j0 with rigorous mathematical proof.

Tables

Table 1. Table 1 : Error profiles and convergence rates for P k subscript 𝑃 𝑘 P_{k} elements with P k + 1 2 superscript subscript 𝑃 𝑘 1 2 P_{k+1}^{2} weak gradient, ( k = 1 , 2 ) 𝑘 1 2 (k=1,2) .

$k$	$h$	$\Vert \left|u_h-Q_{h}u\right|\Vert$	Rate	$\Vert u_0-Q_{0}u\Vert$	Rate
	1/2	8.8205E-01	-	8.0081E-02	-
	1/4	4.3566E-01	1.02	2.5055E-02	1.68
	1/8	2.1565E-01	1.01	6.7091E-03	1.90
1	1/16	1.0747E-01	1.00	1.7092E-03	1.97
	1/32	5.3684E-02	1.00	4.2940E-04	1.99
	1/64	2.6836E-02	1.00	1.0749E-04	2.00
	1/128	1.3417E-02	1.00	2.6871E-05	2.00
	1/2	2.5893E-01	-	1.0511E-02	-
	1/4	6.5338E-02	1.99	1.2889E-03	3.03
	1/8	1.6252E-02	2.00	1.5596E-04	3.05
2	1/16	4.0538E-03	2.09	1.9193E-05	3.02
	1/32	1.0129E-03	2.00	2.3837E-06	3.00
	1/64	2.5321E-04	2.00	2.9713E-07	3.00
	1/128	6.3304E-05	2.00	3.7094E-08	3.00

Table 2. Table 2 : Error profiles and convergence rates for P k subscript 𝑃 𝑘 P_{k} elements with P k + 1 2 superscript subscript 𝑃 𝑘 1 2 P_{k+1}^{2} weak gradient, ( k = 1 , 2 ) 𝑘 1 2 (k=1,2) .

$k$	$N$	$\Vert \left|e_h\right|\Vert$	Rate	$\Vert e_0\Vert$	Rate
	1/2	5.6636E-02	-	5.3064E-03	-
	1/4	2.9670E-02	0.93	1.6876E-03	1.65
	1/8	1.4948E-02	0.99	4.5066E-04	1.90
1	1/16	7.4848E-03	1.00	1.1461E-04	1.98
	1/32	3.7435E-03	1.00	2.8778E-05	1.99
	1/64	1.8719E-03	1.00	7.2024E-06	2.00
	1/128	9.3596E-04	1.00	1.8011E-06	2.00
	1/2	1.3353E-02	-	5.6066E-04	-
	1/4	3.4831E-03	1.94	6.8552E-05	3.03
	1/8	8.7712E-04	1.99	8.2466E-06	3.06
2	1/16	2.1961E-04	2.00	1.0094E-06	3.03
	1/32	5.4967E-05	2.00	1.2482E-07	3.02
	1/64	1.3747E-05	2.00	1.5518E-08	3.00
	1/128	3.4375E-06	2.00	1.9344E-09	3.00

Equations

$$-\Delta u$$

$$u$$

$$(\nabla u,\nabla v)$$

$$V_h$$

$$V_h^0$$

$$({\nabla }_{w}v,\mathbf{q}{)}_T$$

$$(v,w{)}_{{\mathcal{T}}_h}=\sum _{T\in {\mathcal{T}}_h}(v,w{)}_T=\sum _{T\in {\mathcal{T}}_h}{\int }_{T}vwdx,$$

$${⟨v,w⟩}_{\partial {\mathcal{T}}_h}=\sum _{T\in {\mathcal{T}}_h}{⟨v,w⟩}_{\partial T}=\sum _{T\in {\mathcal{T}}_h}{\int }_{\partial T}vwdx.$$

$$A(v,w)$$

$$A(u_h,v)=(f,v_0),\forall v=\{v_0,v_b\}\in V_h^0.$$

$$|||v{|||}^2=\sum _{T\in {\mathcal{T}}_h}({\nabla }_{w}v,{\nabla }_{w}v{)}_T=\sum _{T\in {\mathcal{T}}_h}\Vert {\nabla }_{w}v{\Vert }_T^2.$$

$$\Vert v{\Vert }_{1,h}^2=\sum _{T\in {\mathcal{T}}_h}\left(\Vert \nabla v_0{\Vert }_T^2+h_T^{-1}\Vert v_0-v_b{\Vert }_{\partial T}^2\right).$$

$$\Vert p{\Vert }_{D_1}\le \Vert p{\Vert }_{D_2}\le C{r}^n\Vert p{\Vert }_{D_1},$$

$$\Vert p{\Vert }_{D_1}\le \Vert p{\Vert }_D\le C{r}^n\Vert p{\Vert }_{D_1},\forall p\in P_n(D),$$

$$|e_s|<{\alpha }_0|e_t|;$$

$${\theta }_0<\theta <\pi -{\theta }_0,$$

$$C_1\Vert v{\Vert }_{1,h}\le |||v|||\le C_2\Vert v{\Vert }_{1,h},$$

$$({\nabla }_{w}v,\vec{q}{)}_T=(\nabla v,\vec{q}{)}_T+{⟨v_b-v_0,\vec{q}\cdot \vec{n}⟩}_{\partial T}\forall \vec{q}\in [P_j(T){]}^2.$$

$$(\vec{q},p{)}_T$$

$$\vec{q}\cdot {\vec{n}}_i$$

$${⟨v_b-v_0-\vec{q}\cdot n_0,p⟩}_{e_0}$$

$$\vec{q}={\ell }_2\dots {\ell }_{m-1}{Q}_1\vec{t}=L{Q}_1\vec{t}=Q\vec{t},$$

$${⟨v_b-v_0-{\ell }_2\dots {\ell }_m{Q}_1({\vec{t}}_1\cdot {\vec{n}}_0),p⟩}_{e_0}$$

$$({\ell }_2\dots {\ell }_m{Q}_1,p{)}_T$$

$${⟨{\ell }_2\dots {\ell }_{m-1}{Q}_1,Q_1⟩}_{e_0}=0.$$

$$(y{\ell }_2\dots {\ell }_{m-1}{Q}_2,p{)}_T=0\forall p\in P_{k-1}(T),$$

$$\Vert \vec{q}{\Vert }_T\le {\gamma }_0\Vert \vec{q}\cdot \overrightarrow{n_0}{\Vert }_{e_0},$$

$$(\vec{q}\cdot \overrightarrow{n_0}{)}^2=Q^2{\sin }^2\theta =\vec{q}\cdot \vec{q}{\sin }^2\theta \ge \vec{q}\cdot \vec{q}{\sin }^2{\theta }_0,$$

$$\Vert Q{\Vert }_T\le {\gamma }_0\Vert Q{\Vert }_{e_0},$$

$$e_1=\{(0,y)|0\le y\le 1\},$$

$$d_i$$

$$D_i$$

$$d_i$$

Full text

1 Introduction

In this note, similar to [1], we will consider the following Poisson equation

[TABLE]

as the model problem, where $\Omega$ is a polygonal domain in $\mathbb{R}^{2}$. The variational formulation of the Poisson problem (1)-(2) is to seek $u\in H_{0}^{1}(\Omega)$ such that

[TABLE]

A stabilizer free weak Galerkin finite element method is proposed by Ye and Zhang in [1] as a new method for the solution of Poisson equation on polytopal meshes in 2D or 3D. Let $j$ be the degree of polynomials for the calculation of weak gradient. It has been proved in [1] that there is a $j_{0}>0$ so that the SFWG method converges with optimal order of convergence for any $j\geq j_{0}$. However, when $j$ is too large, the magnitude of the weak gradient can be extremely large, causing numerical instability. In this note, we provide the optimal $j_{0}$ to improve the efficiency and avoid unnecessary numerical difficulties, which has mathematical and practical interests.

2 Weak Galerkin Finite Element Schemes

Let ${\cal T}_{h}$ be a partition of the domain $\Omega$ consisting of polygons in 2D. Suppose that ${\cal T}_{h}$ is shape regular in the sense defined by (2.7)-(2.8) and each $T$ is convex. Let ${\cal E}_{h}$ be the set of all edges in ${\cal T}_{h}$, and let ${\cal E}_{h}^{0}={\cal E}_{h}\setminus\partial\Omega$ be the set of all interior edges. For each element $T\in{\cal T}_{h}$, denote by $h_{T}$ the diameter of $T$, and $h=max_{T\in{\cal T}_{h}}h_{T}$ the mesh size of ${\cal T}_{h}$.

On $T$, let $P_{k}(T)$ be the space of all polynomials with degree no greater than $k$. Let $V_{h}$ be the weak Galerkin finite element space associated with $T\in{\cal T}_{h}$ defined as follows:

[TABLE]

where $k\geq 1$ is a given integer. In this instance, the component $v_{0}$ symbolizes the interior value of $v$, and the component $v_{b}$ symbolizes the edge value of $v$ on each $T$ and $e$, respectively. Let $V^{0}_{h}$ be the subspace of $V_{h}$ defined as:

[TABLE]

For any $v=\{v_{0},v_{b}\}\in V_{h}+H^{1}(\Omega),$ the weak gradient $\nabla_{w}v\in[P_{j}(T)]^{2}$ is defined on $T$ as the unique polynomial satisfying

[TABLE]

where $\mathbf{n}$ is the unit outward normal vector of $\partial T$.

For simplicity, we adopt the following notations,

[TABLE]

[TABLE]

For any $v=\left\{v_{0},v_{b}\right\}$ and $w=\left\{w_{0},w_{b}\right\}$ in $V_{h}$, we define the bilinear forms as follows:

[TABLE]

Stabilizer free Weak Galerkin Algorithm 1**.**

A numerical solution for (1)-(2) can be obtain by finding $u_{h}=\{v_{0},v_{b}\}\in V_{h}^{0}$, such that the following equation holds

[TABLE]

where $A(\cdot,\cdot)$ is defined in (7).

Accordingly, for any $v\in V_{h}^{0}$, we define an energy norm ${|||}\cdot{|||}$ on $V_{h}^{0}$ as:

[TABLE]

An $H^{1}$ norm on $V_{h}^{0}$ is defined as:

[TABLE]

The following well known lemma and corollary will be needed in proving our main result.

Lemma 1**.**

Suppose that $D_{1}\subset D_{2}$ are convex regions in $\mathbb{R}^{d}$ such that $D_{2}$ can be obtain by scaling $D_{1}$ with a factor $r>1$. Then for any $p\in P_{n}(D_{2})$,

[TABLE]

where $C$ depends only on $d\mbox{ and }n$.

Corollary 2.1**.**

Suppose that $D_{1}\mbox{ and }D_{2}$ satisfy the conditions of Lemma 1 and $D_{1}\subseteq D\subseteq D_{2}$. Then

[TABLE]

where $C$ depends only on $d\mbox{ and }n$.

The following Theorem shows that for certain $j_{0}>0$, whenever $j\geq j_{0}$, $\|\cdot\|_{1,h}$ is equivalent to the ${|||}\cdot{|||}$ defined in (9), which is crucial in establishing the feasibility of SFWG. Later on, we will show that $j_{0}$ is optimal.

Theorem 2**.**

(cf.[1], Lemma 3.1) Suppose that $\forall T\in{\cal T}_{h},T$ is convex with at most $m$ edges and satisfies the following regularity conditions: for all edges $e_{t}$ and $e_{s}$ of $T$

[TABLE]

for any two adjacent edges $e_{t}\mbox{ and }e_{s}$ the angle $\theta$ between them satisfies

[TABLE]

where $1\leq\alpha_{0}\mbox{ and }\theta_{0}>0$ are independent of $T$ and $h$. Let $j_{0}=k+m-2\mbox{ or }j_{0}=k+m-3$ when each edge of $T$ is parallel to another edge of $T$. When $deg\nabla_{w}v=j\geq j_{0}$, then there exist two constants $C_{1},C_{2}>0$, such that for each $v=\{v_{0},v_{b}\}\in V_{h}$, the following hold true

[TABLE]

where $C_{1}\mbox{ and }C_{2}$ depend only on $\alpha_{0}\mbox{ and }\theta_{0}$.

Proof.

For simplicity, from now on, all constant are independent of $T$ and $h$, unless otherwise mentioned. They may depend on $k,\alpha_{0},\mbox{ and }\theta_{0}$. Suppose $\nabla_{w}v\in[P_{j}(T)]^{2},v\in P_{k}(T)$. We know that

[TABLE]

Suppose $\partial T=(\cup_{i=0}^{m-1}e_{i})$, and we want to construct a $\vec{q}\in[P_{j_{0}}(T)]^{2}\subseteq[P_{j}(T)]^{2}$, where $j_{0}\leq j$, so that

[TABLE]

where $\vec{n}_{i}$ is the unit outer normal to $e_{i}$. Without loss of generality, we may assume that $e_{0}=\{(x,0)|x\in(0,1)\},\forall(x,y)\in T,y>0$, and one of the vertices of $e_{1}$ is $(0,0)$. Let $\ell_{i}\in P_{1}(T),i=2,\dots,m-1$ be such that $\ell_{i}(e_{i})=0,\ell_{i}(p)\geq 0,\forall p\in T,\mbox{ and }\max_{p\in T}\ell_{i}(p)=1$. Let

[TABLE]

where $\vec{t}$ is a unit tangent vector to $e_{1}$ and $Q_{1}\in P_{j_{0}-m+2}(T)$. It is easy to see that $\vec{q}\cdot\vec{n}_{i}|_{e_{i}}=0,i=1,\dots,m-1$. We want to find $Q_{1}$ so that

[TABLE]

Let $j_{0}-m+2=k$ or $j_{0}=k+m-2$. It is easy to see that if $Q_{1}$ satisfies (17)-(18), then $\vec{q}$ satisfies (13)-(15). In addition, (17)-(18) has $(k+1)(k+2)/2$ equations and the same number of unknowns. To show (17)-(18) has a unique solution $Q_{1}\in P_{k}(T)$ we set $v_{b}-v_{0}=0$. Setting $p=Q_{1}$ in (17) yields

[TABLE]

Since $L=\ell_{2}\dots\ell_{m-1}>0$ inside $e_{0}$, $Q_{1}(x,0)\equiv 0$. Thus $Q_{1}=yQ_{2}$, where $Q_{2}\in P_{k-1}(T)$. Similarly

[TABLE]

implies $Q_{2}=0$ and thus $Q_{1}=0$. Note that for any $j\geq j_{0}$, since $\vec{q}\in[P_{j_{0}}(T)]^{2}\subset[P_{j}(T)]^{2}$, there exists at least one such $\vec{q}$ in $[P_{j}(T)]^{2}$. Note also that when every edge of $\partial T$ is parallel to another one, we can lower the number of $\ell_{i}$ in (16) by one. Thus for such a $T$, we have $j_{0}=m+k-3$. Next, we want to show that the unique solution of (13)-(15) satisfies

[TABLE]

for some $\gamma_{0}>0$. Note that

[TABLE]

where $\theta$ is the angle between $e_{0}$ and $e_{1}$. Thus, (19) is equivalent to

[TABLE]

for some $\gamma_{0}>0$ and we may assume

[TABLE]

because of the assumption (10). Now let $T_{r}$ be the triangle spanned by $e_{0}$ and $e_{1}$ and label the third edge of $T_{r}$ as $e_{t}$. Let $\hat{T}_{r}\subset T_{r}$ be the triangle with vertices $(0,0.25),(0,0.75)$ and $(0.5,0.25)$. Label the edge of $\hat{T}_{r}$ connecting $(0,0.75)$ and $(0.5,0.25)$ as $\hat{e}_{t}$. For $i=2,\dots,m-1$, let

[TABLE]

where $L_{i}$ is the line containing $e_{i}$. Then

[TABLE]

Thus $\forall p\in\hat{T}_{r}$,

[TABLE]

Now we will prove (20) in $3$ steps:

[TABLE]

1. I.

Let $\hat{e}=\{(x,0)|0\leq x\leq 0.75\}$. Then similar to (22), $\forall p\in\hat{e},\ell_{i}(p)\geq\varepsilon,i=2.\dots,m-1$. Thus

[TABLE]

It follows from Corollary 2.1 that

[TABLE] 2. II.

It is easy to see that

[TABLE]

It follows from (22) and Corollary 2.1 that

[TABLE]

Note that

[TABLE]

and

[TABLE]

Since $Ly\geq\dfrac{1}{4}\varepsilon\mbox{ for }(x,y)\in\hat{T}_{r}\mbox{ and }Ly\geq 0,\forall(x,y)\in T$,

[TABLE]

Since $\|\cdot\|_{\hat{T}_{r}}\mbox{ and }\|\cdot\|_{T_{r}}$ are equivalent norms on $P_{k-1}(T_{r})$,

[TABLE]

Let $T^{c}_{e}$ be the circumscribed isosceles-right triangle obtained by scaling $T_{r}$. It is easy to see that

[TABLE]

Since $T_{r}\subseteq T\subseteq T^{c}_{e}$, it follows from (26), (27), and Corollary 2.1 that

[TABLE]

where $\delta_{4}$ depend on (28) and $k$. Letting $p=\bar{Q}$ in (25) yields

[TABLE]

and thus

[TABLE]

Then

[TABLE]

where $N=\max_{(x,y)\in T}y\leq m\alpha_{0}\mbox{ and }M=\max_{(x,y)\in T}x\varepsilon\leq m\alpha_{0}$. It follows from Corollary 2.1 and (26)

[TABLE]

Combining (29) and (24) yields

[TABLE] 3. III.

Applying Corollary 2.1 again yields

[TABLE]

Now we have completed the $3$-step argument.

By a scaling argument, we have

[TABLE]

Plugging $\vec{q}$ into (12) yields

[TABLE]

It follows from (32) that

[TABLE]

Similarly,

[TABLE]

Thus

[TABLE]

By letting $\mathbf{q}=\nabla_{w}v$ in (12), we arrive at

[TABLE]

It follows from the trace inequality and the inverse inequality that

[TABLE]

Thus

[TABLE]

∎

To show that $j_{0}=k+m-2\mbox{ or }j_{0}=k+m-3$ is the lowest possible degree for the weak gradient, we give the following examples.

Example 2.1**.**

Suppose $\Omega=(0,1)\times(0,1)$ and is partitioned as shown in Figure 1.

Let $j=k=1=j_{0}-1$. To see if (8) has a unique solution in $V_{h}^{0}$, we set $f=0$ and $v=u_{h}$. Then $\nabla_{w}u_{h}=0$. By definition,

[TABLE]

Note that the degree of freedom of $[P_{1}({\cal T}_{h})]^{2}$ is $36$ while the degree of freedom of $V_{h}^{0}$ is $39$. Thus, there exist $u_{h}\neq 0$ so that $\nabla_{w}u_{h}=0$. Thus Algorithm $1$ will not work.

Example 2.2**.**

Let $\Omega=(0,1)\times(0,1)$ and is partitioned as shown in Figure 2.

Let $j=k=1=j_{0}-1$. Note that the degree of freedom of $[P_{1}({\cal T}_{h})]^{2}$ is $72$ while the degree of freedom of $V_{h}^{0}$ is $75$. A similar argument as in example 2.1 shows that Algorithm $1$ will not work.

For the completeness, we list the following from [1].

Theorem 3**.**

Let $u_{h}\in V_{h}$ be the weak Galerkin finite element solution of (8). Assume that the exact solution $u\in H^{k+1}(\Omega)$ and $H^{2}$-regularity holds true. Then, there exists a constant $C$ such that

[TABLE]

3 Numerical Experiments

We are devoting this section to verify our theoretical results in previous sections by two numerical examples. The domain in all examples is $\Omega=(0,1)\times(0,1)$. We will implement the SFWG finite element method (8) on triangular meshes. A uniform triangulation of the domain $\Omega$ is used. By refining each triangle for $N=2,4,8,16,32,64,128$, we obtain a sequence of partitions. The first two levels of grids are shown in Figure 3. We apply the SF-WG finite element method with $(P_{k}(T),P_{k}(e),[P_{k+1}(T)]^{2}),k=1,2$ finite element space to find SFWG solution $u_{h}=\{u_{0},u_{b}\}$ in the computation.

Example 3.1**.**

In this example, we solve the Poisson problem (1)-(2) posed on the unit square $\Omega$, and the analytic solution is

[TABLE]

The boundary conditions and the source term $f(x,y)$ are computed accordingly. Table 1 lists errors and convergence rates in $H^{1}$-norm and $L^{2}$-norm.

Example 3.2**.**

In this example, we consider the problem $-\Delta u=f$ with boundary condition and the exact solution is

[TABLE]

The results obtained in Table 2 show the errors of SFWG scheme and the convergence rates in the ${|||}\cdot{|||}$ norm and $\|\cdot\|$ norm. The simulations are conducted on triangular meshes and polynomials of order $k=1,2$. The SFWG scheme with $P_{k}$ element has optimal convergence rate of $O(h^{k})$ in $H^{1}$-norm and $O(h^{k+1})$ in $L^{2}$-norm.

Overall, the numerical experiment’s results shown in the formentioned examples match the theoretical study part in the previous sections.