New error estimates of linear triangle finite elements for the Steklov eigenvalue problem

TL;DR

This paper develops new, optimal error estimates for linear finite element methods applied to the Steklov eigenvalue problem on concave polygons, improving existing theoretical results with supporting numerical experiments.

Contribution

It introduces a novel proof approach leveraging regularity estimates and edge average interpolation, providing sharper error bounds for both conforming and nonconforming finite elements.

Findings

Optimal error estimates in boundary norm for eigenfunctions

Improved theoretical bounds over previous results

Numerical experiments confirm theoretical predictions

Abstract

In this paper we make a further discussion on the finite elements approximation for the Steklov eigenvalue problem on concave polygonal domain. We make full use of the regularity estimate and the characteristic of edge average interpolation operator of nonconforming Crouzeix-Raviart element, which is different from the existing proof argument, and prove a new and optimal error estimate in $\|\cdot\|_{0,\partial\Omega}$ for the eigenfunction of linear conforming finite element and the nonconforming Crouzeix-Raviart element, which is an improvement of the current results. Finally, we present some numerical experiments to support the theoretical analysis.

Equations

$$-div(\alpha \nabla u)+\beta u=0in\Omega ,\alpha \frac{\partial u}{\partial n}=\lambda uon\partial \Omega ,$$

$$a(u,v)=\lambda b(u,v),\forall v\in H^1(\Omega ),$$

$$a(u,v)=\underset{\Omega }{\int }(\alpha \nabla u\cdot \nabla v+\beta uv)dx,b(u,v)=\underset{\partial \Omega }{\int }uvds.$$

$$a(w,v)=b(f,v),\forall v\in H^1(\Omega ).$$

$$\Vert w{\Vert }_{1+\frac{r}{2}}\le C_{\Omega }\Vert f{\Vert }_{0,\partial \Omega };$$

$$\Vert w{\Vert }_{1+r}\le C_{\Omega }\Vert f{\Vert }_{\frac{1}{2},\partial \Omega };$$

$$\Vert w{\Vert }_{\frac{3}{2}+\epsilon }\le C_{\Omega }\Vert f{\Vert }_{\epsilon ,\partial \Omega }.$$

$$a(Af,v)$$

$$Tf=(Af{)}^{\prime },$$

$$Tw=\mu w.$$

$$a_h(u_h,v)={\lambda }_{h}b(u_h,v),\forall v\in S_{nc}^h,$$

$$a_h(u_h,v)=\sum _{K\in {\pi }_h}\underset{K}{\int }(\alpha \nabla u_h\cdot \nabla v+\beta u_{h}v)dx.$$

$$a_h(w_h,v)=b(f,v),\forall v\in S_{nc}^h.$$

$$E_h(w,v)=a_h(w,v)-b(f,v).$$

$$E_h(w,v)\le C{h}^t\Vert w{\Vert }_{1+t}\Vert v{\Vert }_h,\forall v\in S_{nc}^h+H^1(\Omega ),$$

$$\Vert w-w_h{\Vert }_h\le C{h}^t\Vert w{\Vert }_{1+t},$$

$$\Vert w-w_h{\Vert }_{0,\partial \Omega }\le C{h}^{t+\frac{r}{2}}\Vert w{\Vert }_{1+t}.$$

$${\int }_l{I}_{h}uds={\int }_{l}uds\forall l,\forall u\in H^1(\Omega ),$$

$$\Vert u-I_{h}u{\Vert }_0$$

$$\Vert u-I_{h}u{\Vert }_h$$

$$\Vert w-w_h{\Vert }_{-\epsilon ,\partial \Omega }\le C{h}^{t+\frac{1}{2}+\epsilon }\Vert w{\Vert }_{1+t}$$

$$a(v,\phi )=b(g,v),\forall v\in H^1(\Omega ).$$

$$b(f,\phi -{\phi }_h)=b(f,\phi )-b(f,{\phi }_h)$$

$$=a_h(w,\phi )-a_h(w_h,{\phi }_h)$$

$$=a_h(w,\phi )-a_h(w_h,\phi )+a_h(w_h,\phi )-a_h(w_h,{\phi }_h)$$

$$=a_h(w-w_h,\phi )+a_h(w_h-w,\phi -{\phi }_h)+a_h(w,\phi -{\phi }_h).$$

$$a_h(w-w_h,\phi )+a_h(w,\phi -{\phi }_h)$$

$$=a_h(w-w_h,\phi )-b(w-w_h,g)+b(w-w_h,g)$$

