© Wuhan University 2023

0 Introduction

Singularly perturbed problems attract much attention in such applications as optimal control, chemical reactions, fluid dynamics and electrical networks^[1]. The exact solution generally displays boundary layers which causes many numerical difficulties to the traditional finite element method on the quasi-uniform mesh. To obtain satisfactory numerical approximations, several numerical strategies were developed, such as the layer-adapted meshes, the fitted operator methods and the stabilised numerical method, see Refs. [1-3] for a survey.

Discontinuous Galerkin (DG) method is a finite element method whose test function and trail function have possible discontinuous points at the element edges^[4]. As such, the discontinuous finite element space provides much flexibility in solving those problems exhibiting large gradients, boundary layers or even discontinuous interfaces, see Ref. [5] for a survey.

For the singularly perturbed problem whose solution exhibits boundary layers, numerical investigations were performed in Ref. [6] for the DG method. Local behavior was explored on the uniform mesh. Uniform convergence and superconvergence were observed on Shishkin mesh. Along this direction, research on the DG method for singularly perturbed problems was developed in a series of papers, see Refs. [7-10].

However, the above error estimates are often performed on Shishkin mesh, which has a simple structure. Owing to the influence of a logarithmic factor, the convergence rate will be deteriorated as the degree of piecewise polynomials goes larger. It is therefore of much interest to derive optimal convergence rate on general graded layer-adapted meshes. Recently, we studied the local DG method on several layer-adapted meshes. Some uniformly optimal error estimates were established for second order convection-diffusion problem^[11,12].

This paper concerns a uniform convergence of the DG method under a much larger range of the layeradapted meshes. These meshes contain five common Shishkin-type meshes as well as three common Bakhvalov-type meshes^[13,14]. Based on a general analysis framework, we establish optimal convergence for the DG method independent of small perturbation parameter. In particular, we establish a sharper $L^{\mathrm{2}}$-norm convergence rate under the situation that the smooth component of the true solution is a piecewise polynomial. Some numerical results are given to confirm our prediction.

We organize this paper as follows. First, layer-adapted meshes and the DG method are introduced in Section 1. Then we present a local projector as well as its approximation error. Our main result follows in Section 2. Finally, we supplement some numerical results to validate our error estimate.

1 Layer-Adapted Meshes and the DG Method

Consider a one-dimensional model problem

