On the exponent of exponential convergence of the $p$-version FEM spaces

TL;DR

This paper demonstrates that certain non-standard finite element spaces, like serendipity and total degree basis methods, achieve faster exponential convergence rates than traditional tensor product bases for analytic problems, supported by new error bounds and numerical tests.

Contribution

It introduces new $p$-optimal error bounds for specific finite element bases, showing improved exponential convergence rates and revealing that some basis functions are unnecessary for optimal error bounds.

Findings

Serendipity and total degree basis methods have faster exponential convergence than tensor product bases.

New $p$-optimal error bounds are established for projections onto these bases.

Numerical examples confirm the theoretical convergence rates.

Abstract

We study the exponent of the exponential rate of convergence in terms of the number of degrees of freedom for various non-standard {$p$-version} finite element spaces employing reduced cardinality basis. More specifically, we show that serendipity finite element methods and discontinuous Galerkin finite element methods with total degree $\mathcal{P}_p$ basis have a faster exponential convergence with respect to the number of degrees of freedom than their counterparts employing the tensor product $\mathcal{Q}_p$ basis for quadrilateral/hexahedral elements, for piecewise analytic problems under $p$-refinement. The above results are proven by using a new $p$-optimal error bound for the $L^2$-orthogonal projection onto the total degree $\mathcal{P}_p$ basis, and for the $H^1$-projection onto the serendipity finite element space over tensor product elements with dimension $d\geq2$. These new $p$-optimal error bounds lead to a larger exponent of the exponential rate of convergence with respect to the number of degrees of freedom. Moreover, these results show that part of the basis functions in $\mathcal{Q}_p$ basis {plays} no roles in achieving the $hp$-optimal error bound in the Sobolev space. The sharpness of theoretical results is also verified by a series of numerical examples.

Tables

Table 1. Table 1: Example 3: Convergence rate in p 𝑝 p of the FEM(S) and FEM(Q) in H 1 superscript 𝐻 1 H^{1} –seminorm.

$p$	FEM(S)		FEM(Q)		Ratio of Error
	${\left|u-u_h\right|}_{H^1\left(\Omega \right)}$	$p$-rate	${\left|u-u_h\right|}_{H^1\left(\Omega \right)}$	$p$-rate	FEM(S)/FEM(Q)
1	2.09E-01		2.09E-01		1
2	1.25E-01	0.7386	9.62E-02	1.1204	1.303
3	1.20E-01	0.1096	5.99E-02	1.1691	2.0023
4	9.00E-02	0.9971	4.23E-02	1.2087	2.128
5	6.93E-02	1.1703	3.21E-02	1.2372	2.16
10	2.96E-02	1.261	1.32E-02	1.2968	2.2311
15	1.76E-02	1.2921	7.79E-03	1.3143	2.2558
20	1.21E-02	1.306	5.33E-03	1.3215	2.2675
25	9.03E-03	1.3135	3.97E-03	1.3251	2.2741
30	7.10E-03	1.3181	3.12E-03	1.3272	2.2783
35	5.79E-03	1.3211	2.54E-03	1.3285	2.2811
40	4.86E-03	1.3232	2.13E-03	1.3294	2.2832
45	4.15E-03	1.3247	1.82E-03	1.33	2.2847
50	3.61E-03	1.3259	1.58E-03	1.3305	2.2859
55	3.18E-03	1.3268	1.39E-03	1.3309	2.2868
60	2.84E-03	1.3275	1.24E-03	1.3312	2.2876

Equations

$$\Vert u-{\Pi }_{{\mathcal{Q}}_p}u{\Vert }_{L^2(\hat{\kappa })}^2\le C(s)(p+1{)}^{-2s}|u{|}_{H^s(\hat{\kappa })}^2,$$

$$\Vert u-{\Pi }_{{\mathcal{P}}_p}u{\Vert }_{L^2(\hat{\kappa })}^2=\Vert u-{\Pi }_{{\mathcal{Q}}_{\lfloor p/2\rfloor }}u{\Vert }_{L^2(\hat{\kappa })}^2\le \tilde{C}(s)(\lfloor p/2\rfloor +1{)}^{-2s}|u{|}_{H^s(\hat{\kappa })}^2,$$

$$\displaystyle\Phi_{d}(m,n)=\Big{(}\frac{\Gamma(\frac{m-n}{d}+1)}{\Gamma(\frac{m+n}{d}+1)}\Big{)}^{d},$$

$$\Vert u{\Vert }_{V^l(\hat{\kappa })}^2=\sum _{|\alpha |=0}^l|u{|}_{V^l(\hat{\kappa })}^2,\text{and}|u{|}_{V^l(\hat{\kappa })}^2=\sum _{|\alpha |=l}\Vert W^{\alpha }{D}^{\alpha }u{\Vert }_{L^2(\hat{\kappa })}^2.$$

$$a_i={\int }_{\hat{\kappa }}u(x)\prod _{k=1}^d\frac{2i_k+1}{2}{L}_{i_k}(x_k)\mathrm{d}x.$$

$$D^{\alpha }u(x)=\sum _{i\ge \alpha }^{\infty }{a}_i\prod _{k=1}^d{L}_{i_k}^{({\alpha }_k)}(x_k).$$

$${\int }_{-1}^1(1-{\xi }^2{)}^k{L}_i^{(k)}(\xi )L_j^{(k)}(\xi )\mathrm{d}\xi =\frac{2{\delta }_{ij}}{2i+1}\frac{\Gamma (i+k+1)}{\Gamma (i-k+1)},$$

$$\Vert W^{\alpha }{D}^{\alpha }u{\Vert }_{L^2(\hat{\kappa })}^2=\sum _{i\ge \alpha }^{\infty }|a_i{|}^2\prod _{k=1}^d\frac{2}{2i_k+1}\frac{\Gamma (i_k+{\alpha }_k+1)}{\Gamma (i_k-{\alpha }_k+1)}.$$

$$\displaystyle\|{u-\Pi_{{\cal Q}_{p}}u}\|_{L^{2}({\hat{\kappa}})}^{2}\leq\Phi_{1}(p+1,s)\Big{(}\sum_{k=1}^{d}\|{W^{s}_{k}D^{s}_{k}u}\|_{{\hat{\kappa}}}^{2}\Big{)}\leq\Phi_{1}(p+1,s)|u|_{H^{s}(\hat{\kappa})}^{2},$$

$$\sqrt{2\pi }{x}^{x+\frac{1}{2}}{e}^{-x}\le \Gamma (x+1)\le e{x}^{x+\frac{1}{2}}{e}^{-x},$$

$$\displaystyle\Phi_{d}(p+1,s)\leq C(s)\Big{(}\frac{d}{p+1}\Big{)}^{2s},$$

$$F(\xi ,\rho ):=\prod _{k=1}^d\frac{\Gamma ({\rho }_k-{\xi }_k+1)}{\Gamma ({\rho }_k+{\xi }_k+1)}\le {\Phi }_d(M,m).$$

$$\displaystyle\|{u-\Pi_{{\cal P}_{p}}u}\|_{L^{2}({\hat{\kappa}})}^{2}\leq\Phi_{d}(p+1,s)|u|_{V^{s}(\hat{\kappa})}^{2}\leq C(s)\Big{(}\frac{d}{p+1}\Big{)}^{2s}|u|_{H^{s}(\hat{\kappa})}^{2}.$$

$$\Vert u-{\Pi }_{{\mathcal{P}}_p}u{\Vert }_{L^2(\hat{\kappa })}^2$$

$$\displaystyle=\sum_{|\alpha|=s}\hskip 5.69046pt\sum_{|i|=p+1,i\geq\alpha}^{\infty}|a_{i}|^{2}\Big{(}\prod_{k=1}^{d}\frac{2}{2i_{k}+1}\frac{\Gamma(i_{k}+\alpha_{k}+1)}{\Gamma(i_{k}-\alpha_{k}+1)}\Big{)}$$

$$\displaystyle\times\Big{(}\prod_{k=1}^{d}\frac{\Gamma(i_{k}-\alpha_{k}+1)}{\Gamma(i_{k}+\alpha_{k}+1)}\Big{)}$$

$$\le {\Phi }_d(p+1,s)\sum _{|\alpha |=s}\sum _{|i|=p+1,i\ge \alpha }^{\infty }|a_i{|}^2\prod _{k=1}^d\frac{2}{2i_k+1}\frac{\Gamma (i_k+{\alpha }_k+1)}{\Gamma (i_k-{\alpha }_k+1)}$$

$$\le {\Phi }_d(p+1,s)\sum _{|\alpha |=s}\Vert W^{\alpha }{D}^{\alpha }u{\Vert }_{L^2(\hat{\kappa })}^2$$

$$\displaystyle=\Phi_{d}(p+1,s)|u|_{V^{s}(\hat{\kappa})}^{2}\leq C(s)\Big{(}\frac{d}{p+1}\Big{)}^{2s}|u|_{H^{s}(\hat{\kappa})}^{2},$$

$$L(\xi ,\rho ,\mu ,\lambda )=F(\xi ,\rho )+\mu (|\xi |-m)+\lambda (|\rho |-M),$$

$$\frac{\partial L}{\partial {\xi }_k}=-\left(\frac{{\Gamma }^{\prime }({\rho }_k-{\xi }_k+1)}{\Gamma ({\rho }_k-{\xi }_k+1)}+\frac{{\Gamma }^{\prime }({\rho }_k+{\xi }_k+1)}{\Gamma ({\rho }_k+{\xi }_k+1)}\right)F(\xi ,\rho )+\mu =0,$$

$$\frac{\partial L}{\partial {\rho }_k}=\left(\frac{{\Gamma }^{\prime }({\rho }_k-{\xi }_k+1)}{\Gamma ({\rho }_k-{\xi }_k+1)}-\frac{{\Gamma }^{\prime }({\rho }_k+{\xi }_k+1)}{\Gamma ({\rho }_k+{\xi }_k+1)}\right)F(\xi ,\rho )+\lambda =0,$$