$$+a_h(w,\phi -{\phi }_h)-b(f,\phi -{\phi }_h)+b(f,\phi -{\phi }_h)$$

$$=E_h(\phi ,w-w_h)+b(w-w_h,g)+E_h(w,\phi -{\phi }_h)+b(f,\phi -{\phi }_h).$$

$$b(w-w_h,g)=-E_h(\phi ,w-w_h)-E_h(w,\phi -{\phi }_h)-a_h(w_h-w,\phi -{\phi }_h),$$

$$|b(g,w-w_h)|=|b(w-w_h,g)|$$

$$\le |E_h(\phi ,w-w_h)|+|E_h(w,\phi -{\phi }_h)|+|a_h(w_h-w,\phi -{\phi }_h)|.$$

$$|E_h(\phi ,w-w_h)|\le C{h}^{\frac{1}{2}+\epsilon }\Vert \phi {\Vert }_{\frac{3}{2}+\epsilon }\Vert w-w_h{\Vert }_h\le C{h}^{t+\frac{1}{2}+\epsilon }\Vert w{\Vert }_{1+t}\Vert g{\Vert }_{\epsilon ,\partial \Omega },$$

$$|E_h(w,\phi -{\phi }_h)|\le C{h}^t\Vert w{\Vert }_{1+t}\Vert \phi -{\phi }_h{\Vert }_h\le C{h}^{t+\frac{1}{2}+\epsilon }\Vert w{\Vert }_{1+t}\Vert g{\Vert }_{\epsilon ,\partial \Omega },$$

$$|a_h(w_h-w,\phi -{\phi }_h)|\le C\Vert w-w_h{\Vert }_h\Vert \phi -{\phi }_h{\Vert }_h\le C{h}^{t+\frac{1}{2}+\epsilon }\Vert w{\Vert }_{1+t}\Vert g{\Vert }_{\epsilon ,\partial \Omega }.$$

Full text

Key words. Steklov eigenvalue problem, Concave polygonal domain, Linear conforming finite element, Nonconforming Crouzeix-Raviart element, Error estimates.

1 Introduction

Steklov eigenvalue problems have important physical background and many applications. For instance, they appear in the analysis of stability of mechanical oscillators immersed in a viscous fluid (see [12] and the references therein), in the study of surface waves (see [7]), in the study of the vibration modes of a structure in contact with an incompressible fluid (see [6]), in the analysis of the antiplane shearing on a system of collinear faults under slip-dependent friction law (see [10]), etc. Thus the numerical methods for solving these problems have attracted more and more scholars’ attention. Till now, systematical and profound studies on the conforming finite elements approximation for Steklov eigenvalue problems have been made on polygonal domain such as [2, 3, 4, 6, 9, 15, 16, 20, 21]). Recently, the nonconforming finite elements for Steklov problems have also been considered, e.g., see [1, 8, 17, 18, 22]. The aim of this paper is to discuss the error estimates of linear triangle finite elements, including the linear conforming finite element and the nonconforming Crouzeix-Raviart element, approximation for Steklov eigenvalue problems with variable coefficients on concave polygonal domain.

We consider the following Steklov eigenvalue problem

[TABLE]

where $\Omega\subset R^{2}$ is a polygonal domain with $\omega$ being the largest inner angle of $\Omega$, and $\frac{\partial u}{\partial n}$ is the outward normal derivative.

Having in mind that $H^{s}(\Omega)$ denotes the Sobolev space with real order $s$ on $\Omega$, $\|\cdot\|_{s}$ is the norm on $H^{s}(\Omega)$ and $H^{0}(\Omega)=L_{2}(\Omega)$, and $H^{s}(\partial\Omega)$ denotes the Sobolev space with real order $s$ on $\partial\Omega$ with the norm $\|\cdot\|_{s,\partial\Omega}$.

Suppose that the coefficients $\alpha=\alpha(x)$ and $\beta=\beta(x)$ are bounded by above and below by positive constants. We assume that $\alpha\in C^{1}(\bar{\Omega})$.

The weak form of (1.1) is given by: Find $\lambda\in R$, $u\in H^{1}(\Omega),\|u\|_{0,\partial\Omega}=1$, such that

[TABLE]

where

[TABLE]

It is easy to know that $a(\cdot,\cdot)$ is a symmetric, continuous and $H^{1}(\Omega)$-elliptic bilinear form on $H^{1}(\Omega)\times H^{1}(\Omega)$.