$\{\begin{array}{l}-\epsilon q^{\text{'}}+aq=f\mathrm{i}\mathrm{n}\mathrm{\Omega }=(\mathrm{0,1})\\ q(\mathrm{1})=\mathrm{1}\end{array}$(1)

which possesses some essential characteristics of the following second-order problem

$\{\begin{array}{l}-\epsilon u^{″}+a{u}^{\text{'}}+bu=f\mathrm{i}\mathrm{n}\mathrm{\Omega }=(\mathrm{0,1})\\ u^{\text{'}}(\mathrm{1})=\mathrm{1},u(\mathrm{0})=\mathrm{0}\end{array}$(2)

In fact, one can transform (2) into a problem with b=0 if a and b are bothconstants. Then introduce the following transformation

$u^{\text{'}}=q,u(\mathrm{0})=\mathrm{0}$

one obtains the model problem (1).

Assume that $a\ge \alpha >\mathrm{0},$problem (1) exhibits a boundary layer at x =1. Assume that the exact solution of (1) can be expressed as^[15]

$q(x)=q(\mathrm{1}){\mathrm{e}}^{{\int }_{\mathrm{1}}^{x}a\left(s\right)/\epsilon \mathrm{d}s}+{\int }_x^{\mathrm{1}}{\epsilon }^{-\mathrm{1}}f\left(t\right){\mathrm{e}}^{{\int }_t^{x}a\left(s\right)/\epsilon \mathrm{d}s}\mathrm{d}t$

which implies a decomposition $q=\overline{q}+q_{\epsilon }.$ Here $\overline{q}$ is the regular component and $q_{\epsilon }$ is the layer component satisfying

$\left|{\overline{q}}^{\left(j\right)}(x)\right|\le C,\left|q_{\epsilon }^{\left(j\right)}(x)\right|\le C{\epsilon }^{-j}{\mathrm{e}}^{-\alpha (\mathrm{1}-x)/\epsilon },\forall j=\mathrm{0,1},\cdots ,k+\mathrm{1}$(3)

1.1 Layer-Adapted Meshes

Let $\phi$ be a monotonically increasing, continuous and piecewisely differentiable function. Let $\phi (\mathrm{0})=\mathrm{0}$. Introduce a mesh transition parameter

$\tau :=\mathrm{m}\mathrm{i}\mathrm{n}\left\{\frac{\mathrm{1}}{\mathrm{2}},\frac{\sigma \epsilon }{\alpha }\phi \left(\frac{\mathrm{1}}{\mathrm{2}}\right)\right\}$

where $\sigma >\mathrm{0}$ is a constant. Assume that $\tau =\sigma \epsilon \phi (\mathrm{1}/\mathrm{2})/\alpha$ for a small $\epsilon$. Otherwise, the problem is non-singularly perturbed and one can carry out the error analysis in a classical framework. Let $N\ge \mathrm{2}$ be an even integer. Denote $\mathrm{\Omega }={\mathrm{\Omega }}^c\bigcup {\mathrm{\Omega }}^f$, where ${\mathrm{\Omega }}^c=[\mathrm{0,1}-\tau ]$ is the rough region and ${\mathrm{\Omega }}^f=[\mathrm{1}-\tau ,\mathrm{1}]$ is the refined region; they have both equalelements. The mesh points are given by

$x_j=\{\begin{array}{ll}\frac{\mathrm{2}j}{N}(\mathrm{1}-\tau ),& j=\mathrm{0,1},\cdots ,\frac{N}{\mathrm{2}}-\mathrm{1}\\ \mathrm{1}-\frac{\sigma \epsilon }{\alpha }\phi (\mathrm{1}-\frac{j}{N}),& j=\frac{N}{\mathrm{2}},\frac{N}{\mathrm{2}}+\mathrm{1},\cdots ,N\end{array}$(4)

In Table 1, we list eight common layer-adapted meshes, which are simplified as S-, BS-, mBS-,VS-, pS-, B-, mVB- and RS-meshes. Here $\psi ={\mathrm{e}}^{-\phi }$ is mesh characterizing function. See Ref. [13] for more details.

Set ${\mathrm{\Omega }}_N=\{I_j{\}}_{{}_{j=\mathrm{1}}}^N$, where each element $I_j=(x_{j-\mathrm{1}},x_j)$ has the mesh size $h_j=x_j-x_{j-\mathrm{1}}$. Assume that $\epsilon \le N^{-\mathrm{1}}$, then one has $\psi (\mathrm{1}/\mathrm{2})\le N^{-\mathrm{1}}$ for each mesh in Table 1. This property will be frequently used in the following analysis.

Lemma 1^[12] Define

${{G}}_j=\mathrm{m}\mathrm{i}\mathrm{n}\left\{\frac{h_j}{\epsilon },\mathrm{1}\right\}{\mathrm{e}}^{-\alpha (\mathrm{1}-x_j)/\sigma \epsilon },j=N/\mathrm{2}+\mathrm{1},\cdots ,N$

Then one has

$\underset{N/\mathrm{2}+\mathrm{1}\le j\le N}{\mathrm{m}\mathrm{a}\mathrm{x}}{{G}}_j\le C{N}^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}|{\psi }^{\text{'}}|$

${\displaystyle \sum _{j=N/\mathrm{2}+\mathrm{1}}^N}{{G}}_j\le C$

Here $C>\mathrm{0}$ does not depend on $\epsilon$ and N.

1.2 The DG Method

Define a finite element space as

$V_N=\left\{z\in L^{\mathrm{2}}(\mathrm{\Omega }):z{|}_{I_j}\in P^k(I_j),\forall I_j\in {\mathrm{\Omega }}_N\right\}$

where $P^k(I_j)$ is a space of polynomial with degree no larger than k. The functions in the above discontinuous finite element space have possible discontinuous points at the cell ends. Define $z_j^{\pm }=\underset{x\to x_j^{\pm }}{\mathrm{l}\mathrm{i}\mathrm{m}}z(x)$ and the jumps as

$\begin{array}{c}{⟦z⟧}_{\mathrm{0}}=z_{\mathrm{0}}^{+},{⟦z⟧}_j=z_j^{+}-z_j^{-},j=\mathrm{1},\cdots ,N-\mathrm{1},\\ {⟦z⟧}_N=-z_N^{-}.\end{array}$

The DG method reads: Find $Q\in V_N$ such that

${\int }_{I_j}Q(\epsilon r^{\text{'}}+ar)\mathrm{d}x-\epsilon {\widehat{Q}}_j{r}_j^{-}+\epsilon {\widehat{Q}}_{j-\mathrm{1}}{r}_{j-\mathrm{1}}^{+}={\int }_{I_j}fr\mathrm{d}x$(5)

holds for any $r\in V_N$ and $I_j(j=\mathrm{1,2},\cdots ,N)$, where

${\widehat{Q}}_j=\{\begin{array}{l}{Q}_j^{+},j=\mathrm{0,1},\cdots ,N-\mathrm{1}\\ q(\mathrm{1}),j=N\end{array}$

Denote $〈\phi ,\varphi 〉:={\displaystyle \sum _{j=\mathrm{1}}^N}{〈\phi ,\varphi 〉}_{I_j}:={\displaystyle \sum _{j=\mathrm{1}}^N}{\int }_{I_j}\phi \varphi \mathrm{d}x.$ Rewrite the scheme (5) into a compact form: Find $Q\in V_N$ such that

$B(Q;r)=〈f,r〉+\epsilon q(\mathrm{1})r_N^{-},\forall r\in V_N$(6)

where

$B(Q;r)=〈Q,\epsilon r^{\text{'}}+ar〉+\epsilon {\displaystyle \sum _{j=\mathrm{1}}^{N-\mathrm{1}}}{Q}_j^{+}{⟦r⟧}_j+\epsilon Q_{\mathrm{0}}^{+}{r}_{\mathrm{0}}^{+}$

One obtains an energy norm

${\Vert Q\Vert }_{\mathrm{E}}^{\mathrm{2}}:=B(Q;Q)={\Vert a^{\mathrm{1}/\mathrm{2}}Q\Vert }^{\mathrm{2}}+\frac{\epsilon }{\mathrm{2}}{\displaystyle \sum _{j=\mathrm{0}}^N}{⟦Q⟧}_j^{\mathrm{2}}$

Let $f=q(\mathrm{1})=\mathrm{0}$ in (6) then one has $B(Q;Q)={\Vert Q\Vert }_{\mathrm{E}}^{\mathrm{2}}=\mathrm{0}$ and $Q=\mathrm{0},$ which implies the uniquely existence of the computed solution determined by the DG method (6).

2 Convergence Analysis

Divide the error $e=q-Q$ as follows

$e=\left(q-{\pi }^{+}q\right)-\left(Q-{\pi }^{+}q\right):=\eta -\xi$

where ${\pi }^{+}:H^{\mathrm{1}}({\mathrm{\Omega }}_N)\to V_N$ is the local Gauss-Radau projecter such that for any function $z\in H^{\mathrm{1}}({\mathrm{\Omega }}_N)$ and each element $I_j$,

$({\pi }^{+}{z)}_{j-\mathrm{1}}^{+}=z_{j-\mathrm{1}}^{+},{〈{\pi }^{+}z,v〉}_{I_j}={〈z,v〉}_{I_j},\hspace{1em}\forall v\in P^{k-\mathrm{1}}(I_j)$(7)

From Ref.[16], one can verify the well-posedness of the above projection. Furthermore,

${\Vert {\pi }^{+}z\Vert }_{I_j}\le C\left[{\Vert z\Vert }_{I_j}+h_j^{\mathrm{1}/\mathrm{2}}|z_{j-\mathrm{1}}^{+}|\right]$(8)

${\Vert {\pi }^{+}z\Vert }_{L^{\mathrm{\infty }}(I_j)}\le C{\Vert z\Vert }_{L^{\mathrm{\infty }}(I_j)}$(9)

${\Vert z-{\pi }^{+}z\Vert }_{L^{{l}}(I_j)}\le C{h}_j^{k+\mathrm{1}}{\Vert z^{(k+\mathrm{1})}\Vert }_{L^{{l}}(I_j)},{l}=\mathrm{2},\mathrm{\infty }$(10)

Lemma 2   Assume $\sigma \ge k+\mathrm{1.5}$ in the definition of the layer-adapted mesh (4) . Let $z=\overline{z}+z_{\epsilon }$ with each component $\overline{z}$ and $z_{\epsilon }$ satisfying (3), then one has

$\Vert z-{\pi }^{+}z\Vert \le C[N^{-(k+\mathrm{1})}+\sqrt[]{\epsilon }(N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}|{\psi }^{\text{'}}{|)}^{k+\mathrm{1}}]$(11)

${\left({\displaystyle \sum _{j=\mathrm{0}}^N}{⟦z-{\pi }^{+}z⟧}_j^{\mathrm{2}}\right)}^{\mathrm{1}/\mathrm{2}}\le C(N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}|{\psi }^{\text{'}}{|)}^{k+\mathrm{1}/\mathrm{2}}$(12)

Furthermore, if $\overline{z}\in V_N$, one has

$\Vert z-{\pi }^{+}z\Vert \le C\sqrt[]{\epsilon }(N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}|{\psi }^{\text{'}}{|)}^{k+\mathrm{1}}$(13)

Proof   Denote ${\eta }_{\varsigma }=\varsigma -{\pi }^{+}\varsigma$ for $\varsigma =\overline{z},z_{\epsilon }$. From (3) and (10), one has

$\Vert {\eta }_{\overline{z}}\Vert \le C{N}^{-(k+\mathrm{1})}$(14)

For the monotonic increasing function ${\mathrm{e}}^{-\alpha (\mathrm{1}-x)/\epsilon },x\in [\mathrm{0,1}],$ using the stability (8), one has

$\begin{array}{l}{\displaystyle \sum _{j=\mathrm{1}}^{N/\mathrm{2}}}{\Vert {\eta }_{z_{\epsilon }}\Vert }_{I_j}^{\mathrm{2}}\le C{\displaystyle \sum _{j=\mathrm{1}}^{N/\mathrm{2}}}\left[{\Vert z_{\epsilon }\Vert }_{I_j}^{\mathrm{2}}+h_j|z_{\epsilon }(x_{j-\mathrm{1}}){|}^{\mathrm{2}}\right]\\ \le C{\displaystyle \sum _{j=\mathrm{1}}^{N/\mathrm{2}}}\left[{\Vert {\mathrm{e}}^{-\alpha (\mathrm{1}-x)/\epsilon }\Vert }_{I_j}^{\mathrm{2}}+h_j{\mathrm{e}}^{-\mathrm{2}\alpha (\mathrm{1}-x_{j-\mathrm{1}})/\epsilon }\right]\\ \le C{\int }_{\mathrm{0}}^{\mathrm{1}-\tau }{\mathrm{e}}^{-\mathrm{2}\alpha (\mathrm{1}-x)/\epsilon }\mathrm{d}x\\ \le C\epsilon N^{-\mathrm{2}\sigma }\end{array}$(15)

On the refined domain, if $\sigma \ge k+\mathrm{1.5}$, one obtains from Lemma 1 that

$\begin{array}{l}{\displaystyle \sum _{j=N/\mathrm{2}+\mathrm{1}}^N}{\Vert {\eta }_{z_{\epsilon }}\Vert }_{I_j}^{\mathrm{2}}\\ \le C{\displaystyle \sum _{j=N/\mathrm{2}+\mathrm{1}}^N}min\{h_j^{\mathrm{2}(k+\mathrm{1})}{\Vert {z_{\epsilon }}^{(k+\mathrm{1})}\Vert }_{I_j}^{\mathrm{2}},h_j|z_{\epsilon }(x_{j-\mathrm{1}}){|}^{\mathrm{2}}+{\Vert z_{\epsilon }\Vert }_{I_j}^{\mathrm{2}}\}\\ \le C{\displaystyle \sum _{j=N/\mathrm{2}+\mathrm{1}}^N}min\left\{{(\frac{h_j}{\epsilon })}^{\mathrm{2}(k+\mathrm{1})},\mathrm{1}\right\}{\Vert {\mathrm{e}}^{-\alpha (\mathrm{1}-x)/\epsilon }\Vert }_{I_j}^{\mathrm{2}}\\ \le C{\displaystyle \sum _{j=N/\mathrm{2}+\mathrm{1}}^N}\epsilon {\left\{\mathrm{m}\mathrm{i}\mathrm{n}\{\frac{h_j}{\epsilon },\mathrm{1}\}{\mathrm{e}}^{-\alpha (\mathrm{1}-x_j)/\sigma \epsilon }\right\}}^{\mathrm{2}(k+\mathrm{3}/\mathrm{2})}\\ \le C\epsilon \mathrm{m}\mathrm{a}{\mathrm{x}}_{N/\mathrm{2}+\mathrm{1}\le j\le N}{{{G}}_j}^{\mathrm{2}(k+\mathrm{1})}{\displaystyle \sum _{j=N/\mathrm{2}+\mathrm{1}}^N}{{G}}_j\\ \le C\epsilon (N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}|{\psi }^{\text{'}}{|)}^{\mathrm{2}(k+\mathrm{1})}\end{array}$(16)

Consequently, (11) follows from (14)-(16). Note that if $\overline{z}\in V_N$, (13) follows from (15),(16).

By the $L^{\mathrm{\infty }}$ approximation property (10), one has

$\begin{array}{l}{\displaystyle \sum _{j=\mathrm{0}}^N}{⟦{\eta }_{\overline{z}}⟧}_j^{\mathrm{2}}\le C{\displaystyle \sum _{j=\mathrm{1}}^N}{\Vert {\eta }_{\overline{z}}\Vert }_{L^{\mathrm{\infty }}(I_j)}^{\mathrm{2}}\\ \le C{\displaystyle \sum _{j=\mathrm{1}}^N}{N}^{-\mathrm{2}(k+\mathrm{1})}\le C{N}^{-(\mathrm{2}k+\mathrm{1})}\end{array}$

By Lemma 1 and $\sigma \ge k+\mathrm{1}$ one obtains

$\begin{array}{l}{\displaystyle \sum _{j=\mathrm{0}}^N}{⟦{\eta }_{z_{\epsilon }}⟧}_j^{\mathrm{2}}\le C{\displaystyle \sum _{j=N/\mathrm{2}+\mathrm{1}}^N}{\Vert {\eta }_{z_{\epsilon }}\Vert }_{L^{\mathrm{\infty }}(I_j)}^{\mathrm{2}}+C{\displaystyle \sum _{j=\mathrm{1}}^{N/\mathrm{2}}}{\Vert {\eta }_{z_{\epsilon }}\Vert }_{L^{\mathrm{\infty }}(I_j)}^{\mathrm{2}}\\ \le C{\displaystyle \sum _{j=N/\mathrm{2}+\mathrm{1}}^N}min\left\{h_j^{\mathrm{2}(k+\mathrm{1})}{\Vert z_{\epsilon }^{(k+\mathrm{1})}\Vert }_{L^{\mathrm{\infty }}(I_j)}^{\mathrm{2}},{\Vert z_{\epsilon }\Vert }_{L^{\mathrm{\infty }}(I_j)}^{\mathrm{2}}\right\}\\ +C{\displaystyle \sum _{j=\mathrm{1}}^{N/\mathrm{2}}}{\Vert z_{\epsilon }\Vert }_{L^{\mathrm{\infty }}(I_j)}^{\mathrm{2}}\\ \le C\mathrm{m}\mathrm{a}{\mathrm{x}}_{N/\mathrm{2}+\mathrm{1}\le j\le N}{{{G}}_j}^{\mathrm{2}k+\mathrm{1}}{\displaystyle \sum _{j=N/\mathrm{2}+\mathrm{1}}^N}{{G}}_j+C{N}^{-\mathrm{2}\sigma +\mathrm{1}}\\ \le C(N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}|{\psi }^{\text{'}}{|)}^{\mathrm{2}k+\mathrm{1}}\end{array}$

which leads to (12).

Theorem 1   Let q and $Q\in V_N$ are respectively the solutions of (1) and (5) . Then one has

$\Vert q-Q\Vert \le C\left[N^{-\left(k+\mathrm{1}\right)}+\sqrt[]{\epsilon }{\left(N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}\left|{\psi }^{\text{'}}\right|\right)}^{k+\mathrm{1}}\right]$(17)

Furthermore, if $\overline{q}\in V_N$, one has

$\Vert q-Q\Vert \le C\sqrt[]{\epsilon }{\left(N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}\left|{\psi }^{\text{'}}\right|\right)}^{k+\mathrm{1}}$(18)

Proof   For any$r\in V_N$, one has Galerkin orthogonality $B(q-Q;r)=\mathrm{0}$ which leads to

$\begin{array}{l}{\Vert \xi \Vert }_{\mathrm{E}}^{\mathrm{2}}=B\left(\xi ;\xi \right)=B\left(\eta ;\xi \right)\\ =\epsilon {\displaystyle \sum _{j=\mathrm{0}}^{N-\mathrm{1}}}{\eta }_j^{+}{⟦\xi ⟧}_j+〈\eta ,\epsilon {\xi }^{\text{'}}+a\xi 〉\\ =〈\eta ,a\xi 〉\le C\Vert \eta \Vert \Vert a^{\mathrm{1}/\mathrm{2}}\xi \Vert \\ \le C\Vert \eta \Vert {\Vert \xi \Vert }_{\mathrm{E}}\end{array}$

Then one gets from (11) that

$\begin{array}{l}\Vert \xi \Vert \le {\Vert \xi \Vert }_{\mathrm{E}}\le C\Vert \eta \Vert \\ \le C\left[N^{-(k+\mathrm{1})}+\sqrt[]{\epsilon }{\left(N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}\left|{\psi }^{\text{'}}\right|\right)}^{k+\mathrm{1}}\right]\end{array}$

Hence,

$\begin{array}{l}\Vert e\Vert \le \Vert \eta \Vert +\Vert \xi \Vert \le C\Vert \eta \Vert \\ \le C\left[N^{-(k+\mathrm{1})}+\sqrt[]{\epsilon }{\left(N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}\left|{\psi }^{\text{'}}\right|\right)}^{k+\mathrm{1}}\right]\end{array}$

If $\overline{q}\in V_N$, by (13), one has

$\Vert e\Vert \le C\Vert \eta \Vert \le C\sqrt[]{\epsilon }(N^{-\mathrm{1}}{\mathrm{m}\mathrm{a}\mathrm{x}\left|{\psi }^{\text{'}}\right|)}^{k+\mathrm{1}}$

Remark 1   Theorem 1 presents a convergence rate $O(N^{-\left(k+\mathrm{1}\right)}+\sqrt[]{\epsilon }{\left(N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}\left|{\psi }^{\text{'}}\right|\right)}^{k+\mathrm{1}})$ for an equal number of mesh elements in the rough and refined domains. However, it is possible to explore different number of mesh elements in rough and refined domains (denote $N_{\mathrm{1}}$ and $N_{\mathrm{2}}$ respectively). Following the similar line, one derives a convergence rate $O({N_{\mathrm{1}}}^{-\left(k+\mathrm{1}\right)}+\sqrt[]{\epsilon }{\left({N_{\mathrm{2}}}^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}\left|{\psi }^{\text{'}}\right|\right)}^{k+\mathrm{1}})$. Then it is possible to use larger $N_{\mathrm{2}}$ to balance the influence of $\mathrm{m}\mathrm{a}\mathrm{x}\left|{\psi }^{\text{'}}\right|$ on the convergence rate and arrive at a convergence $O({N_{\mathrm{1}}}^{-\left(k+\mathrm{1}\right)})$.

3 Numerical Experiments

We perform the DG method (5) on the eight meshes given by Table 1. Set $\sigma =k+\mathrm{1.5}.$ For the pS-mesh, take $m=\mathrm{2}.$ For the RS-mesh, take $l=\mathrm{3}.$ Compute the convergence rate by the formulae

$\begin{array}{l}{r}_{\mathrm{2}}=\frac{\mathrm{l}\mathrm{o}\mathrm{g}{e}^N-\mathrm{l}\mathrm{o}\mathrm{g}{e}^{\mathrm{2}N}}{\mathrm{l}\mathrm{o}\mathrm{g}\mathrm{2}},\\ r_s=\frac{\mathrm{l}\mathrm{o}\mathrm{g}{e}^N-\mathrm{l}\mathrm{o}\mathrm{g}{e}^{\mathrm{2}N}}{\mathrm{l}\mathrm{o}\mathrm{g}(\mathrm{2}\mathrm{l}\mathrm{n}N/\mathrm{l}\mathrm{n}\mathrm{2}N)}\end{array}$

The quantities $r_{\mathrm{2}}$ and $r_s$ are used to reflect the convergence rates from the error bounds of the forms $C{N}^{-r}$and $C(N^{-\mathrm{1}}{\mathrm{l}\mathrm{n}N)}^r$, respectively.

Example 1 Consider the problem (1) with $b=\mathrm{1}$ and the true solution is set as

$q(x)=\mathrm{c}\mathrm{o}\mathrm{s}(\mathrm{1}-x){\mathrm{e}}^{-(\mathrm{1}-x)/\epsilon }$

As mentioned before, for a large $\epsilon$, one can use arbitrary mesh with maximum mesh size h and expect a convergence rate $O(h^{k+\mathrm{1}})$. In Table 2, we list the numerical results on four meshes. One observes that the convergence rates on uniform mesh and B-mesh are both $O(N^{-(k+\mathrm{1})})$, the convergence rate on S-mesh is $O((N^{-\mathrm{1}}{\mathrm{l}\mathrm{n}N)}^{k+\mathrm{1}})$, while the convergence rate on BS-mesh is $O(\mathrm{m}\mathrm{a}\mathrm{x}\{\epsilon ,N^{-\mathrm{1}}{\}}^{k+\mathrm{1}})$ because its maximum mesh size is bounded by $\mathrm{m}\mathrm{a}\mathrm{x}\{\epsilon ,N^{-\mathrm{1}}\}.$

Table 3 and Table 4 list the $L^{\mathrm{2}}$-error of the DG method on the eight meshes for $\epsilon ={\mathrm{10}}^{-\mathrm{4}}$ and $\epsilon ={\mathrm{10}}^{-\mathrm{8}}$, respectively. One observes for the S-mesh a convergence rate $O((N^{-\mathrm{1}}{\mathrm{l}\mathrm{n}N)}^{k+\mathrm{1}})$, while for the BS-, mBS-, VS-, B- and mVB- meshes a general convergence rate $O(N^{-(k+\mathrm{1})})$. For the pS-mesh, the convergence rate is a little smaller than k+1 because of the influence of the logarithmic factor (here $\mathrm{m}\mathrm{a}\mathrm{x}\left|{\psi }^{\text{'}}\right|={\left(\mathrm{l}\mathrm{n}N\right)}^{\mathrm{1}/\mathrm{2}}$). For the RS-mesh, the convergence rate behaviors as $O(N^{-\mathrm{2}(k+\mathrm{1})/\mathrm{3}})$ which agrees with theoretical prediction (17) in view of l=3 and $\mathrm{m}\mathrm{a}\mathrm{x}|{\psi }^{\text{'}}|=N^{\mathrm{1}/\mathrm{3}}.$ Thus, the numerical convergence rate is generally $O((N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}|{\psi }^{\text{'}}{|)}^{k+\mathrm{1}})$ for every meshes, which confirms our theoretical convergence rate in Theorem 1.

Furthermore, Figure 1 and Table 5 show that the $L^{\mathrm{2}}$-error is uniform regarding the singular perturbation parameter in the case that $k=\mathrm{2}$ and N=64 except the S-mesh and RS-mesh which are slightly influenced by the factor $\sqrt[]{\epsilon }$ as $\epsilon$ is suitably large.

Fig. 1 ${\mathit{L}}^{\mathrm{2}}$-error on the singular perturbation parameter for Example 1

Example 2 We continue example 1 but employ a different true solution

$q(x)=\mathrm{2}-{\mathrm{e}}^{-(\mathrm{1}-x)/\epsilon }$

Now the regular component of q is in the finite element space. Let $\epsilon ={\mathrm{10}}^{-\mathrm{4}}$,${\mathrm{10}}^{-\mathrm{8}}$. Tables 6-7 show a general convergence rate $O((N^{-\mathrm{1}}\mathrm{m}\mathrm{a}\mathrm{x}|{\psi }^{\text{'}}{|)}^{k+\mathrm{1}})$ for each mesh. Furthermore, Figure 2 and Table 8 demonstrate the influence of the factor $\sqrt[]{\epsilon }$ on the upper bound of these $L^{\mathrm{2}}$-errors, which confirms our prediction (18).

Fig. 2 ${\mathit{L}}^{\mathrm{2}}$-error on the singular perturbation parameter for Example 2

References

All Tables

All Figures

	Fig. 1 ${\mathit{L}}^{\mathrm{2}}$-error on the singular perturbation parameter for Example 1
In the text

	Fig. 2 ${\mathit{L}}^{\mathrm{2}}$-error on the singular perturbation parameter for Example 2
In the text