$$\frac{{\Gamma }^{\prime }({\rho }_k-{\xi }_k+1)}{\Gamma ({\rho }_k-{\xi }_k+1)}=\frac{\mu -\lambda }{2F(\xi ,\rho )}\text{and}\frac{{\Gamma }^{\prime }({\rho }_k+{\xi }_k+1)}{\Gamma ({\rho }_k+{\xi }_k+1)}=\frac{\mu +\lambda }{2F(\xi ,\rho )},$$

$$\phi(z+1)=-\gamma+\sum_{n=1}^{\infty}\frac{z}{n(n+z)}=-\gamma+\sum_{n=1}^{\infty}\Big{(}\frac{1}{n}-\frac{1}{n+z}\Big{)},\quad z\neq-1,-2,\dots,$$

$$F(\tilde{\xi },\tilde{\rho })={\Phi }_d(M,m).$$

$${\Phi }_n(M,m)\le {\Phi }_d(M,m),n=1,\dots ,d-1.$$

$$\displaystyle\frac{\Phi_{d}(M,m)}{\Phi_{n}(M,m)}\geq\Big{(}\frac{\sqrt{2\pi}}{e}\Big{)}^{d+n}\Big{(}\frac{d}{n}\Big{)}^{2\delta M}\Big{(}\frac{1-\delta}{1+\delta}\Big{)}^{\frac{d-n}{2}}.$$

$$\displaystyle\frac{\Phi_{d}(M,m)}{\Phi_{n}(M,m)}=\frac{(\Gamma(\frac{2M}{n}+1))^{n}}{(\Gamma(\frac{2M}{d}+1))^{d}}\geq\frac{(\sqrt{2\pi})^{n}}{e^{d}}\Big{(}\frac{d}{2M}\Big{)}^{\frac{d-n}{2}}\Big{(}\frac{d}{n}\Big{)}^{2M+\frac{n}{2}}.$$

$$F({\xi }^b,{\rho }^b)={\Phi }_1(M,m)\le {\Phi }_2(M,m)=F(\tilde{\xi },\tilde{\rho }).$$

$$F({\xi }^b,{\rho }^b)={\Phi }_2(M,m)\le {\Phi }_3(M,m)=F(\tilde{\xi },\tilde{\rho }).$$

$${\Phi }_{d-1}(M,m)\le {\Phi }_d(M,m),$$

$${\mathcal{S}}_p(\hat{\kappa }):={\mathcal{P}}_p(\hat{\kappa })+\text{span}\{x_1^p{x}_2,x_1{x}_2^p\},p\ge 1.$$

$${\pi }_{{\mathcal{Q}}_p}u(x_1,x_2):=$$

Full text

1 Introduction

Polynomial approximation on tensor product domains plays an important role in deriving the exponential rate of convergence with respect to the number of degrees of freedom for $hp$-version finite element methods (FEMs) gui1986h ; guo1986hp1 ; guo1986hp2 ; MR946376 ; MR954788 ; MR1301344 ; MR3351174 ; schwab ; melenk1999hp and $hp$-version discontinuous Galerkin finite element methods (DGFEMs) hss ; newpaper ; Wihler2002 ; schotzau2013hp1 ; schotzau2013hp2 ; schotzau2016hp . In general, the proof of the exponential rate of convergence usually depends on the $hp$-approximation results for some suitable projection operators onto a local polynomial space consisting of polynomials with degree less or equal than $p$ in each variable (known as ${\cal Q}_{p}$ basis) over a tensor product element (quadrilateral/hexahedral elements), for dimension $d\geq 2$.

The key reason for using the ${\cal Q}_{p}$ basis over a tensor product element is because $hp$-optimal approximation results for the multi-dimensional projection operators can be derived by using the stability and approximation results of the one-dimensional projections via tensor product arguments. On the other hand, the $hp$-approximation results for $L^{2}$-orthogonal projections onto polynomial basis with total degree less or equal than $p$ (${\cal P}_{p}$ basis) and $H^{1}$-projections onto serendipity basis (${\cal S}_{p}$ basis) have not been fully explored. Typically, $hp$-error bounds for projections onto the ${\cal P}_{p}$ or ${\cal S}_{p}$ basis have been derived using the fact that there exists a $q\leq p$ such that the bases ${\cal P}_{p}$ or ${\cal S}_{p}$ contain ${\cal Q}_{q}$ as a subset, together with the help of the $hp$-optimal approximation results for the projections onto the basis ${\cal Q}_{q}$, see Corollary 4.52 inschwab .

For instance, we consider the $L^{2}$–norm error bound of the two-dimensional $L^{2}$-orthogonal projection $\Pi_{{\cal Q}_{p}}$ onto the ${\cal Q}_{p}$ basis as an example, cf. canuto1982approximation ; newpaper . Let $\hat{\kappa}=(-1,1)^{2}$ and $u\in H^{l}(\hat{\kappa})$, $l$ is an integer with $l\geq 0$. Then, the following estimate holds,

[TABLE]

where the constant $C(s)$ is independent of $p$ and $0\leq s\leq\min\{p+1,l\}$. It is straightforward to see that the above error bound is sharp in the sense that it is $p$-optimal in both Sobolev regularity index $l$ and polynomial approximation order $p$.

Next, we consider the $L^{2}$–norm error bound of $L^{2}$-orthogonal projection $\Pi_{{\cal P}_{p}}$ onto the ${\cal P}_{p}$ basis. Following the Lemma 6 in MR1974174 , we define $\Pi_{{\cal P}_{p}}=\Pi_{{\cal Q}_{\lfloor p/2\rfloor}}$, with $\lfloor p/2\rfloor$ denoting the largest integer which is less than or equal to $p/2$. Then, the following bound holds:

[TABLE]

where the constant $\tilde{C}(s)$ is independent of $p$ and $0\leq s\leq\min\{\lfloor p/2\rfloor+1,l\}$. We emphasize that *for function $u\in H^{l}(\hat{\kappa})$, with $p$ sufficiently large, the above error bound is $p$-optimal because $s=l$. However, if function $u$ is sufficiently smooth or even analytic, then the above error bound is $p$-suboptimal by at least ${p/2}$ orders because $s={\lfloor p/2\rfloor+1}$. * The similar $p$-suboptimal error bound holds for $H^{1}$-projections onto ${\cal S}_{p}$ basis.

Using the $p$-suboptimal error bound for $L^{2}$-orthogonal projections onto the ${\cal P}_{p}$ basis and $H^{1}$-projections onto the ${\cal S}_{p}$ basis, it is possible to derive an exponential rate of convergence for $hp$-FEMs employing the ${\cal S}_{p}$ basis and $hp$-DGFEMs employing the ${\cal P}_{p}$ basis, but the resulting exponent is much smaller with respect to the number of degrees of freedom than the exponent of FEMs and DGFEMs employing the ${\cal Q}_{p}$ basis. This contradicts the numerical observation in work cangiani2013hp ; cangiani2015hp ; cangiani2016hp ; Zhaonan , where it is observed that the error with respect to the number of degrees of freedom for DGFEMs with the ${\cal P}_{p}$ basis on tensor product elements has a steeper exponential convergence compared to DGFEMs with the ${\cal Q}_{p}$ basis, for sufficiently smooth problems. This situation has been numerically tested on many different examples. We also observed numerically that the ratio of the slope of the exponential error decay for DGFEMs with the ${\cal P}_{p}$ basis compared to that of the ${\cal Q}_{p}$ basis depends only on the space dimension. The same phenomenon is also observed when comparing conforming FEMs with the ${\cal S}_{p}$ basis and the ${\cal Q}_{p}$ basis.

The disagreement between the numerical observations and theoretical results implies that the error bound (2) is not a sharp bound for ${\cal P}_{p}$ and ${\cal S}_{p}$ basis. To address this, in this work, we derive an $hp$-optimal error bound for the $L^{2}$-orthogonal projection onto the ${\cal P}_{p}$ basis in the $L^{2}$–norm, and for the $H^{1}$-projection onto the ${\cal S}_{p}$ basis in the $L^{2}$–norm and $H^{1}$–seminorm.

The technique for proving the new error bounds is different from the existing techniques for $hp$-approximation with the ${\cal Q}_{p}$ basis, due to the lack of a tensor product structure in the ${\cal P}_{p}$ and ${\cal S}_{p}$ bases, thereby hindering the use of the usual tensor product arguments. The key tools used in this work are: a multi-dimensional orthogonal polynomial expansion and the careful selection of basis functions. To the best author’s knowledge, the new error bounds for both projections never appeared in the literatures. The resulting bounds are $hp$-optimal with respect to both Sobolev regularity and polynomial approximation order. Moreover, it also shows that the ${\cal Q}_{p}$ basis contains in a sense “extra” basis functions that are unnecessary for optimal convergence. These basis functions do not increase the order in $p$ of the error bound, but instead only reduce its “constant”.

By using the new $hp$-optimal error bound for the $L^{2}$-orthogonal projection onto the ${\cal P}_{p}$ basis and the $H^{1}$-projection onto the ${\cal S}_{p}$ basis, we can prove that methods using ${\cal P}_{p}$ and ${\cal S}_{p}$ bases offer exponential convergence with a larger exponent with respect to the number of degrees of freedom than comparable methods using ${\cal Q}_{p}$ basis for piecewise analytic problem under $p$-refinement. Furthermore, the approximation results also show that there are a lot of basis functions in ${\cal Q}_{p}$ basis with no roles in improving the $hp$-optimal error bound, which can be generalized to other FEM with the local polynomial space employing reduced cardinality basis. Finally, we emphasize that we are using DGFEM employing ${\cal P}_{p}$ basis for quadrilateral and hexahedral elements also and this is the key novelty of the approach, since this is possible for DGFEM and essentially for serendipity spaces.