In the existing literatures, the error estimate of linear triangle elements eigenfunction, including conforming element and nonconforming Crouzeix-Raviart element (hereafter termed C-R element for simplicity), in $\|\cdot\|_{0,\partial\Omega}$ is all $O(h^{r+\frac{r}{2}})$ where $r$ is the regularity exponent of the eigenfunction (see Lemma 2.1). It is obvious that this estimate is not optimal since it doesn’t achieve the order of interpolation error. In this paper, we improve this estimate when eigenfunctions are singular (i.e., $r<1$) and prove that in this case the error estimate of linear triangle elements eigenfunction can achieve $O(h^{r+\frac{1}{2}})$. Comparing the proof arguments of existing estimates (e.g., see [3, 9, 17, 22]), we make full use of the regularity estimate and the characteristic of edge average interpolation operator of C-R element, especially in the analysis for conforming finite elements, and obtain the improved error estimates (2.25) and (3.6) which are optimal.

Throughout this paper, $C$ denotes a positive constant independent of $h$, which may not be the same constant in different places.

2 The nonconforming Crouzeix-Raviart element approximation for the Steklov eigenvalue problem

Consider the source problem (2.1) associated with (1.1): Find $w\in H^{1}(\Omega)$, such that

[TABLE]

As for the source problem (2.1), there hold the following regularity results.

Lemma 2.1**.**

If $f\in L_{2}(\partial\Omega)$, then $w\in H^{1+\frac{r}{2}}(\Omega)$ and

[TABLE]

if $f\in H^{\frac{1}{2}}(\partial\Omega)$, then $w\in H^{1+r}(\Omega)$ and

[TABLE]

if $f\in H^{\varepsilon}(\partial\Omega)$, $\varepsilon\in(0,r-1/2)$, then $w\in H^{\frac{3}{2}+\varepsilon}(\Omega)$ and

[TABLE]

Here $r=1$ when $\omega<\pi$, and $r<\frac{\pi}{\omega}$ which can be arbitrarily close to $\frac{\pi}{\omega}$ when $\omega>\pi$, and $C_{\Omega}$ is a priori constant.

Proof. See [14].    $\Box$

Note that $a(\cdot,\cdot)$ is coercive, using the source problem (2.1) associated with (1.2) we can define the operator $A:L_{2}(\partial\Omega)\to H^{\frac{3}{2}}(\Omega)\subset H^{1}(\Omega)$, satisfying

[TABLE]

Define the operator $T:L_{2}(\partial\Omega)\to H^{1}(\partial\Omega)$ satisfying

[TABLE]

where ′ denotes the restriction to $\partial\Omega$.

Bramble and Osborn [9] proved that (1.2) has the operator form:

[TABLE]

Namely, if $(\mu,w)$ is an eigenpair of (2.5), then $(\lambda,Aw)$ is an eigenpair of (1.2), $\lambda=\frac{1}{\mu}$. Conversely, if $(\lambda,u)$ is an eigenpair of (1.2), then $(\mu,u^{\prime})$ is an eigenpair of (2.5), $\mu=\frac{1}{\lambda}$.

Let $\lambda$ be the $j$-th eigenvalue of $T$. We arrange eigenvalues by the increasing order with each eigenvalue counted according to its algebraic multiplicity. And let $M(\lambda)$ denote the space spanned by eigenfunctions of (1.2) corresponding to the eigenvalue $\lambda$.

Let $\pi_{h}=\{K\}$ be a regular triangulation of $\Omega$ in the sense of the minimal internal angle condition (see [11], pp. 131). We denote $h=\max_{K\in\pi_{h}}h_{K}$ where $h_{K}$ is the diameter of element $K$. Let $S^{h}_{nc}$ be the C-R element space (see [13]) defined on $\pi_{h}$:

$S^{h}_{nc}=\{v\in L_{2}(\Omega):v\mid_{K}\in span\{1,x_{1},x_{2}\}$, $v$ is continuous at the midpoints of the edges of elements$\}$.

The C-R element approximation of (1.2) is: Find $\lambda_{h}\in R$, $u_{h}\in S^{h}_{nc},\|u_{h}\|_{0,\partial\Omega}=1$, such that

[TABLE]

where

[TABLE]