The remainder of this work is structured as follows. In Section 2, we introduce the required notation and the weighted Sobolev spaces together with some properties about the orthogonal polynomials. Then, the $p$-optimal error bound for the $L^{2}$-orthogonal projection onto the ${\cal P}_{p}$ basis in $L^{2}$–norm is proved in Section 3. In Section 4, we derive the $p$-optimal error bound for $H^{1}$-projection onto the ${\cal S}_{p}$ basis in both $L^{2}$–norm and $H^{1}$–seminorm. Section 5 is devoted to deriving the exponential rate of convergence for the $L^{2}$- and the $H^{1}$-projections employing different local polynomial bases. The sharpness of the approximation results is verified through a series of numerical examples in Section 6.

2 Preliminaries

2.1 Notation

We employ the multi-indices ${i}=(i_{1},\dots,i_{d})$, and ${\alpha}=(\alpha_{1},\dots,\alpha_{d})$, where each component is non-negative. We denote by $|\cdot|$ the $l_{1}$–norm of the multi-index ${i}$, with $|{i}|=\sum_{j=1}^{d}|i_{k}|$. Further, for multi-indices, the relation ${i}\geq{\alpha}$ means that $i_{k}\geq\alpha_{k}$ for all $k=1,\dots,d$.

Next, we define the following shorthand notation for the summations of indices. For multi-indices $i$ and $\alpha$ satisfying $i\geq\alpha$ we define $\sum_{i\geq\alpha}^{\infty}:=\sum_{i_{1}=\alpha_{1}}^{\infty}\dots\sum_{i_{d}=\alpha_{d}}^{\infty}$ and the summation for multi-indices $i$ satisfying $|i|\geq p$ is defined as $\sum_{|i|=p}^{\infty}.$ Moreover, we also define a summation for a multi-index $i$ satisfying multiple conditions, e.g. multi-index $i$ satisfying the condition $i\geq\alpha$ and the condition $|i|\geq p$ is defined as $\sum_{|i|=p,i\geq\alpha}^{\infty}.$

We introduce a function $\Phi_{d}(m,n)$ which will be used frequently in this work, given by

[TABLE]

where $\Gamma$ is the Gamma function satisfying $\Gamma(n+1)=n!$ for integer $n\geq 0$.

2.2 Weighted Sobolev spaces

For the reference element $\hat{\kappa}:=(-1,1)^{d}$, let $W^{\alpha}({x})=\prod_{k=1}^{d}W_{k}({x}_{k})^{\alpha_{k}},$ where the weight function $W_{k}(x_{k}):=(1-x^{2}_{k})^{1/2},$ for $k=1,\dots,d$, and $\alpha_{k}\geq 0$ are integers.

Next, we define the weighted Sobolev spaces $V^{l}(\hat{\kappa})$ as a closure of $C^{\infty}(\hat{\kappa})$ in the norm with the weights $W^{\alpha}$, defined by

[TABLE]

It is easy to see that $|u|_{V^{l}(\hat{\kappa})}\leq|u|_{H^{l}(\hat{\kappa})}$, $\forall u\in H^{l}(\hat{\kappa})$, with some integer $l\geq 0$. We note that the above definition for weighted Sobolev spaces can be extended to the fractional order weighted Sobolev spaces and weighted Besov spaces by using the real interpolation techniques, cf. MR0482275 .

For $u\in L^{2}(\hat{\kappa})$, we introduce the Legendre polynomial expansion over the reference element $\hat{\kappa}$, given by $u(x)=\sum_{|i|=0}^{\infty}a_{i}\prod_{k=1}^{d}L_{i_{k}}(x_{k}),$ where $x=(x_{1},\dots,x_{d})$, and $L_{i_{k}}({x}_{k})$ denotes the Legendre polynomial with order $i_{k}$ over the variable $x_{k}$. The coefficients $a_{i}$ are defined by

[TABLE]

The derivatives of the function $u$ can be expressed as

[TABLE]

The derivatives of the Legendre polynomials satisfy the orthogonality property

[TABLE]

see (schwab, , Lemma 3.10). By employing (7), we have

[TABLE]

Identity (8) establishes a link between the derivatives of the functions in the weighted $L^{2}$–norms and their Legendre polynomial expansions.

Remark 1

The weighted Sobolev space in the above definition is a special case of the general Jacobi-weighted Sobolev spaces introduced in babuska2002direct_part1 . The key reason to introduce the Jacobi-weighted Sobolev spaces is to deal with the loss of orthogonality suffered by orthogonal polynomials in standard Sobolev spaces; the $L^{2}$-orthogonality is preserved in Jacobi-weighted Sobolev spaces. As we shall see in the forthcoming analysis, orthogonality plays a key role in deriving optimal error bounds in the polynomial order $p$.

3 The $L^{2}$-orthogonal projection operator onto the ${\cal P}_{p}$ basis

In this section, we derive an $hp$-optimal error bound for the $L^{2}$-orthogonal projection over the reference element $\hat{\kappa}:=(-1,1)^{d}$.

3.1 The $L^{2}$-orthogonal projection operator

For the reference element $\hat{\kappa}$, we define ${\cal P}_{p}(\hat{\kappa})$ and ${\cal Q}_{p}(\hat{\kappa})$ be the space of all polynomials with total degree less than or equal to $p$ and with seperate degree less than or equal to $p$, respectively.

In order to distinguish the same projections onto spaces with different polynomial basis, we use subscripts to signify the basis type: we use $\Pi_{{\cal Q}_{p}}:=\Pi_{p}^{(1)}\Pi_{p}^{(2)}\dots\Pi_{p}^{(d)}$ to denote the $L^{2}$-projection onto ${\cal Q}_{p}$, which is constructed by using tensor product arguments together with the one dimensional $L^{2}$-projection with respect to variable $x_{k}$, given by $\Pi_{p}^{(k)}$. On the other hand, the $L^{2}$-projection onto ${\cal P}_{p}$ is denoted by $\Pi_{{\cal P}_{p}}$.

First, we have the following $hp$-optimal approximation result for the $L^{2}$-orthogonal projection $\Pi_{{\cal Q}_{p}}$, c.f. (newpaper, , Lemma 3.4).

Lemma 1

Let $\hat{\kappa}=(-1,1)^{d}$. Suppose that $u\in H^{l}(\hat{\kappa})$, for some interger $l\geq 0$. Let $\Pi_{{\cal Q}_{p}}u$ be the $L^{2}$-projection of $u$ onto $\mathcal{Q}_{p}(\hat{\kappa})$ with $p\geq 0$. Then, for any integer $s$, with $0\leq s\leq\min\{p+1,l\}$, and $W_{k}=W_{k}({x}_{k})$, we have:

[TABLE]

where $\Phi_{1}(p+1,s)$ is defined in (3).

Proof

The result is proved by modifying the proof of Lemma 3.4 in newpaper . Instead of using triangle inequality, we use the orthogonality and stability of the one-dimensional $L^{2}$-orthogonal projection, which leads to the error bound (9).

We remark on the asymptotic behaviour of the Gamma function. Making use of sharp double side inequalities for the Gamma function, see Theorem 1.6. in Batir2008 , for all positive real numbers $x\geq 1$, we have

[TABLE]

and it follows

[TABLE]

with $0\leq s\leq\min\{p+1,l\}$ and $C(s)$ depending on the constant $s$ only. This implies that the error bound (9) is optimal in $p$ with respect to both the Sobolev regularity index $l$ and polynomial order $p$. In fact, by modifying the proof of Theorem 6.2 in MR3556398 , it is can be shown that the constant $C(s)=(\frac{e}{2})^{2s}$.

Next, we introduce a useful lemma which is the key tool in proving the optimal error bounds in $p$. The proof of the lemma is postponed until Section 3.2.

Lemma 2

*Let $\xi=(\xi_{1},\xi_{2},\dots,\xi_{d})$ and $\rho=(\rho_{1},\rho_{2},\dots,\rho_{d})$ be two non-negative integer valued vectors with $\rho\geq\xi$, satisfying $|\rho|=M$, $|\xi|=m$ for $M,m\in\mathbb{N}$. Then, we have the (global) upper bound *

[TABLE]

Furthermore, the maximum value of $F(\xi,\rho)$ under the above constraints on $\xi$ and $\rho$ is attained at $\xi_{k}=m/d$, $\rho_{k}=M/d$, $k=1,\dots,d$.

Theorem 3.1

Let $\hat{\kappa}=(-1,1)^{d}$. Suppose that $u\in H^{l}(\hat{\kappa})$, for some integer $l\geq 0$. Let $\Pi_{{\cal P}_{p}}u$ be the $L^{2}$-projection of $u$ onto $\mathcal{P}_{p}(\hat{\kappa})$ with $p\geq 0$. Then, for any integer $s$, $0\leq s\leq\min\{p+1,l\}$, we have:

[TABLE]

where $\Phi_{d}(p+1,s)$ is defined in (3).

Proof

Using the relation (7) , for any integer $s$, $0\leq s\leq\min\{p+1,l\}$, we have

[TABLE]

where in step one, the index set is enlarged; indeed, some of the terms with multi-index $|i|\geq p+1$ have been used more than once; in step three, we use Lemma 2, taking $\xi_{k}=\alpha_{k}\geq 0$, $\rho_{k}=i_{k}\geq 0$, $M=p+1$, $m=s$, together with the restriction $0\leq s\leq\min\{p+1,l\}$; in step four, we used (8) and in the last step the bound holds from (11).

Remark 2

We point out that the above proof for the $L^{2}$-orthogonal projection $\Pi_{{\cal P}_{p}}$ on $d$-dimensional reference element is a natural extension of the proof for one-dimensional result, see schwab for details.

By comparing the $L^{2}$–norm bound (9) for the projection $\Pi_{{\cal Q}_{p}}$ and (13) for the projection $\Pi_{{\cal P}_{p}}$, it is easy to see that both bounds are $p$-optimal with respect to Sobolev regularity index $l$ and also for polynomial order $p$. Moreover, we can see that the bound in (13) will have a larger constant compared to the bound in (9), and this constant depends on the dimension $d$. This result will play a key role in deriving the exponential convergence for the ${\cal P}_{p}$ basis.

3.2 The Proof of Lemma 2

The proof will be split into three steps.

Step 1: The proof follows a constrained optimization procedure. We set,

[TABLE]

and we calculate the stationary points. We consider the partial derivatives with respect to $\xi_{k}$ and $\rho_{k}$, $k=1,\dots,d$,

[TABLE]

and

[TABLE]

which satisfy the equations

[TABLE]

with $k=1,\dots,d$, by using the fact that $F(\xi,\rho)>0$. The right-hand sides of the two equations in (15) are independent of the index $k$. Moreover, the function $\phi(z)=\Gamma(z)^{\prime}/\Gamma(z)$ is the so-called Digamma function with the property that (see MR0167642 , (6.3.16))

[TABLE]

where $\gamma$ is the Euler-Mascheroni constant. For $z\geq 0$, the function $\phi(z+1)$ is a continuous monotonically increasing function, which shows that (15) have only one solution. This solution is $\tilde{\xi}_{k}=m/d$ and $\tilde{\rho}_{k}=M/d$, $k=1,\dots,d$, and the $F(\xi,\rho)$ will have the extreme value at this stationary point, given by

[TABLE]

Step 2: In order to find the global maximum, we need to prove the following asymptotic relationship:

[TABLE]

This is proven by considering three different cases. We first consider the special case $m=0$. In this case, (17) holds trivially. Next, we consider the case $m=\delta M$, with $0<\delta<1$. By using the property (10) of Gamma functions, we have the following bound:

[TABLE]

By recalling that $0<\delta<1$ and $n=1,\dots,d-1$, we have that $0<\frac{1-\delta}{1+\delta}<1$ and the function $(\frac{d}{n})^{2\delta M}$ is monotonically increasing with respect to $M$. For $M\geq\left((d+n)\log({\frac{e}{\sqrt{2\pi}}})+\frac{d-n}{2}\log({\frac{1+\delta}{1-\delta}})\right)\left(2\delta\log(\frac{d}{n})\right)^{-1}$, the above quotient formula is greater than $1$ and therefore (17) holds.

Finally, we consider the case $m=M$. Using the same techniques used to derive (18) together with the fact that $\Gamma(1)=1$, we have

[TABLE]

By using the fact that exponentially increasing functions grow faster than polynomials, we know that for sufficiently large $M$ the right hand side of (19) is greater than $1$ and therefore (17) holds.

Step 3: Finally, we need to show that the extreme value (16) is the global maximum value of $F({\xi},{\rho})$ under the constraints $|\xi|=m$ and $|\rho|=M$.

First, we can see that the function $F({\xi},{\rho})$ is symmetric and continuous with respect to $\xi$ and $\rho$. The constraints $|\xi|=m$ and $|\rho|=M$ restrict the domain of $\xi$ and $\rho$ to be a $(d-1)$-dimensional simplex, which is convex and compact. So the maximum value of the function $F({\xi},{\rho})$ over the domain will be obtained only at the boundary of the domain or the stationary point of $F(\xi,\rho)$. We have calculated the function value at the stationary point in (16) already, so now we just need to check the function values on the boundary of the domain.

This may be proved by induction. We start with the case $d=2$: the domain of $\xi$ and $\rho$ satisfying the constrains are two straight lines, ${\rho_{1}}+{\rho_{2}}=M$ and ${\xi_{1}}+{\xi_{2}}=m$. Here, the stationary point is the mid-point of each of the two lines $\tilde{\xi}=(m/2,m/2)$, $\tilde{\rho}=(M/2,M/2)$, and the boundary of the domain consist of the points $\xi^{b}=(0,m)$, $\rho^{b}=(0,M)$ or $\xi^{b}=(m,0)$, $\rho^{b}=(M,0)$, due to the constraints $\rho\geq\xi$. Using the symmetry of the function and of the domain, we know that at the two boundary points of the domain, $F(\xi,\rho)$ will attain the same value, with $F(\xi^{b},\rho^{b})=\Phi_{1}(M,m)$. By using the asymptotic relation (17), we find

[TABLE]

The above relation shows that the extreme value (16) is the global maximum value under the constraints for $d=2$.

Next, we consider the case $d=3$, where the domain of each of $\xi$ and $\rho$ will be a triangle. In this case, the stationary point of $F(\xi,\rho)$ is when $\xi$ and $\rho$ are located at the barycenter of their respective triangle. The boundary of each domain consists of $3$ straight lines. We need to calculate the maximum value of $F(\xi,\rho)$ on the boundary of the domain. By using the symmetry of $F(\xi,\rho)$, and that fact that $|\xi|=m$ and $|\rho|=M$, we only need to consider one part of domain boundary where $\xi_{3}$ = 0 and $\rho_{3}=0$. Then, the maximum of $F(\xi,\rho)$ on the domain boundary can be viewed as exactly the same problem with the same constraints as in the case $d=2$. Consequently, the maximum value of $F(\xi,\rho)$ along the boundary of the domain is $F(\xi^{b},\rho^{b})=\Phi_{2}(M,m)$. Again, by using the same techniques as for $d=2$, we deduce that

[TABLE]

The above relation shows that the extreme value (16) is the global maximum value under the constraints for $d=3$. For the general $d$-dimensional case, the proof can be carried out in a similar way. The key observation is that the maximum value of $F(\xi,\rho)$ on the boundary of $d$-dimensional domain will be at the stationary point of $F({\xi},{\rho})$ on the $(d-1)$-dimensional domain. By using the relation

[TABLE]

the proof is complete.

4 The $H^{1}$-projection operator onto the ${\cal S}_{p}$ basis

In this section, we shall consider the $H^{1}$-projection over the reference element $\hat{\kappa}:=(-1,1)^{d}$ with $d=2,3$. Since the three dimensional results depend on the two dimensional results, we start with the two dimensional case.

4.1 The $H^{1}$-projection operator on the reference square

First, we introduce the two-dimensional serendipity finite element space, cf. schwab

[TABLE]

We can see in Figure 1 that the serendipity space ${\cal S}_{p}$ contains two more basis functions than the ${\cal P}_{p}$ basis for $p\geq 2$. Another way to define the serendipity basis is to consider the decomposition of the $C^{0}$ finite element space with ${\cal Q}_{p}$ basis over $\hat{\kappa}$. For polynomial order $p$, the ${\cal S}_{p}$ basis has the same number of nodal basis functions and edge basis functions as the ${\cal Q}_{p}$ basis, but the ${\cal S}_{p}$ basis only has internal moment basis functions (those with zero value along the element boundary $\partial\hat{\kappa}$) whose total degree is less than or equal $p$, cf. schwab ; MR1164869 . We note that serendipity FEMs can be defined in a dimension-independent fashion, see arnold2011serendipity .

Similarly to the case of the $L^{2}$-projection, we use $\pi_{{\cal Q}_{p}}:=\pi_{p}^{(1)}\pi_{p}^{(2)}$ to denote the $H^{1}$-projection onto the ${\cal Q}_{p}$ basis, which can be constructed via a tensor product of one dimensional $H^{1}$-projection with respect to variable $x_{k}$, given by $\pi_{p}^{(k)}$. Similarly, the $H^{1}$-projection onto the ${\cal S}_{p}$ basis is denoted by $\pi_{{\cal S}_{p}}$, which is defined in (4.1).

Now, we construct the two-dimensional $H^{1}$-projection explicitly by using the one-dimensional $H^{1}$-projection and tensor product arguments, see schwab ; hss . For $u\in H^{l}(\hat{\kappa})$, $l\geq 2$, the projection $\pi_{{\cal Q}_{p}}u\in{\cal Q}_{p}(\hat{\kappa})$, $p\geq 1$, is defined by

[TABLE]

the projection $\Pi_{{\cal Q}_{p-1}}$ and $\Pi_{p-1}^{(k)}$ are the two dimensional and one dimensional $L^{2}$-orthogonal projections, respectively, the coefficients $a_{{i_{1}}{i_{2}}}$, $b_{i_{1}}$ and $c_{i_{2}}$ are given by:

[TABLE]

and the polynomial function $\psi_{j}(z)=\int_{-1}^{z}L_{j}(z)\,\mathrm{d}z$ with degree $j+1$, and satisfies $\psi_{j}(\pm 1)=0$ for $j\geq 1$. Moreover, for $j\geq 1$, $\psi_{j}(z)=-\frac{1}{j(j+1)}(1-z^{2})L_{j}^{\prime}(z)$ has the following properties, cf. schwab ,

[TABLE]

Next, we rearrange the relation (4.1) by separating the internal moment basis functions:

[TABLE]

so that the first double summation in (4.1) only contains the internal moment basis functions. From the definition of ${\cal S}_{p}$, $\pi_{{\cal S}_{p}}$ can be constructed by removing the internal moment basis functions with polynomial order greater than $p$ in $\pi_{{\cal Q}_{p}}$. More specifically, $\pi_{{\cal S}_{p}}u\in{\cal S}_{p}(\hat{\kappa})$, $p\geq 4$, is defined by

[TABLE]

For $1\leq p\leq 3$, the first term in (4.1) will vanish, because there are no internal moment basis functions for the serendipity basis in that case. In this work, we focus on the high order polynomial cases, so we only consider the $H^{1}$-projection $\pi_{{\cal S}_{p}}$ for $p\geq 4$.

Next, we recall the following approximation lemma for $\pi_{{\cal Q}_{p}}$ from hss .

Lemma 3

Let $\hat{\kappa}=(-1,1)^{2}$. Suppose that $u\in H^{l+1}(\hat{\kappa})$, for some $l\geq 1$. Let $\pi_{{\cal Q}_{p}}u$ be the $H^{1}$-projection of $u$ onto $\mathcal{Q}_{p}(\hat{\kappa})$ with $p\geq 1$. Then, we have

[TABLE]

and the following error estimates hold:

[TABLE]

and

[TABLE]

for any integer $s$, $1\leq s\leq\min\{p,l\}$.