Define $\|v\|_{h}=(\sum\limits_{K\in\pi_{h}}\|v\|_{1,K}^{2})^{\frac{1}{2}}$, $\|v\|_{1,K}^{2}=\int_{K}(|\frac{\partial v}{\partial x_{1}}|^{2}+|\frac{\partial v}{\partial x_{2}}|^{2}+|v|^{2})dx$. Evidently, $\|\cdot\|_{h}$ is the norm on $S^{h}_{nc}$ and it is simple to show that $a_{h}(\cdot,\cdot)$ is uniformly $S^{h}_{nc}$-elliptic.

The C-R element approximation of (2.1) is: Find $w_{h}\in S^{h}_{nc}$, such that

[TABLE]

Denote the consistency term of the C-R element by

[TABLE]

And based on the standard method (see, for example [1, 8, 17]), the following consistency error estimate can be proved.

Lemma 2.2**.**

Suppose that $w\in H^{1+t}(\Omega)$ with $t\in[\frac{r}{2},1]$ is the weak solution of (2.1), then

[TABLE]

Proof. See, e.g., Lemma 2.2 in [8], or Theorem 2.1 in [17].    $\Box$

Define the interpolation operator $I_{h}:H^{1}(\Omega)\to S^{h}_{nc}$:

[TABLE]

where $l$ is an edge of arbitrary element in $\pi_{h}$.

According to the interpolation theory (see [11]), we have

[TABLE]

Theorem 2.3**.**

Let $w\in H^{1+t}(\Omega)$ with $t\in(\frac{1}{2},1]$ be the weak solution of (2.1), then

[TABLE]

where $t=\frac{1}{2}+\varepsilon$ if $f\in H^{\varepsilon}(\partial\Omega)$, $0<\varepsilon<r-\frac{1}{2}$, and $t=r$ if $f\in H^{\frac{1}{2}}(\partial\Omega)$.

Proof. For each $g\in H^{\varepsilon}(\partial\Omega)$, let $\varphi$ be the unique solution of the following variational problem:

[TABLE]

From (2.4) we know that $\varphi\in H^{\frac{3}{2}+\varepsilon}(\Omega)$. Let $\varphi_{h}=I_{h}\varphi\in S^{h}$ be the interpolation of $\varphi$, then

[TABLE]

By the definition of consistency term we have

[TABLE]

Combining the above two relationships, we get

[TABLE]

then

[TABLE]

From (2.9), (2.3), (2.4), (2.10) and the error estimate of interpolation, we can deduce that

[TABLE]

And substituting the above three estimates into (2), we obtain

[TABLE]

By the definition of negative norm, we have

[TABLE]

namely, (2.15) is true.   $\Box$

Let us denote by $\delta S^{h}_{nc}$ the functions defined on $\partial\Omega$, which are restriction of functions in $S^{h}_{nc}$ to $\partial\Omega$. From [1], pp.189 we know that

[TABLE]

Since $a_{h}(\cdot,\cdot)$ is uniformly elliptic with respect to $h$, the approximate source problem (2.7) associated with (2.6) is uniquely solvable. Thus, we can define the operator $A_{h}:L_{2}(\partial\Omega)\to S^{h}_{nc}$, satisfying

[TABLE]

Define $T_{h}:L_{2}(\partial\Omega)\to\delta S^{h}_{nc}\subset H^{\varepsilon}(\partial\Omega)\subset L_{2}(\partial\Omega)$, satisfying

[TABLE]

[22] proved that (2.6) has the operator form:

[TABLE]

Namely, if $(\mu_{h},w_{h})$ is an eigenpair of (2.19), then $(\lambda_{h},A_{h}w_{h})$ is an eigenpair of (2.6), $\lambda_{h}=\frac{1}{\mu_{h}}$. Conversely, if $(\lambda_{h},u_{h})$ is an eigenpair of (2.6), then $(\mu_{h},u_{h}^{\prime})$ is an eigenpair of (2.19), $\mu_{h}=\frac{1}{\lambda_{h}}$.

We prove the following interpolation estimates.

Lemma 2.4**.**

Let $u\in H^{1+r}(\Omega)$, then the following estimates hold:

[TABLE]

Proof. Let $l\subset\partial\Omega$ be the edge of the element $K$, then by the trace inequality ( see Lemma 7.1.1 in [19]) we have

[TABLE]

thus (2.20) is valid.

For any $g\in H^{\varepsilon}(\partial\Omega)~{}(\varepsilon\in(0,r-\frac{1}{2}))$, let $I_{0}g$ be the piecewise constant interpolation of $g$ on $\partial\Omega$. From the definition of $I_{h}$ and the interpolation estimates we have