Then, we derive the $hp$-error bound for the $H^{1}$-projection $\pi_{{\cal S}_{p}}$ for $p\geq 4$.

Theorem 4.1

Let $\hat{\kappa}=(-1,1)^{2}$. Suppose that $u\in H^{l+1}(\hat{\kappa})$, for some $l\geq 1$. Let $\pi_{{\cal S}_{p}}u$ be the $H^{1}$ projection of $u$ onto $\mathcal{S}_{p}(\hat{\kappa})$ with $p\geq 4$. Then, we have

[TABLE]

and for any integer $s$, $1\leq s\leq\min\{p,l\}$, $p$ sufficiently large, the following error estimates hold:

[TABLE]

and

[TABLE]

Proof

The key observation is the fact that the serendipity basis $\mathcal{S}_{p}$ differs from the $\mathcal{Q}_{p}$ basis only at the internal moment basis functions which vanish along the boundary of $\hat{\kappa}$. Indeed, using (4.1) and (4.1), we have

[TABLE]

Using the fact that $\psi_{j}(\pm 1)=0$ for $j\geq 1$, we deduce that $(\pi_{{\cal Q}_{p}}u-\pi_{{\cal S}_{p}}u)|_{\partial{\hat{\kappa}}}=0$. Thus, (29) is proved.

Next, we derive (4.1). The first step is the use of the triangle inequality,

[TABLE]

Thus, we only need to consider the error from the second term in the above bound. By using the orthogonality relation (23) of $\psi_{j}(x)$ for $j\geq 1$ and $1\leq s\leq\min\{p,l\}$, we have

[TABLE]

where in step two, we enlarged the summation index sets by adding the high order internal moment basis functions with coefficients $a_{i_{1}i_{2}}$, $i_{k}\geq 1$ for $k=1,2$ and $|i|\geq p-1$. Thus, we have

[TABLE]

Where we used $\frac{1}{i_{k}(i_{k}+1)}\leq\frac{6}{(i_{k}+\alpha_{k}+1)(i_{k}+\alpha_{k}+2)}$, since $i_{k}\geq\alpha_{k}$ and $i_{k}\geq 1$. Now, we have

[TABLE]

in step one, we enlarge the index set by adding functions with coefficients $a_{i_{1}i_{2}}$ whose index satisfying the relation $|i|\geq p-1$, $\prod_{k=1}^{2}i_{k}=0$, while in step two, we use Lemma 2, with $\xi_{1}=\alpha_{1}+1\geq 1$, $\xi_{2}=\alpha_{2}+1\geq 1$, $\rho_{1}=i_{1}+1\geq 1$, $\rho_{2}=i_{2}+1\geq 1$, $M=p+1$, and $m=s+1$, together with the restriction $1\leq s\leq\min\{p,l\}$; in step three, we use (8) and (4.1) to build up the link between the derivatives of $u$ and coefficients $a_{i_{1}i_{2}}$ and in the last step, we use (11).

Using the same techniques, we can derive the error estimate for the $H^{1}$–seminorm. We have

[TABLE]

In the last step, we enlarge the summation index sets by adding the high order internal moment basis functions with coefficients $a_{i_{1}i_{2}}$, $i_{k}\geq 1$ for $k=1,2$ and $|i|\geq p-1$. Thus, we have

[TABLE]

where in step two, we enlarge the index set again; in step three we use Lemma 2, taking $\xi_{1}=\alpha_{1}\geq 0$, $\xi_{2}=\alpha_{2}+1\geq 1$, $\rho_{1}=i_{1}\geq 0$, $\rho_{2}=i_{2}+1\geq 1$, $M=p$, and $m=s$, together with the restriction $1\leq s\leq\min\{p,l\}$.

Therefore, we have the bound

[TABLE]

Finally, using (4.1), (39) and Lemma 3, the bounds (4.1) and (4.1) follow.

4.2 The $H^{1}$-projection operator on the reference cube

In this section, we shall consider the $H^{1}$-projection operator over the reference cube $\hat{\kappa}:=(-1,1)^{3}$. First, we introduce the $3$D serendipity finite element space.

A simple way to define the serendipity basis is to consider a decomposition of the $C^{0}$ finite element space with ${\cal Q}_{p}$ basis over $\hat{\kappa}$. For polynomial order $p$, the ${\cal S}_{p}$ basis has the same number of nodal basis functions and edge basis functions as the ${\cal Q}_{p}$ basis, but the ${\cal S}_{p}$ basis only has face basis functions (those with zero value on twelve edges and eight vertices) and internal moment basis functions (those with zero value along the element boundary $\partial\hat{\kappa}$) whose total degree is less than or equal $p$. The number of basis functions of ${\cal S}_{p}$ basis is calculated in the following way

[TABLE]

here, we note that for $p=1$, the serendipity basis only contains $8$ nodal basis functions and ${\cal S}_{1}(\hat{\kappa}):={\cal Q}_{1}(\hat{\kappa})$. For $p\geq 2$, the serendipity basis contains $(p-1)$ edge basis functions for each of the $12$ edges. For $p\geq 4$, the serendipity basis contains $(p-2)(p-3)/{2}$ face basis functions for each of the $6$ faces. For $p\geq 6$, the serendipity basis contains ${(p-3)(p-4)(p-5)}/{6}$ internal moment basis functions.

Similarly to the $2$D case, we use $\pi_{{\cal Q}_{p}}:=\pi_{p}^{(1)}\pi_{p}^{(2)}\pi_{p}^{(3)}$ to denote the $H^{1}$-projection onto the ${\cal Q}_{p}$ basis. The $H^{1}$-projection onto the ${\cal S}_{p}$ basis is denoted by $\pi_{{\cal S}_{p}}$. Additionally, we introduce some new notation for the forthcoming analysis. The projection $\pi_{{\cal S}_{p}}^{(1,2)}$ shall denote the $H^{1}$-projection onto the serendipity spaces ${\cal S}_{p}$ with variables $(x_{1},x_{2})$ only, and the projections $\pi_{{\cal S}_{p}}^{(1,3)}$ and $\pi_{{\cal S}_{p}}^{(2,3)}$ are defined in an analogous manner.

First, we explicitly construct the three-dimensional projection $\pi_{{\cal Q}_{p}}=\pi_{p}^{(1)}\pi_{p}^{(2)}\pi_{p}^{(3)}$. For $u\in H^{l}(\hat{\kappa})$, $l\geq 3$, the projection $\pi_{{\cal Q}_{p}}u\in{\cal Q}_{p}(\hat{\kappa})$, $p\geq 1$, is defined by

[TABLE]

Then, the following Legendre polynomial expansion holds:

[TABLE]

with coefficients $a_{{i_{1}}{i_{2}}{i_{3}}}$, $b_{i_{1}i_{2}}$, $c_{i_{1}i_{3}}$, $d_{i_{2}i_{3}}$, give by

[TABLE]

together with $e_{i_{1}}$, $f_{i_{2}}$ and $g_{i_{3}}$

[TABLE]

Now, we separate the face basis functions and internal moment basis functions from (4.2).

[TABLE]

Here, the first triple summation terms contains all the internal moment basis functions only. Three double summation terms contain all the face basis functions. The edge basis functions and nodal basis functions will not be written explicitly because they play no role in the analysis.

From the definition of ${\cal S}_{p}$, $\pi_{{\cal S}_{p}}u$ can be constructed by removing the face basis functions and internal moment basis functions with polynomial order greater than $p$ in $\pi_{{\cal Q}_{p}}u$. More specifically, $\pi_{{\cal S}_{p}}u\in{\cal S}_{p}(\hat{\kappa})$, $p\geq 6$, is defined by

[TABLE]

For $1\leq p\leq 3$, both face basis functions and internal moment basis functions in (4.2) will vanish. For $4\leq p\leq 5$, internal moment basis functions in (4.2) will vanish. Similar to the $2$D case, we only consider the $H^{1}$-projection $\pi_{{\cal S}_{p}}$ for $p\geq 6$.

Next, by using the stability and approximation results for one dimensional $H^{1}$-projection in hss , we can derive the following approximation results for $\pi_{{\cal Q}_{p}}$.

Lemma 4

Let $\hat{\kappa}=(-1,1)^{3}$. Suppose that $u\in H^{l+1}(\hat{\kappa})$, for some $l\geq 2$. Let $\pi_{{\cal Q}_{p}}u$ be the $H^{1}$-projection of $u$ onto $\mathcal{Q}_{p}(\hat{\kappa})$ with $p\geq 1$. Then, we have

[TABLE]

*and the following error estimates hold: *

[TABLE]

*and *

[TABLE]

for any integer $s$, $2\leq s\leq\min\{p,l\}$.

Then, we derive the $hp$-error bound for the $H^{1}$-projection $\pi_{{\cal S}_{p}}$ for $p\geq 6$.

Theorem 4.2

Let $\hat{\kappa}=(-1,1)^{3}$. Suppose that $u\in H^{l+1}(\hat{\kappa})$, for some $l\geq 2$. Let $\pi_{{\cal S}_{p}}u$ be the $H^{1}$ projection of $u$ onto $\mathcal{S}_{p}(\hat{\kappa})$ with $p\geq 6$. Then, we have

[TABLE]

and for any integer $s$, $2\leq s\leq\min\{p,l\}$, $p$ sufficiently large, the following error estimates hold:

[TABLE]

and

[TABLE]

Here, $C_{1}$ and $C_{2}$ are positive constants independent of $p$, $l$ and $s$.

Proof

The proof of (49) is similar to the two dimensional case, by using the fact that the serendipity basis $\mathcal{S}_{p}$ differs from $\mathcal{Q}_{p}$ only in the face basis functions and internal basis functions which have zero values at each vertex.

Next, we begin to prove relation (4.2). Using (4.2) and (4.2), we have

[TABLE]

We note that the term $T_{1}$ only contains the internal moment basis functions, and the three other terms only contain the face basis functions. By using the orthogonality relation (23) of $\psi_{j}(x)$ for $j\geq 1$ and $2\leq s\leq\min\{p,l\}$, we have