[TABLE]

and using the definition of negative norm we know that (2.21) holds.   $\Box$

Lemma 2.5**.**

Suppose that $M(\lambda)\subset H^{1+r}(\Omega)$ and $\lambda$ be the $j$-th eigenvalue of (1.2). Let $\lambda_{h}$ be the $j$-th eigenvalue of (2.6) and $u_{h}$ be an eigenfunction corresponding to $\lambda_{h}$ with $\|u_{h}\|_{0,\partial\Omega}=1$. Then there exists $u\in M(\lambda)$ with $\|u\|_{0,\partial\Omega}=1$, such that

[TABLE]

Proof. See [1, 17, 22].   $\Box$

Lemma 2.5 is an existing conclusion. Next we will improve the estimate (2.24).

Theorem 2.6**.**

Under the conditions of Lemma 2.5, further assume that $\pi_{h}$ is a quasi-uniform mesh (see pp.135 in [11]), then

[TABLE]

Proof. Since $Au$ and $A_{h}u$ are solutions of (2.1) and (2.7) with $f=u$, respectively, then from (2.15) we know that

[TABLE]

Using (2.21) we obtain

[TABLE]

From (2.26) and (2.27), we have

[TABLE]

By the definition of negative norm and the inverse estimates, we have

[TABLE]

thus

[TABLE]

By using (2.29) and (2.20), we get

[TABLE]

It has been proved in [17, 8] that $\|T-T_{h}\|_{0,\partial\Omega}\to 0(h\to 0)$, thus, from Theorem 7.4 in [5] we get

[TABLE]

Substituting (2.30) into (2.31), we obtain (2.25).    $\Box$

Remark 2.1. If $r\in(\frac{1}{2},1)$, i.e., $\Omega$ is concave, it is clear that the estimate (2.25) is better than (2.24).

3 The conforming element approximation for the Steklov eigenvalue problem

Let $S^{h}_{c}\subset C(\overline{\Omega})$ be a space of piecewise linear polynomials defined on $\pi_{h}$. The conforming element approximation of (1.2) is: Find $\lambda_{h}\in R$, $u_{h}\in S^{h}_{c}$ with $\|u_{h}\|_{0,\partial\Omega}=1$, such that

[TABLE]

As for the conforming finite element approximation (3.1), the following results are valid (see [3, 9]).

Lemma 3.1**.**

Let $(\lambda_{h},u_{h})$ be the $j$-th eigenpair of (3.1), $\lambda$ be the $j$-th eigenvalue of (1.2), and $M(\lambda)\subset H^{1+r}(\Omega)$. Then there exists $u\in M(\lambda)$ such that

[TABLE]

where the principle to determine $r$ see Lemma 2.1.

Now, let $P_{h}:H^{1}(\Omega)\to S^{h}_{c}$ be the Ritz projection defined by

[TABLE]

We can define the operator $A_{h}:L_{2}(\partial\Omega)\to S^{h}_{c}$, satisfying

[TABLE]

It is easy to know that $A_{h}=P_{h}A$.

Let $\delta S^{h}_{c}$ be the space of functions defined on $\partial\Omega$, which are restriction of functions in $S^{h}_{c}$ to $\partial\Omega$. Define $T_{h}:L_{2}(\partial\Omega)\to\delta S^{h}_{c}\subset H^{1}(\partial\Omega)$, satisfying

[TABLE]

It has been proved in [5, 9] that $\|T-T_{h}\|_{0,\partial\Omega}\to 0~{}(h\to 0)$, and (3.1) has the operator form:

[TABLE]

Next we will give a new error estimate for the conforming finite element.

Theorem 3.2**.**

Under the conditions of Lemma 3.1, further assume that $\pi_{h}$ is quasi-uniform mesh, then

[TABLE]

Proof. For each $g\in H^{\varepsilon}(\partial\Omega)$, let $\varphi$ be the unique solution of the following variational problem:

[TABLE]

From (2.4) we know that $Ag=\varphi\in H^{\frac{3}{2}+\varepsilon}(\Omega)$, and

[TABLE]

thus, by the definition of negative norm, we have

[TABLE]

Let $I_{h}u$ be the interpolation of $u$ defined by (2.12). By using the inverse estimates, (3.7), (2.21) and (2.20), we get

[TABLE]

By using the spectral approximation theory, we get

[TABLE]

Substituting (3) into (3.9), we obtain (3.6).   $\Box$