[TABLE]

where in step two, we enlarged the summation index sets by adding the high order internal moment basis functions with coefficients $a_{i_{1}i_{2}i_{3}}$, $i_{k}\geq 1$ for $k=1,2,3$ and $|i|\geq p-2$; in the last step, we used the relation $\frac{1}{i_{k}(i_{k}+1)}\leq\frac{6}{(i_{k}+\alpha_{k}+1)(i_{k}+\alpha_{k}+2)}$ for internal moment basis functions with $i_{k}\geq 1$ and $i_{k}\geq\alpha_{k}$, $k=1,2,3$. Thus, we have

[TABLE]

where in step one, we enlarged the summation index set by adding functions with coefficients $a_{i_{1}i_{2}i_{3}}$ whose index satisfying the relation $|i|\geq p-2$, $\prod_{k=1}^{3}i_{k}=0$; in step two, we used Lemma 2, with $\xi_{k}=\alpha_{k}+1\geq 1$, $\rho_{k}=i_{k}+1\geq 1$, $k=1,2,3$, together with $M=p+1$, $m=s+1$, and the restriction $2\leq s\leq\min\{p,l\}$ and in the last step, we used the relation (11).

Next, we will derive the error bound for the term $T_{2}$. We first rewrite $T_{2}$ in the following way by adding and subtracting the same terms

[TABLE]

The key observation is that $T_{2,a}=\pi_{p}^{(1)}\pi_{p}^{(2)}u-\pi_{{\cal S}_{p}}^{(1,2)}u$, see Appendix. By using the $2$D approximation results (4.1), we have the following bound

[TABLE]

We note that $T_{2,b}$ only contains the internal moment basis functions, so using (LABEL:weight) and the orthogonality relation (23) of $\psi_{j}(x)$ for $j\geq 1$, we have

[TABLE]

where in step two, we used the fact that the multi-index set for $T_{2,b}$ is a subset of multi-index $i$ with $|i|\geq p-2$ and $i_{k}\geq 1$, $k=1,2,3$ and then, we enlarged the index set of the summations by adding the high order internal moment basis functions with coefficients $a_{i_{1}i_{2}i_{3}}$, $i_{k}\geq 1$ for $k=1,2,3$ and $|i|\geq p-2$. By using and (4.2), the following bound holds

[TABLE]

Combining (4.2) and (58) together with the asymptotic relation (17), the following bound for $T_{2}$ holds

[TABLE]

It is easy to verify that $T_{3}$ and $T_{4}$ satisfy the same error bound as $T_{2}$, and thus the $L^{2}$–norm error bound (4.2) is proved.

The error estimate for the $H^{1}$–seminorm (4.2) can be derived by using the same techniques used in deriving (4.2). As such, we have

[TABLE]

Now, we have

[TABLE]

where in step one, we enlarged the summation index set by adding functions with coefficients $a_{i_{1}i_{2}i_{3}}$ whose index satisfying the relation $|i|\geq p-2$, $\prod_{k=1}^{3}i_{k}=0$; in step two we use Lemma 2, taking $\xi_{1}=\alpha_{1}\geq 0$, $\xi_{2}=\alpha_{2}+1\geq 1$, $\xi_{3}=\alpha_{3}+1\geq 1$, $\rho_{1}=i_{1}\geq 0$, $\rho_{2}=i_{2}+1\geq 1$, $\rho_{3}=i_{3}+1\geq 1$, $M=p$, and $m=s$, together with the restriction $2\leq s\leq\min\{p,l\}$.

Next, we consider the error bound for term $T_{2}$. By using (4.1), the following bound holds for $T_{2,a}$

[TABLE]

With the help of (4.2) and (4.2), the following bound holds

[TABLE]

Combing (4.2), (4.2) together with the asymptotic relation (17), then the following error bound for term $T_{2}$.

[TABLE]

The above error bound for $T_{2}$ term is also the error bound for terms $T_{3}$ and $T_{4}$. By using (4.2) and (4.2), it is to see that the following relation holds

[TABLE]

where $C^{*}$ is a positive constant independent of $p$, $s$ and $l$. We note that above $L^{2}$–norm error bound also holds for the rest two partial derivatives. So the $H^{1}$–seminorm bound (4.2) is proved.

Remark 3

We again make a comparison between the bounds in the $L^{2}$–norm and $H^{1}$–seminorm, given in Lemma 3 for $d=2$ and Lemma 4 for $d=3$ respectively for $\pi_{{\cal Q}_{p}}$, and Theorem 4.1 for $d=2$ and Theorem 4.2 for $d=3$ respectively for $\pi_{{\cal S}_{p}}$. Similarly to the comparisons for the $L^{2}$-projection onto ${\cal P}_{p}$ and ${\cal Q}_{p}$, both bounds are $p$-optimal in both Sobolev regularity and polynomial order. We can also see that the bounds for $\pi_{{\cal S}_{p}}$ have a larger constant than those for $\pi_{{\cal Q}_{p}}$, and this constant depends on dimension $d$. Moreover, we point out that the optimal approximation results for the $H^{1}$-projection with ${\cal S}_{p}$ basis in Theorem 4.1 and Theorem 4.2 directly imply the $hp$-optimal error bound for the $L^{2}$–norm on the trace of $\hat{\kappa}$ for $\pi_{{\cal S}_{p}}$.

Remark 4

We note that in the Theorem 4.1 and 4.2, the minimum Sobolev regularity requirements for defining $H^{1}$-projection is $u\in H^{d}(\hat{\kappa})$ for the reference element. In fact, this regularity requirement can be relaxed by using the the tensor product Sobolev spaces, cf. schwab ; thesis . In this work, we do not consider the minimum regularity assumptions because we only consider the standard Sobolev spaces.

4.3 The $H^{1}$-projection operator onto the ${\cal P}_{p}$ basis

Finally, we present the error bound for $\pi_{{\cal P}_{p}}$ which we shall define now. The key observation is that the ${\cal P}_{p}$ basis with polynomial order $p$ contains the ${\cal S}_{p+1-d}$ basis for $p\geq d$, see arnold2011serendipity . Then, we can simply define $\pi_{{\cal P}_{p}}=\pi_{{\cal S}_{p+1-d}}$ for $d=2,3$.

Corollary 1

Let $\hat{\kappa}=(-1,1)^{d}$, $d=2,3$. Suppose that $u\in H^{l+1}(\hat{\kappa})$, for some $l\geq d-1$. Let $\pi_{{\cal P}_{p}}u:=\pi_{{\cal S}_{p+1-d}}u$ be the $H^{1}$ projection of $u$ onto $\mathcal{P}_{p}(\hat{\kappa})$ with $p\geq 3d-1$. Then, we have:

[TABLE]

and the following error estimates hold:

[TABLE]

and

[TABLE]

for any integer $s$, $d-1\leq s\leq\min\{p+1-d,l\}$, $p$ sufficiently large.

Remark 5

We emphasize that * the above error bound for the $\pi_{{\cal P}_{p}}$ projection is $p$-suboptimal by one order for $d=2$ and two orders for $d=3$ for sufficiently smooth functions, but it is $p$-optimal for functions with finite Sobolev regularity in the case $l\leq p+1-d$. * However, sub-optimality by one or two orders in $p$ is better than using the $\pi^{\cal Q}_{\lfloor p/d\rfloor}$ projection, as suggested by schwab (see Corollary 4.52 on p190), which is sub-optimal in $p$ by at least $p/2$ orders for sufficiently smooth functions for $d=2$. Moreover, the one or two order sub-optimality in $p$ for analytic functions does not influence the exponent of the exponential rate of convergence, as we shall see below.

5 Exponential convergence for analytic solutions

We shall be concerned with the proof of exponential convergence for serendipity FEMs and DGFEMs with the ${\cal P}_{p}$ basis over tensor product elements. For simplicity, we only consider the case when the given problem is piecewise analytic over the whole computational domain. Exponential convergence is then achieved by fixing the computational mesh, and increasing the polynomial order $p$. Only parallelepiped meshes are considered, which are the affine family obtained from the reference element $\hat{\kappa}=(-1,1)^{d}$. The analysis of FEMs and DGFEMs with a general $hp$-refinement strategy is beyond the scope of this analysis (see MR3351174 ; schotzau2013hp1 ; schotzau2013hp2 ; schotzau2016hp for the analysis for both methods employing the ${\cal Q}_{p}$ basis).

The proof of exponential convergence for FEMs and DGFEMs depends on proving exponential convergence of $L^{2}$- and $H^{1}$-projections for piecewise analytic functions under $p$-refinement. The $H^{1}$-projection $\pi_{{\cal S}_{p}}$ onto ${\cal S}_{p}$ can be directly applied to $p$-FEMs for second order elliptic problems with the same optimal rate as the $H^{1}$ projection $\pi_{{\cal Q}_{p}}$, see schwab for details. For deriving error bounds of DGFEMs using the $L^{2}$- and $H^{1}$- projections onto ${\cal Q}_{p}$, we refer to hss ; newpaper ; thesis . Following similar techniques, we can prove the corresponding $hp$-bounds for DGFEMs employing the ${\cal P}_{p}$ basis, albeit with sub-optimal rate in $p$. The sub-optimality in $p$ is due to the fact that the $p$-optimal bound for $L^{2}$-projection onto ${\cal P}_{p}$ basis over the trace of the tensor product elements is still open. Additionally, the $H^{1}$-projection onto the ${\cal P}_{p}$ basis is suboptimal in $p$ by $d-1$ orders for sufficiently smooth functions. However, we point out that the suboptimality in $p$ by $d-1$ order, with $d=2,3$, does not influence the exponent of the exponential rate of convergence.

Next, we focus on deriving the exponential convergence for the $L^{2}$-projections in the $L^{2}$–norm and $H^{1}$-projections in the $L^{2}$–norm and $H^{1}$–seminorm on analytic problems under $p$-refinement on shape-regular $d$-parallelepiped meshes. The extension to anisotropic meshes will be consider in the future.

Let $\kappa$ be a parallelepiped element. For a function $u$ having an analytic extension into an open neighbourhood of $\bar{\kappa}$, we have:

[TABLE]

where $|\kappa|$ denotes the measure of element $\kappa$, cf. (davis1975interpolation, , Theorem 1.9.3).

Lemma 5

Let $u:\kappa\rightarrow\mathbb{R}$ have an analytic extension to an open neighbourhood of $\bar{\kappa}$. Also let $p_{\kappa}\geq 0$ and $0\leq s_{\kappa}\leq p_{\kappa}+1$ be two positive numbers such that $s_{\kappa}=\epsilon(p_{\kappa}+1)$, $0\leq\epsilon\leq 1$ and $d=2,3$. Then the following bounds hold:

[TABLE]

and

[TABLE]

Here, $C$ and $C(u)$ are positive constants depending elemental shape regularity, and $C(u)$ also depends on $u$. $F_{1}(R_{\kappa},\epsilon)=\frac{(1-\epsilon)^{1-\epsilon}}{(1+\epsilon)^{1+\epsilon}}(\epsilon R_{k})^{2\epsilon}$, $\epsilon_{\min}={1}/{\sqrt{1+R_{\kappa}^{2}}}$, $b_{1,\kappa}:=\frac{1}{2}|\log F_{1}(R_{\kappa},\epsilon_{\min})|+\epsilon_{\min}|\log{h_{\kappa}}|$ and $b_{2,\kappa}:=b_{1,\kappa}-\epsilon_{\min}\log d$.

Proof

Using standard scaling arguments for $\kappa$ together with Lemma 1 and Theorem 3.1, we have the approximation results for the $L^{2}$-projection over $\kappa$. For brevity, we set $q_{\kappa}=p_{\kappa}+1$. By employing the relation (11) and the fact $|u|_{V^{l}({\kappa})}\leq|u|_{H^{l}({\kappa})}$, we have the bounds:

[TABLE]

where

[TABLE]

Recalling (69), we have $R_{\kappa}>0$,

[TABLE]

Thus, we have

[TABLE]

Therefore, we have the exponential convergence for the $L^{2}$-projection $\Pi_{{\cal Q}_{p_{\kappa}}}$, via

[TABLE]

with $b_{1,\kappa}:=\frac{1}{2}|\log F_{1}(R_{\kappa},\epsilon_{\min})|+\epsilon_{\min}|\log{h_{\kappa}}|$. Similarly, for the $L^{2}$-projection $\Pi_{{\cal P}_{p_{\kappa}}}$, Stirling’s formula implies

[TABLE]

where,

[TABLE]

with the minimum,

[TABLE]

In order to compare with the slope of projection $\Pi_{{\cal Q}_{p_{\kappa}}}$, here we will use the same $\epsilon_{\min}$. We have

[TABLE]

Thus, we have

[TABLE]

with slope $b_{2,\kappa}:=\frac{1}{2}|\log F_{1}(R_{\kappa},\epsilon_{\min})|+\epsilon_{\min}(|\log{h_{\kappa}}|-\log d)$.

Next, we begin to derive the exponential convergence for $H^{1}$-projections.

Lemma 6

Let $u:\kappa\rightarrow\mathbb{R}$ have an analytic extension to an open neighbourhood of $\bar{\kappa}$. Also let $p_{\kappa}\geq 2d$ and $(d-1)\leq s_{\kappa}\leq p_{\kappa}$ be two positive numbers such that $s_{\kappa}=\epsilon p_{\kappa}$, $0<\epsilon\leq 1$ and $d=2,3$. Then the following bounds hold:

[TABLE]

[TABLE]

and

[TABLE]

[TABLE]

Here, $C$ and $C(u)$ are positive constants depending elemental shape regularity, and $C(u)$ also depends on $u$. $F_{1}(R_{\kappa},\epsilon)=\frac{(1-\epsilon)^{1-\epsilon}}{(1+\epsilon)^{1+\epsilon}}(\epsilon R_{k})^{2\epsilon}$, $\epsilon_{\min}={1}/{\sqrt{1+R_{\kappa}^{2}}}$, $b_{1,\kappa}:=\frac{1}{2}|\log F_{1}(R_{\kappa},\epsilon_{\min})|+\epsilon_{\min}|\log{h_{\kappa}}|$ and $b_{2,\kappa}:=b_{\kappa}^{1}-\epsilon_{\min}\log d$.

Proof

The proof follows by the same techniques used in Lemma 5.

In the above Lemma 5 and Lemma 6, we can see that the $L^{2}$–norm error for both $L^{2}$-projections $\Pi_{{\cal Q}_{p_{\kappa}}}$ and $\Pi_{{\cal P}_{p_{\kappa}}}$, and the $L^{2}$–norm and $H^{1}$–seminorm error for the $H^{1}$-projections $\pi_{{\cal S}_{p_{\kappa}}}$ and $\pi_{{\cal Q}_{p_{\kappa}}}$ decay exponentially for analytic functions under $p$-refinement. If we measure the error against $p$, the exponent $b_{1,\kappa}$ for the ${\cal Q}_{p}$ basis is slightly greater than the exponent $b_{2,\kappa}$ for the ${\cal P}_{p}$ basis and ${\cal S}_{p}$ basis by a small factor of $(\log d)/{\sqrt{1+R_{\kappa}^{2}}}$. By using Lemma 5 and Lemma 6, we can also derive the following theorem.

Theorem 5.1

Let $u$ be an analytic function as defined in (69), and exponent $b_{1,\kappa}$ and $b_{2,\kappa}$ defined in Lemma 5 Then, there exists $C>0$ such that following bounds hold:

[TABLE]

[TABLE]

and

[TABLE]

[TABLE]

and

[TABLE]

[TABLE]

Proof

By recalling the relationship between degrees of freedom and polynomial order $p$ for both the ${\cal Q}_{p}$ and ${\cal P}_{p}$ bases, we have

[TABLE]

and

[TABLE]

Then, (74) and (75) follow from Lemma 5.

By using relations (20) and (40), we have the asymptotic relation

[TABLE]

The relations (76), (77), (78) and (79) follow from the Lemma 6.

For $d=2,3$, if the following condition

[TABLE]

holds, then we have $b_{2,\kappa}\approx b_{1,\kappa}$. It is easy to see that for small $R_{\kappa}$ or small mesh size $h$, the condition (83) will be satisfied. Moreover, we point out that an analytic function having sufficiently small $R_{\kappa}$ is equivalent to the function having an analytic continuation into a sufficiently large open neighbourhood of $\bar{\kappa}$, see davis1975interpolation for details.

Now, if we consider the error in terms of $\sqrt[d]{Dof}$ for the above bounds, the exponent for the exponential convergence rate of the ${\cal P}_{p}$ basis and the ${\cal S}_{p}$ basis are larger than the exponent for the ${\cal Q}_{p}$ basis by a fixed factor of $\sqrt[d]{d!}$.

We have observed a steeper slope in error against $\sqrt[d]{Dof}$ for FEMs with ${\cal S}_{p}$ basis and DGFEMs with ${\cal P}_{p}$ basis. For $d=2$, this suggests a typical ratio between convergence slopes of DGFEMs with ${\cal P}_{p}$ and ${\cal Q}_{p}$ basis, FEMs with ${\cal S}_{p}$ and ${\cal Q}_{p}$ basis, to be $\sqrt{2!}\approx 1.414$. For $d=3$, this ratio is $\sqrt[3]{3!}\approx 1.817$. The numerical examples in Section 6 show that the ratio is slightly worse than the ideal ratio, but it is not far from the ideal ratio.

6 Numerical examples

We present some numerical examples to confirm the theoretical analysis in the previous sections. All the numerical examples are computed by Matlab on the High Performance Computing facility ALICE of the University of Leicester. For simplicity of presentation, we use DGFEM(P) and DGFEM(Q) to denote the DGFEMs with local polynomial basis consisting of either ${\cal P}_{p}$ or ${\cal Q}_{p}$ polynomials and use FEM(S) and FEM(Q) to denote the FEMs with local polynomial basis consisting of either ${\cal S}_{p}$ or ${\cal Q}_{p}$ polynomials.

The comparisons are mainly made between the slope of FEM(S) and FEM(Q) over square meshes for $d=2$ and hexahedral meshes for $d=3$ under $p$-refinement. The slopes of the convergence lines are calculated by taking the average of the last two slopes of the line segments of each convergence line. We will also present an example comparing DGFEM(P) and DGFEM(Q). For more numerical examples for DGFEMs, see D17 ; Zhaonan .

6.1 Example 1

In the first example, we investigate the computational efficiency of DGFEM(P) and DGFEM(Q) schemes. To this end, we consider a partial differential equation with nonnegative characteristic form of mixed type. Let $\Omega=(-1,1)^{2}$, and consider the PDE problem:

[TABLE]

with exact solution:

[TABLE]

This problem is hyperbolic in the region $y\leq 0$ and parabolic for $y>0$. In order to ensure continuity of the normal flux across $y=0$, where the partial differential equation changes type, the exact solution has a discontinuity across the line $y=0$, cf. cangiani2015hp ; thesis .

By following cangiani2015hp , we use the symmetric interior penalty DGFEMs employing a special class of quadrilateral meshes for which the discontinuity in the exact solution lies on element interfaces. In this setting, we modify the discontinuity-penalization parameter $\sigma$, so that $\sigma$ vanishes on edges which form part of the interface $y=0$; this ensures that the (physical) discontinuity present in the exact solution is not penalized within by the numerical scheme.