Remark 3.1. Comparing (3.4) and (3.6), we can see that when eigenfunctions are singular, i.e., $r<1$, the error estimate in $\|\cdot\|_{0,\partial\Omega}$ is improved.

When we prove the improved estimates (2.25) and (3.6), we make full use of the regularity estimate (2.4) to analyze the negative norm estimate, then use the negative norm estimate and the interpolation of C-R element, especially in the analysis for conforming elements, to obtain the optimal estimates in $L^{2}(\partial\Omega)$; while the existing work is to analyze directly the error in $L^{2}(\partial\Omega)$ by using (2.2) which leads to the lost of error order.

Remark 3.2. We prove the estimates (2.25) and (3.6) under the condition that $\pi_{h}$ is quasi-uniform. In fact, this condition is not a restriction. Since when $\tilde{\pi}_{h}$ is a regular partition derived from $\pi_{h}$ by local refinement, the approximate eigenfunction $\tilde{u}_{h}$ computed on $\tilde{\pi}_{h}$ generally satisfies $\|\tilde{u}_{h}-u\|_{0,\partial\Omega}\leq C\|u_{h}-u\|_{0,\partial\Omega}$, then, for such regular meshes (2.25) and (3.6) are still valid.

4 Numerical Experiments

Consider the problem (1.1), where $\alpha(x)=\beta(x)=1$, $\Omega\subset R^{2}$, $\Omega=[0,1]\times[0,\frac{1}{2}]\bigcup[0,\frac{1}{2}]\times[\frac{1}{2},1]$ is a L-shaped domain with the largest inner angle $\omega=\frac{3}{2}\pi$, or $\Omega=[0,1]\times[0,1]\setminus\{x=(x_{1},x_{2}):\frac{1}{2}\leq x_{1}\leq 1,x_{2}=\frac{1}{2}\}$ is the unit square with a slit which the largest inner angle $\omega=2\pi$.

We adopt a uniform isosceles right triangulation $\pi_{h}$. We use the formula $ratio(\lambda_{h})=lg(\frac{\lambda_{h}-\lambda}{\lambda_{h/2}-\lambda})/lg2$ and $ratio(u_{h})=lg(\frac{\|u_{h}-u\|_{0,\partial\Omega}}{\|u_{h/2}-u\|_{0,\partial\Omega}})/lg2$ to compute the convergence order of approximations of linear conforming element to validate our analysis.

By calculation we find that the eigenfunction associated with $\lambda_{2}$ is singular. So in our numerical experiments we compute the approximation of the second eigenvalue $\lambda_{2,h}$ and the corresponding eigenfunction $u_{2,h}$. Since the exact eigenpairs of the problem (1.1) are unknown, we use the adaptive method to compute a high-precision approximation $\lambda_{2}\in[0.89364476,0.89364690]$ for the L-shaped domain and $\lambda_{2}\in[0.734554376,0.73455822]$ for the unit square with a slit, and use them as the exact values, and the corresponding eigenfunction $u$ is taken as the approximation computed on the uniform mesh with the mesh diameter $h=\frac{\sqrt{2}}{1024}$. The numerical results on the L-shaped domain and the slit domain are listed in Table 1 and Table 2, respectively.

Table 1: The results by using linear conforming element on the L-shaped domain

[TABLE]

Table 2: The results by using linear conforming element on the unit square with a slit

[TABLE]

For the L-shaped domain $r=\frac{2}{3}$, $2r=\frac{4}{3}$. From Table 1 we can see that the convergence order of $\lambda_{2,h}$ is approximately equal to $2r=\frac{4}{3}\approx 1.333333$. It also can be seen from Table 1 that the convergence order of $u_{2,h}$ is very close to $r+\frac{1}{2}=\frac{7}{6}\approx 1.166667$, which is coincide with the theoretical result (3.6); while the convergence order of $u_{2,h}$ according to the previous conclusion (3.4) should be $\frac{3r}{2}=1$.

For the unit square with a slit $r=\frac{1}{2}$. From Table 2 we can see that the convergence order of $\lambda_{2,h}$ is approximately equal to $2r=1$. We can also see from Table 2 that the convergence order of $u_{2,h}$ is very close to $r+\frac{1}{2}=1$, which is coincide with the theoretical result (3.6); while the previous conclusion (3.4) states that the convergence order of $u_{2,h}$ is $\frac{3r}{2}=0.75$.

Acknowledgments. This work was supported by the National Natural Science Foundation of China (Grant Nos. 11201093, 10761003).