In this case, the exact solution is piecewise analytic on the two parts of the domain. In Figure 2, we observe that the DG–norm $|\kern-0.7113pt\|{u-u_{h}}|\kern-0.7113pt\|_{\rm DG}$ decays exponentially for both DGFEM(P) and DGFEM(Q) under $p$-refinement on $64$, $4096$ and $16384$ uniform square elements. The definition of DG–norm $|\kern-0.7113pt\|{\cdot}|\kern-0.7113pt\|_{\rm DG}$ can be found in cangiani2015hp . Moreover, the slope of the convergence line for the DGFEM(P) is greater than the line of DGFEM(Q) in error against $\sqrt{Dof}$. The ratio between the two slopes is about $1.39$ on coarse meshes and fine meshes. The numerical observation confirms the theoretical results in Theorem 5.1.

6.2 Example 2

In the second example, we investigate the computational efficiency of FEM(S) and FEM(Q) on standard tensor-product elements (quadrilaterals in 2D and hexahedra in 3D).

Firstly, we consider the following two–dimensional Poisson problem: let $\Omega=(0,1)^{2}$ and select $f=2\pi^{2}\sin(\pi x)\sin(\pi y)$, so that the exact solution is given by $u=\sin(\pi x)\sin(\pi y)$.

In this case, the exact solution is piecewise analytic on the domain. In Figure 3, we observe that the $H^{1}$–seminorm $|u-u_{h}|_{H^{1}(\Omega)}$ decays exponentially for both FEM(S) and FEM(Q) under $p$-refinement on $64$, $4096$ and $16384$ uniform square elements. Again, we observe that the slope of the convergence line for the FEM(S) is greater than the line of FEM(Q) in error against $\sqrt{Dof}$. The ratio between the two slopes is about $1.39$ on coarse meshes and fine meshes.

We now consider the three–dimensional variant of the above problem. Let $\Omega=(0,1)^{3}$ and select $f=3\pi^{2}\sin(\pi x)\sin(\pi y)\sin(\pi z)$, so that the exact solution is given by $u=\sin(\pi x)\sin(\pi y)\sin(\pi z)$.

In Figure 4, we observe that the $H^{1}$–seminorm $|u-u_{h}|_{H^{1}(\Omega)}$ decays exponentially for both FEM(S) and FEM(Q) under $p$-refinement on $64$, $4096$ and $32768$ uniform hexahedral elements. Moreover, we observe that the slope of the convergence line for the FEM(S) is greater than the line of FEM(Q) in error against $\sqrt[3]{Dof}$. The ratio between the two slopes is about $1.62$ on coarse meshes and $1.73$ on fine meshes. The numerical observation confirms the theoretical results in Theorem 5.1.

6.3 Example 3

In the third example, we investigate the convergence behaviour of the FEM(S) and FEM(Q) approaches for the Poisson problem on a non-smooth domain with fixed computational meshes under $p$-refinement. To this end, we let $\Omega$ be the L-shaped domain $(-1,1)^{2}\setminus[0,1)\times(-1,0]$. Uniform square meshes consisting of $12$ elements are used. Then, writing $(r,\varphi)$ to denote the system of polar coordinates, we impose an appropriate inhomogeneous boundary condition for $u$ so that

[TABLE]

cf. Wihler2002 . We note that $u$ is analytic in $\overline{\Omega}\setminus\{{\bf 0}\}$, but $\nabla u$ is singular at the origin; indeed, here $u\not\in H^{2}(\Omega)$. This example reflects the typical (singular) behaviour that solutions of elliptic boundary value problems exhibit in the vicinity of reentrant corners in the computational domain.

In fact, $u\in H^{\frac{5}{3}-\epsilon}(\Omega)$, for any $\epsilon>0$. We investigate the convergence rate of the FEM(S) and FEM(Q) under $p$-refinement for this problem. In Table 1, we list the $H^{1}$–seminorm error and also the convergence rate of FEM(S) and FEM(Q) with polynomial order $p=1,\dots,60$. We point out that due to the singularity at the origin, geometrically graded quadrature points towards the origin are used in order to get the desired accuracy (see chernov2011exponential ).

As we can see, the convergence rate in $p$ for both FEM(S) and FEM(Q) are approximately $\mathcal{O}(p^{4/3})$. The convergence rate in $p$ is double the theoretical rate with respect to the Sobolev regularity of $u$. This is the doubling order convergence in the $p$-version finite element, see Babuska-Suri:SINUM:1987 ; schwab for details. The reason of this doubling order convergence in $p$ is related to the fact that standard Sobolev space can not optimally characterize the singularity of $r^{\gamma}\log^{\nu}r$ type, $\gamma\in\mathbb{R}^{+}$, $\nu\in\mathbb{N}$; indeed from babuska2002direct_part1 ; babuska2002direct_part2 ; MR1754720 , we know that the modified Jacobi-weighted Besov spaces provide a sharper function space setting to characterize such singular functions. By using the results in babuska2002direct_part2 , FEM(Q) has the following sharp error bound under $p$-refinement for this problem

[TABLE]

where constants $C$ and $c$ are independent of $p$. Moreover, we observe that the FEM(S) error is greater than FEM(Q) error by a factor about $2.29$ for fixed $p$. By noting that $2^{4/3}\approx 2.52$, we suppose that the error bound for FEM(S) satisfying the following relation

[TABLE]

The proof the above intuitive optimal $p$-version error bound for FEM(S) is beyond the scope of this work, which depends on the approximation theory for orthogonal projections onto the ${\cal P}_{p}$ basis and the ${\cal S}_{p}$ basis in the modified Jacobi-weighted Besov spaces.

Now, let’s consider the relation of the error against $\sqrt{Dof}$. By using the asymptotic relation (80), we have the following relation for FEM(Q)

[TABLE]

Then, we derive the following relation for FEM(S) by using the relation (82),

[TABLE]

Hence, both FEM(S) and FEM(Q) have the same algebraic convergence rate in degrees of freedom under $p$-refinement. However, the error bound of FEM(S) seems to be larger than the error bound of FEM(Q) by a constant. This result is observed in the Figure 5, where we can see that the convergence line of FEM(S) and FEM(Q) have the same slope in error against $p$ (left) and error against ${Dof}^{1/2}$ (right). But the error of FEM(S) is always larger than the error of FEM(Q) in both graphs.

7 Conclusions and remarks

In this work, we derived new $hp$-optimal approximation results for the $L^{2}$-orthogonal projection onto the total degree ${\cal P}_{p}$ basis, and the $H^{1}$-projection onto the serendipity basis ${\cal S}_{p}$ over tensor product elements. With these results, we proved that the exponent of the exponential rate of convergence with respect to the number of degrees of freedom for the DGFEM employing ${\cal P}_{p}$ basis and the serendipity FEM are greater than the exponent of the exponential rate of convergence for their counterparts employing ${\cal Q}_{p}$ basis for analytic functions under $p$-refinement. Moreover, the exponent for the ${\cal P}_{p}$ and ${\cal S}_{p}$ bases are larger than that of ${\cal Q}_{p}$ basis by a constant only depending on dimension. The sharpness of the theoretical results has been verified by numerical examples.

Finally, we remark on some applications and potential extensions of results in this work. First, we note that the $hp$-optimal error bounds for the $L^{2}$-orthogonal projections onto the ${\cal P}_{p}$ basis and $H^{1}$-projections onto the ${\cal S}_{p}$ basis can be used to improve the $hp$-error bounds for mixed-FEMs employing the ${\bf BDFM}$-elements and virtual element methods, see MR1974174 and MR3509090 . Next, as we have already observed in the numerical example $7.2$ in cangiani2013hp , the $hp$-adaptive DGFEM employing the ${\cal P}_{p}$ basis gives a greater exponent of the exponential rate of convergence in terms of number of degrees of freedom than the exponent of DGFEM employing the ${\cal Q}_{p}$ basis for solving the L-shaped domain problem. We are very interested in applying the results in this work for $hp$-version FEMs employing the serendipity basis and DGFEMs employing the ${\cal P}_{p}$ basis in a general $hp$-refinement setting.

Acknowledgements.

The author wishes to express his gratitude to Emmanuil Georgoulis (University of Leicester & National Technical University of Athens) and Andrea Cangiani (University of Leicester) for their helpful comments. Z. D. was supported by the Leverhulme Trust (grant no. RPG-2015-306).

Appendix. Technical results

We present the details for proving that $T_{2,a}=\pi_{p}^{(1)}\pi_{p}^{(2)}u-\pi_{{\cal S}_{p}}^{(1,2)}u$ in relation (4.2). The idea is that both of the $H^{1}$-projections $\pi_{p}^{(1)}\pi_{p}^{(2)}u$ and $\pi_{{\cal S}_{p}}^{(1,2)}u$ are independent of the variable $x_{3}$ for $u\in H^{3}(\hat{\kappa})$, $\hat{\kappa}=(-1,1)^{3}$. By recalling the definition of the $2$D $H^{1}$-projection $\pi_{p}^{(1)}\pi_{p}^{(2)}u$ in (4.1), then we have the following relation

[TABLE]

with $\tilde{a}_{{i_{1}}{i_{2}}}(x_{3})$, $\tilde{b}_{i_{1}}(x_{3})$ and $\tilde{c}_{i_{2}}(x_{3})$ given by:

[TABLE]

The above coefficients in (Appendix. Technical results) all depend the variable $x_{3}$. We note that for function $u\in H^{3}(\hat{\kappa})$, the above expansion can be rewritten in following way,

[TABLE]

Next, we expand above partial derivative of coefficients $\tilde{a}_{{i_{1}}{i_{2}}}(x_{3})$, $\tilde{b}_{i_{1}}(x_{3})$ and $\tilde{c}_{i_{2}}(x_{3})$ into Legendre polynomials $L_{i_{3}}(x_{3})$, and compare coefficients in (Appendix. Technical results), (4.2) and (4.2), the following relation holds

[TABLE]

Similarly, by using the same techniques, in conjunction with the definition of the $2$D $H^{1}$-projection $\pi_{{\cal S}_{p}}^{(1,2)}u$ in (4.1). The following relation holds for $\pi_{{\cal S}_{p}}^{(1,2)}u$

[TABLE]

By combining above two relation (Appendix. Technical results) and (Appendix. Technical results), we derive the desired relation,

[TABLE]

The proof is complete.