Well-conditioned frames for high order finite element methods

TL;DR

This paper introduces a novel frame-based representation for high order $C^0$ finite element spaces on simplicial meshes, ensuring bounded $L^2$ condition numbers regardless of polynomial degree, which improves numerical stability.

Contribution

It presents the first known construction of frame representations with degree-independent condition numbers for high order finite element methods.

Findings

Constructed frames with bounded $L^2$ condition number independent of polynomial degree

Utilized bubble transform and Jacobi polynomial properties for the construction

Discussed implications for preconditioned iterative methods in finite element analysis

Abstract

The purpose of this paper is to discuss representations of high order $C^0$ finite element spaces on simplicial meshes in any dimension. When computing with high order piecewise polynomials the conditioning of the basis is likely to be important. The main result of this paper is a construction of representations by frames such that the associated $L^2$ condition number is bounded independently of the polynomial degree. To our knowledge, such a representation has not been presented earlier. The main tools we will use for the construction is the bubble transform, introduced previously in [Falk and Winther, Found Comput Math (2016) 16: 297], and properties of Jacobi polynomials on simplexes in higher dimensions. We also include a brief discussion of preconditioned iterative methods for the finite element systems in the setting of representations by frames.

Figures

Tables

Table 1. Table 1. Results for Test 1.

r	3	5	7	9
${\lambda }_{\max }$	6.623	6.819	6.893	6.930
${\lambda }_{\min }^{\ast }$	0.427	0.457	0.472	0.481
cond. num.	15.503	14.922	14.600	14.414
dim. of frame	33	98	199	336
rank of frame	19	46	85	136

Table 2. Table 2. Results for Test 2.

r	3	5	7	9
${\lambda }_{\max }$	6.624	6.819	6.893	6.930
${\lambda }_{\min }^{\ast }$	0.333	0.383	0.414	0.435
cond. num.	19.869	17.809	16.647	15.929
dim. of frame	33	98	199	336
rank of frame	19	46	85	136

Table 3. Table 3. Results for Test 3.

r	3	5	7	9
${\lambda }_{\max }$	6.624	6.819	6.893	6.930
${\lambda }_{\min }^{\ast }$	0.348	0.396	0.427	0.446
cond. num.	19.049	17.203	16.162	15.525
dim. of frame	73	222	455	772
rank of frame	43	106	197	316

Equations

$$a(u,v)=f(v),v\in H^1(\Omega ),$$

$$a(u_h,v)=f(v),v\in V_h.$$

$${\mu }_{h}f\cdot c=\sum _{i=1}^n⟨f,{\varphi }_i⟩c_i=⟨f,{\tau }_h(c)⟩,$$

$${\mathbb{A}}_{h}c={\mu }_h(f_h)\equiv {\tau }_h^{\ast }(f_h),$$

$$\begin{diagram}$$

$$⟨{\mathcal{I}}_{h}f,v{⟩}_{L^2}=⟨f,v⟩,f\in V_h^{\ast },v\in V_h.$$

$$a(u,v)=\lambda ⟨u,v{⟩}_{L^2},u,v\in V_h.$$

$$\kappa ({\mathbb{A}}_h)=\kappa ({\tau }_h^{\ast }{\mathcal{A}}_h{\tau }_h)\ll \kappa ({\mathcal{I}}_h{\mathcal{A}}_h).$$

$$H^1(\Omega ):=\{u\in L^2(\Omega ),\mathrm{grad}u\in L^2(\Omega {)}^d\},$$

$$\text{mathaccent}23H^1(\Omega ):=\{u\in H^1(\Omega ):{\mathrm{tr}}_{\partial \Omega }u=0\}.$$

$${\mathcal{P}}_r({\mathcal{T}}_h):=\{u\in C^0(\Omega ):u{|}_T\in {\mathcal{P}}_r(T),\forall T\in {\mathcal{T}}_h\}.$$

$${\mathbb{A}}_h={\tau }_h^{\ast }{\mathcal{I}}_h^{-1}{\mathcal{I}}_h{\mathcal{A}}_h{\tau }_h={\mathbb{M}}_h({\tau }_h^{-1}{\mathcal{I}}_h{\mathcal{A}}_h{\tau }_h).$$

$$\kappa ({\mathbb{A}}_h)\le \kappa ({\mathcal{I}}_h{\mathcal{A}}_h)\kappa ({\mathbb{M}}_h).$$

$${\lambda }_{\max }({\mathbb{A}}_h)\le {\lambda }_{\max }({\mathcal{I}}_h{\mathcal{A}}_h){\lambda }_{\max }({\mathbb{M}}_h),\text{and }{\lambda }_{\min }({\mathbb{A}}_h)\ge {\lambda }_{\min }({\mathcal{I}}_h{\mathcal{A}}_h){\lambda }_{\min }({\mathbb{M}}_h),$$

$${\lambda }_{\max }({\mathbb{A}}_h)$$

$$\le {\lambda }_{\max }({\mathcal{I}}_h{\mathcal{A}}_h){\lambda }_{\max }({\mathbb{M}}_h),$$

$$a(u,v)={\int }_{-1}^1{u}^{\prime }(x)v^{\prime }(x)dx,$$

$${\int }_{-1}^1{J}_s^{2,2}(x)\cdot J_t^{2,2}(x)(1-x{)}^2(1+x{)}^{2}dx={\delta }_{st},$$

$$S_m^c:=\{\lambda =({\lambda }_0,\cdots ,{\lambda }_m)\in {\mathbb{R}}^{m+1}:\sum _{j=0}^m{\lambda }_j\le 1,{\lambda }_j\ge 0,\forall 0\le j\le m\}.$$

$$S_m:=\{\lambda \in S_m^c:\sum _{j=0}^m{\lambda }_j=1\},$$

$${\Omega }_f^e=\bigcup _{j=0}^m{\Omega }_{x_j}.$$

$${\lambda }_f^{\ast }v(x)=v({\lambda }_f(x)),$$

$$Q_r\left(S_m^c\right):=\left\{v\in {\mathcal{P}}_r\left(S_m^c\right):{\mathrm{tr}}_{\partial S_m^c\backslash S_m}v=0\right\}.$$

$$Q_{f,r}^{\ast }:={\lambda }_f^{\ast }\left(Q_r\left(S_{\dim f}^c\right)\right),m<d.$$

$$\mathfrak{B}:H^1(\Omega )\mapsto \prod _{f\in \Delta }\text{mathaccent}23H^1({\Omega }_f).$$

$$u=\sum _{f\in \Delta }{B}_{f}u,\forall u\in H^1(\Omega ),$$

$$B_f:H^1({\Omega }_f^e)\mapsto \text{mathaccent}23H^1({\Omega }_f)$$

$$B_{f}u:={\lambda }_f^{\ast }\circ K_m\circ A_f{u}^m,\forall f\in {\Delta }_m,$$

$$u^m:=\left(u-\sum _{g\in {\Delta }_j,j<m}{B}_{g}u\right),$$

$$A_{f}u(\lambda )=|{\Omega }_f{|}^{-1}{\int }_{{\Omega }_f}u(G_m(\lambda ,y))dy,$$

$$G_m(\lambda ,y)=y+\sum _{j=0}^m{\lambda }_j(x_j-y).$$

Full text

1. Introduction

The discussion in this paper is motivated by finite element discretizations of second order elliptic equations, where $C^{0}$ piecewise polynomial spaces of high polynomial degree are used as the finite dimensional space. As the polynomial degree increases the choice of basis can have a substantial effect on the conditioning of the linear systems to be solved. The purpose of this paper is to discuss how to obtain representations of the finite element spaces which are uniformly well-conditioned with respect to the polynomial degree. Here the conditioning of the representation is measured by the $L^{2}$ condition number. Furthermore, we will explain how this influences the conditioning of the corresponding discrete systems. Since our main goal is to discuss dependence with respect to the polynomial degree we will consider the mesh $\mathcal{T}_{h}$ to be fixed throughout the discussion below.

To motivate the discussion below, we consider a second order elliptic equation, defined on a bounded domain $\Omega\in\mathbb{R}^{d}$, which admits a weak formulation of the form:

Find $u\in H^{1}(\Omega)$ such that

[TABLE]

where $H^{1}(\Omega)$ denotes the Sobolev space of all functions in $L^{2}$ which also have all first order partial derivates in $L^{2}$. Furthermore, $f$ is a bounded linear functional, and $a$ is a symmetric, bounded, and coercive bilinear form on $H^{1}(\Omega)$. The formulation above reflects that we are considering an elliptic problem with natural boundary condition. If we instead consider problems with an essential boundary condition on parts of the boundary, we will obtain a weak formulation with respect to a corresponding subspace of $H^{1}(\Omega)$. However, the effect of such modifications of (1.1) will have minor effects on the discussion below. Therefore, we will restrict the discussion to problems of the form (1.1) throughout this paper.

A discretization of the problem (1.1) can be derived from a finite dimensional subspace $V_{h}$ of $H^{1}(\Omega)$. In the finite element method $V_{h}$ is typically a space of piecewise polynomials with respect to a partition, or a mesh, $\mathcal{T}_{h}$, with global $C^{0}$ continuity, and where the mesh parameter $h$ indicates the size of the cells of the partition. The corresponding discrete solution is defined by:

Find $u_{h}\in V_{h}$ such that

[TABLE]

This system can alternatively be written as a linear system of the form ${\mathcal{A}_{h}}u_{h}=f_{h}$, where $f_{h}\in V_{h}^{*}$, and where the operator ${\mathcal{A}_{h}}:V_{h}\to V_{h}^{*}$ is defined by ${\mathcal{A}_{h}}u(v)=a(u,v),$ for all $u,v\in V_{h}$. Hence, ${\mathcal{A}_{h}}$ is symmetric in the sense that for all $u,v\in V_{h}$, $\langle{\mathcal{A}_{h}}u,v\rangle=\langle{\mathcal{A}_{h}}v,u\rangle,$ where $\langle\cdot,\cdot\rangle$ is the duality pairing between $V_{h}^{*}$ and $V_{h}$. To turn the discrete system (1.2) into a system of linear equations, written in a matrix/vector form, we need to introduce a basis $\{\phi_{j}\}_{j=1}^{n}$ for the space $V_{h}$. This means that any element $v\in V_{h}$ can be written uniquely on the form $v=\sum_{j}c_{j}\phi_{j}$. We denote the map from $\mathbb{R}^{n}$ to $V_{h}$ given by $c\mapsto v$ for $\tau_{h}$. In a corresponding manner we define $\mu_{h}:V_{h}^{*}\to\mathbb{R}^{n}$ by $(\mu_{h}f)_{i}=\langle f,\phi_{i}\rangle$. We note that if $f\in V_{h}^{*}$ and $c\in\mathbb{R}^{n}$ then

[TABLE]

where $\mathbb{R}^{n}$ is equipped with the standard Euclidean inner product, and where we adopt the standard “dot notation” for this inner product. Hence, $\mu_{h}:V_{h}^{*}\to\mathbb{R}^{n}$ can be identified as $\tau_{h}^{*}$. If $c$ is the unknown vector, $c=\tau_{h}^{-1}u_{h}$, then the system (1.2) is equivalent to the linear system

[TABLE]

where ${{{\mathbb{A}}}_{h}}$ corresponds to the $n\times n$ matrix representing the operator $\tau_{h}^{*}{\mathcal{A}_{h}}\tau_{h}:\mathbb{R}^{n}\to\mathbb{R}^{n}$. The matrix ${{{\mathbb{A}}}_{h}}$ is usually referred to as the stiffness matrix, and the element $({{{\mathbb{A}}}_{h}})_{i,j}$ is given as $a(\phi_{i},\phi_{j})$. Furthermore, we note that the diagram

[TABLE]

commutes. However, there is a striking difference between the operator ${\mathcal{A}_{h}}:V_{h}\to V_{h}^{*}$ and its matrix representation ${{{\mathbb{A}}}_{h}}$. The stiffness matrix ${{{\mathbb{A}}}_{h}}$ depends strongly on the choice of basis, while the operator ${\mathcal{A}_{h}}$ only depends on the bilinear form $a$ and the space $V_{h}$.

For piecewise polynomial spaces of high order the choice of basis can have dramatic effect on the conditioning of the stiffness matrix ${\mathbb{A}}_{h}$. Therefore, there are a number of contributions in the literature discussing how to choose proper bases for $C^{0}$ piecewise polynomial spaces of high order. The purpose of these constructions is usually motivated by the desire to control specific condition numbers or to control the sparsity of the resulting matrices. Let ${\mathcal{I}_{h}}:V_{h}^{*}\to V_{h}$ be the Riesz map given by

[TABLE]

The operator ${\mathcal{I}_{h}}{\mathcal{A}_{h}}:V_{h}\to V_{h}$ is symmetric and positive definite in the $L^{2}$ inner product. This operator is also basis independent, and its eigenvalues are given by the generalized eigenvalue problem

[TABLE]

If $r$ is the polynomial degree then its spectral condition number, $\kappa({\mathcal{I}_{h}}{\mathcal{A}_{h}})$, generally may grow like $r^{4}$, cf. [1, 7, 16]. Therefore, one possiblity is to design a special basis such that the spectral condition number of the stiffness matrix is much smaller than that of the operator ${\mathcal{I}_{h}}{\mathcal{A}_{h}}$, i.e.,

[TABLE]

The obvious constructions which will lead to this is to consider bases which are close to orthonormal with respect to the bilinear form $a$. This approach is taken by Schwab in [16], where integrated Legendre polynomials are used to construct a basis in one space dimension. However, the generalization of this approach to higher dimensions is not obvious. Babuška and Szabó [17] used these basis functions in a tensor product setting to obtain bases for cubical meshes in higher dimensions. In [7] it is established that for such bases one has an estimate for the condition number of the stiffness matrix of the form $\kappa({\mathbb{A}}_{h})\lesssim r^{4(d-1)}$, while $\kappa({\mathbb{M}}_{h})\lesssim r^{4d}$, where ${\mathbb{M}}_{h}$ is the mass matrix given by $({\mathbb{M}}_{h})_{i,j}=\langle\phi_{i},\phi_{j}\rangle_{L^{2}}$ and $d$ is the space dimension. In particular, these estimates indicate that in higher dimensions it may occur that $\kappa({\mathbb{A}}_{h})$ is much larger than its basis independent counterpart, $\kappa({\mathcal{I}_{h}}{\mathcal{A}_{h}})$. Similar constructions on triangular and tetrahedral meshes are based on orthogonal polynomials with respect to triangles and tetrahedrons constructed by the Duffy transform, i.e., mappings between simplexes and cubes, cf. [2, 6, 8, 11, 22]. In a slightly different direction, Xin and Cai [19] used multivariate orthogonal polynomials on simplexes to design $L^{2}$ hierarchical bases for the discontinuous Galerkin method.

Our aim in this paper is to construct representations of $C^{0}$ finite element spaces with $L^{2}$ condition numbers which are bounded independently of the polynomial degree. If the spaces are represented by a basis this quantity can also be characterized as the spectral condition number of the mass matrix. We will argue below that if $\kappa({\mathbb{M}}_{h})$ is well-behaved then the condition number of the stiffness matrix is basically controlled by the basis independent quantity $\kappa(\mathcal{I}_{h}\mathcal{A}_{h})$, cf. (2.2) below. In fact, we will not restrict the discussion below to bases, but instead also allow representations of finite element spaces by frames, i.e., generating sets with redundancy. This more general set up allows us to identify a construction which will lead to representations with $L^{2}$ condition numbers which are independent of the polynomial degree $r$. The key tools we will use for this construction include the properties of the bubble transform, cf. [5]. By combining proper results for the bubble transform with general results for frames based on space decompositions, the construction of frames with $L^{2}$ condition numbers bounded independently of the polynomial degree is reduced to a pure polynomial problem. More precisely, we need to construct Jacobi polynomials on standard simplexes in higher dimensions, and these constructions are well known, cf. [4].

The rest of this paper will be organised as follows. In Section 2, we will present some notation and preliminaries needed for our discussions. In particular, we present a simple bound that relates the spectral condition number of the stiffness matrix and the mass matrix, and we present some elementary numerical examples to illustrate the sharpness of this bound. Furthermore, we will recall the construction of the bubble transform and its properties. Section 3 is devoted to a discussion of frames obtained by a space decomposition, where we focus on general estimates for the appropriate condition numbers. The construction of the specific frames are completed in Section 4, where we explain how to utilize well-conditioned bases on local subdomains. We focus on the construction of $L^{2}$ orthogonal local bases in Section 4 and present explicit forms of the local bases based on Jacobi polynomials on simplexes. We present numerical experiments based on these local constructions which give results that are consistent with the theory. In particular, the results of the experiments indicate that the behavior of the frame condition numbers are robust with respect to perturbations of the mesh, even when the mesh is nearly degenerate. Sections 5 is devoted to the discussion of preconditioned Krylov space methods for the associated frame systems, which in general will be singular. In particular, we make the observation that under the assumption of standard representations of the discrete elliptic operator and the preconditioner, the conditioning of the preconditioned system is, in a proper sense, independent of the choice of basis or frame. On the other hand, the representation will substantially effect the individual conditioning of the stiffness matrix and the matrix representation of the preconditioner. However, we will argue that as long as the $L^{2}$ condition number of the representation stays bounded, then these matrices will roughly behave like their basis independent counterparts. Finally, some concluding remarks are given in §6.

2. Notation and preliminaries

We assume that $\Omega$ is a bounded polyhedral domain in $\mathbb{R}^{d}$. We recall the definition of the Sobolev spaces

[TABLE]

and the corresponding subspace of functions with vanishing trace:

[TABLE]

We assume that $\mathcal{T}_{h}$ is a simplicial mesh on $\Omega$. If $S$ is a subset of $\mathbb{R}^{d}$ we let $\mathcal{P}_{r}(S)$ denote the polynomials with degrees less than or equal to $r$ on $S$, while the corresponding space of $C^{0}$-piecewise polynomials with respect to the mesh $\mathcal{T}_{h}$ is denoted ${\mathcal{P}_{r}}(\mathcal{T}_{h})$, i.e.,

[TABLE]

2.1. Representation of discrete operators

Consider a finite element system of the form (1.3), where the space $V_{h}={\mathcal{P}_{r}}(\mathcal{T}_{h})$ for a suitable $r\geq 1$, and let $\{\phi_{j}\}_{j=1}^{n}$ be any basis for $V_{h}$. We recall that the $n\times n$ stiffness matrix ${\mathbb{A}}_{h}$ is given as $\tau_{h}^{*}\mathcal{A}_{h}\tau_{h}$, where $\mathcal{A}_{h}:V_{h}\to V_{h}^{*}$ is the discrete elliptic operator defined by the variational problem (1.2). The corresponding mass matrix is the $n\times n$ matrix with elements $\langle\phi_{i},\phi_{j}\rangle_{L^{2}}$. Alternatively, ${\mathbb{M}}_{h}=\tau_{h}^{*}\mathcal{I}_{h}^{-1}\tau_{h}$, where we recall that $\mathcal{I}_{h}$ is the Riesz map, mapping $V_{h}^{*}$ to $V_{h}$. The matrices ${\mathbb{A}}_{h}$ and ${\mathbb{M}}_{h}$ are related by the relation

[TABLE]

We note that the operator $\tau_{h}^{-1}\mathcal{I}_{h}\mathcal{A}_{h}\tau_{h}$ is similar to the basis independent operator $\mathcal{I}_{h}\mathcal{A}_{h}$. A direct consequence of the identity (2.1), using the characterization of the extreme eigenvalues by the Raleigh quotient, is the inequality

[TABLE]

where the spectral condition number $\kappa({\mathcal{I}_{h}}{\mathcal{A}_{h}})$ is basis independent, while $\kappa({\mathbb{M}}_{h})$ is independent of the underlying elliptic operator. In fact, (2.2) follows from the stronger properties

[TABLE]

where $\lambda_{\min}$ and $\lambda_{\max}$ denote the extreme eigenvalues. To see this just observe that

[TABLE]

and a similar argument establishes the corresponding inequality for $\lambda_{\min}$.

We can therefore conclude that if the basis is chosen such that $\kappa({\mathbb{M}}_{h})$ is properly bounded, then the conditioning of the stiffness matrix ${\mathbb{A}}_{h}$ is no worse than its basis independent analog. To illustrate the effect on the conditioning of the stiffness matrix ${\mathbb{A}}_{h}$, by controlling the $L^{2}$ condition number of the basis, we present some simple numerical examples in one space dimension. In other words, we are testing the sharpness of the bound (2.2) in the simplest possible setting.

Example. We consider the Laplace problem in one space dimension, with homogeneous Dirichlet boundary conditions, i.e., the bilinear form $a$ is given by

[TABLE]

and we use a mesh consisting of one interval. Therefore, the corresponding spaces $V_{h}$ will be given as $V_{h}=\mathaccent 23{\mathcal{P}}_{r}(\Omega)$, where $\Omega$ is the interval $(-1,1)$ and $\mathaccent 23{\mathcal{P}}_{r}(\Omega)$ is the space of polynomials of degree on $\Omega$ which vanish on the boundary. We will investigate the effect of choosing three different bases $\{\phi_{1},\phi_{2},\cdots,\phi_{r-1}\}$ for the spaces $V_{h}$, by computing the condition numbers of the mass matrix $\mathbb{M}_{h}$, the stiffmess matrix $\mathbb{A}_{h}$ and the condition number of the basis independent operator $\mathcal{I}_{h}\mathcal{A}_{h}$. In fact, for any basis the latter is equal to $\kappa(\mathbb{M}_{h}^{-1}\mathbb{A}_{h})$.

Our first test is based on an $L^{2}$ orthonormal basis. We consider the polynomials $(1-x)(1+x)J^{2,2}_{r}(x),~{}r=0,1,\cdots$, where $J_{r}^{2,2}(x)$ is the orthonormal Jacobi polynomials on $[-1,1]$ with respect to the weight $(1-x)^{2}(1+x)^{2}$, i.e.,

[TABLE]

cf. Appendix A below. Since these polynomials form an orthonormal basis for $V_{h}$ $\kappa({\mathbb{M}}_{h})=1$ for all $r$ and $\kappa(\mathbb{M}_{h}^{-1}\mathbb{A}_{h})=\kappa({\mathbb{A}}_{h})$. The logarithms of the latter, for increasing values of $r$, are shown in Figure 2, while $\log(\kappa(\mathbb{A}_{h}))$ are compared to $\log(r)$ in Figure 2. The growth depicted here is consistent with the asympotic upper bound, $\kappa(\mathcal{I}_{h}\mathcal{A}_{h})\lesssim r^{4}$, cf., [1].

Next we consider the Bernstein basis. More precisely, for any $r\geq 2$ consider the functions $b_{s,r}\left(\frac{x+1}{2}\right)$, where $b_{s,r}(x)={r\choose s}x^{s}(1-x)^{r-s},\quad 1\leq s\leq r-1.$ The condition numbers of the mass and stiffness matrices are shown in Figure 4, and the results indicate that both $\kappa({\mathbb{M}}_{h})$ and $\kappa({\mathbb{A}}_{h})$ grow exponentially with $r$. This is consistent with the explicit formula, $\kappa(M_{h})=\sqrt{2r+1\choose r}$, which holds in the case of no boundary condidtions [3, 12]). Furthermore, we observe that $\kappa({\mathbb{A}}_{h})$ is several magnitudes larger than the basis independent quantity $\kappa(\mathcal{I}_{h}\mathcal{A}_{h}).$

Finally, we consider the corresponding power basis $(1+x)(1-x)x^{r-2},~{}~{}r=2,3,\cdots$. The condition numbers of the mass and stiffness matrices are in this case even much larger than for the Bernstein basis, cf. Figure 4. Due to the extremely bad condition numbers, the computations are only reliable for small values of $r$.

The experiments just presented illustrate an effect of the bound (2.2). For the Jacobi basis the condition number of the mass matrix is well controlled, and as a consequence the growth of $\kappa({\mathbb{A}}_{h})$ is moderate. On the other hand, for the two other bases $\kappa({\mathbb{M}}_{h})$ and $\kappa({\mathbb{A}}_{h})$ both grows much faster than $\kappa(\mathcal{I}_{h}\mathcal{A}_{h}).$

2.2. The local spaces $Q_{f,r}^{*}$

The rest of this paper will mostly be devoted to the construction of representations for the spaces ${\mathcal{P}_{r}}(\mathcal{T}_{h})$ which admit $L^{2}$ condition numbers which are independent of the polynomial degree $r$. To obtain such a result we will not restrict the discussion to representations by bases, but we will allow more general representations by frames. In particular, our construction will rely on results for the bubble transform derived in [5]. To present these results, and to explain how they will be used here, we will first introduce some additional notation.

If $\mathcal{T}_{h}$ is simplicial triangulation of $\Omega$, we let $\Delta_{m}=\Delta_{m}(\mathcal{T}_{h})$ be the set of all the subsimplexes of $\mathcal{T}_{h}$ of dimension $m$, while $\Delta=\cup_{m=0}^{d}\Delta_{m}$ contains all the subsimplexes. If $T\in\mathcal{T}_{h}$ we let $\Delta(T)$ be the set of all subsimplexes of $T$. For $f\in\Delta_{m}$, the local patch, or macroelement, $\Omega_{f}=\cup\{T\in\Delta\,:\,f\in\Delta(T)\}$ is the union of all the elements of the mesh which contains $f$. Furthermore, $\mathcal{T}_{f,h}$ is the partition $\mathcal{T}_{h}$ restricted to $\Omega_{f}$.

When $x_{j}\in\Delta_{0}$ is a vertex, we use $\lambda_{j}(x)$ to denote the piecewise linear function which equals one at $x_{j}$, and equals zero at other vertices. From another point of view, $\lambda_{j}$ is the barycentric coordinate associated to $x_{j}$, extended by zero outside the macroelement $\Omega_{x_{j}}$. For $m<d$ and $f=[x_{0},x_{1},\dots,x_{m}]\in\Delta_{m}$, i.e. the convex hull of vertices $x_{0},x_{1},\cdots,x_{m}$, we use $\lambda_{f}$ to denote the vector field $(\lambda_{0},\lambda_{1},\ldots\lambda_{m})$. Following the approach taken in [5] we will consider $\lambda_{f}$ as a mapping from the domain $\Omega$ to the standard simplex $S_{m}^{c}$ in $\mathbb{R}^{m+1}$ given by

[TABLE]

We also define $S_{m}$ to be the face of $S_{m}^{c}$ opposite the origin, i.e,

[TABLE]

such that $S_{m}^{c}=\left[0,S_{m}\right]$, i.e., $S_{m}^{c}$ is the set of all convex combinations of the origin and elements of $S_{m}$. The mapping $\lambda_{f}$, restricted to $\Omega_{f}$, is surjective but not injective (see Figure 5 for the case $m=1$).

If $f=[x_{0},x_{1},\ldots,x_{m}]\in\Delta_{m}$ then the associated macroelement $\Omega_{f}$ can be characterized as $\Omega_{f}=\cap_{j=0}^{m}\Omega_{x_{j}}$, while the corresponding extended macroelement, $\Omega_{f}^{e}\supset\Omega_{f}$, is defined by

[TABLE]

The pull back of the extended barycentric coordinates, $\lambda_{f}^{\ast}$, given by

[TABLE]

maps functions on $S_{m}^{c}$ to functions on $\Omega$ which are constant and equal to $v(0)$ outside $\Omega_{f}^{e}$. Furthermore, if ${\operatorname{tr}}_{\partial S_{m}^{c}\setminus S_{m}}v=0$ then $\lambda_{f}^{*}v$ vanishes on the boundary of $\Omega_{f}$.

The space of polynomials of degree less than equal to $r$ which vanish on $\partial S_{m}^{c}\setminus S_{m}$ will be denoted $Q_{r}\left(S_{m}^{c}\right)$, i.e.,

[TABLE]

By applying the pullback, $\lambda_{f}^{*}$, to this polynomial space we obtain

[TABLE]

The elements of this space are polynomials in the variables $\lambda_{0}(x),\lambda_{1}(x),\ldots,\lambda_{m}(x)$, and they vanish on the boundary of $\Omega_{f}$. In other words, the space ${Q_{f,r}^{\ast}}$ can be identified with a subspace of $\mathaccent 23{\mathcal{P}_{r}}(\mathcal{T}_{f,h})$, the subspace of ${\mathcal{P}_{r}}(\mathcal{T}_{f,h})$ which consists of functions which vanish on $\partial\Omega_{f}$. In fact, in the special case when $m=d$, i.e., when $f\in\Delta_{d}=\mathcal{T}_{h}$ we define the ${Q_{f,r}^{\ast}}$ to be equal to $\mathaccent 23{\mathcal{P}_{r}}(f)$. Alternatively, if we define $Q_{r}\left(S_{m}^{c}\right)$ to be $\mathaccent 23{\mathcal{P}_{r}}(S_{m})$ for $m=d$, then the identification (2.3) also holds in this case. The local spaces ${Q_{f,r}^{\ast}}$ will act as key building blocks in our construction below.

2.3. The bubble transform

The bubble transform is a map that depends on the mesh $\mathcal{T}_{h}$, but no piecewise polynomial space occurs in the definition. In particular, it does not depend on a degree parameter $r$. In [5] the construction of the bubble transform was partly motivated by the desire to design local projections onto the piecewise polynomial spaces $\mathcal{P}_{r}(\mathcal{T}_{h})$ with proper bounds independent of $r$. In this paper, we will utilize the properties of the bubble transform to construct frames for the spaces $\mathcal{P}_{r}(\mathcal{T}_{h})$ that admit $L^{2}$ condition numbers which are independent of the degree parameter $r$.

The bubble transform is a map $\mathfrak{B}=\mathfrak{B}_{\mathcal{T}_{h}}$ of the form

[TABLE]

It is a tool to decompose an $H^{1}$ function defined on $\Omega$ into components $B_{f}u$ with local support in $\Omega_{f}$. More precisely,

[TABLE]

where

[TABLE]

gives the component of $u$ which is supported on $\Omega_{f}$. In particular, we observe that the operator $B_{f}$ is local in the sense that $B_{f}u$ depends on $u|_{\Omega_{f}^{e}}$. Another key property of the map $\mathcal{B}$ is that it is invariant with respect to the piecewise polynomial spaces $\mathcal{P}_{r}(\mathcal{T}_{h})$, i.e., if $u\in\mathcal{P}_{r}(\mathcal{T}_{h})$ then $B_{f}u\in\mathaccent 23{\mathcal{P}}_{r}(\mathcal{T}_{f,h})$.

The bubble transform has a recursive definition. We briefly recall its construction, but for more details we refer to [5]. For $m=0,1,\dots,d-1$, $B_{f}u$ is of the form

[TABLE]

where

[TABLE]

while $B_{T}u=u^{d}|_{T}$ if $T\in\Delta_{d}=\mathcal{T}_{h}$. The pull back $\lambda_{f}^{*}$, mapping functions on $S_{m}^{c}$ to functions on $\Omega$, is discussed above. The operator $A_{f}$ is an average operator, while $K_{m}$ is refereed to as a cut-off operator. If $f=[x_{0},x_{1},\ldots,x_{m}]$ then for any $\lambda\in S_{m}^{c}$

[TABLE]

where $G_{m}:S_{m}^{c}\times\Omega_{f}\to\Omega_{f}$ is given by

[TABLE]

The operator $A_{f}$ maps a function $u$ defined on $\Omega_{f}$ to a function $A_{f}u$ defined on $S_{m}^{c}$, and it is a smoothing operator in the sense that for any $u\in L^{2}(\Omega)$ the function $A_{f}u$ will have point values away from the simplex $S_{m}$. On the other hand, ${\operatorname{tr}}_{S_{m}}A_{f}u$ corresponds exactly to ${\operatorname{tr}}_{f}u$. Furthermore, the operator $A_{f}$ has the property that it maps piecewise polynomials into polynomials. More precisely, if $u\in\mathcal{P}_{r}(\mathcal{T}_{f,h})$ then $A_{f}u\in\mathcal{P}_{r}(S_{m}^{c})$.

The cut-off operator $K_{m}$ maps the set of functions defined on $S_{m}^{c}$ to itself. Its key property is that it preserves the trace on $S_{m}$, i.e., ${\operatorname{tr}}_{S_{m}}K_{m}v={\operatorname{tr}}_{S_{m}}v$, while it kills the trace on the rest of the boundary of $\partial S_{m}^{c}$, i.e, ${\operatorname{tr}}_{\partial S_{m}^{c}\setminus S_{m}}K_{m}v=0$. In addition, $K_{m}$ is polynomial preserving in the sense that if $v\in\mathcal{P}_{r}(S_{m}^{c})$, with ${\operatorname{tr}}_{\partial S_{m}}v=0$, then $K_{m}v\in\mathcal{P}_{r}(S_{m}^{c})$.

The key properties of the bubble transform are stated in [5, Section 4]. For the convenience of the readers, and since most of these properties are essential for the discussion below, we summerize the main results here.

1. (1)

The construction using barycentric coordinates: For $f\in\Delta_{m}$ and $0\leq m<d$

[TABLE]

where $u^{m}$ is defined by (2.4), while $B_{f}u=u^{d}|_{f}$ if $f\in\Delta_{d}$. 2. (2)

The boundedness in $L^{2}$ and $H^{1}$:

[TABLE]

where $b$ is a generic constant not depending on the function $u$. 3. (3)

Partition of unity: If $u\in L^{2}(\Omega)$ then $\sum_{f\in\Delta}B_{f}u=u,$ 4. (4)

Local support: If $u\in H^{1}(\Omega)$ then $B_{f}u\in\mathaccent 23{H}^{1}(\Omega_{f})$. 5. (5)

Polynomial preserving: If $u\in{\mathcal{P}_{r}}(\mathcal{T}_{h})$ then $B_{f}u\in{Q_{f,r}^{\ast}}.$

We also note that if $u=\sum_{f\in\Delta}u_{f}$, where $\mathrm{supp}(u_{f})\subset\Omega_{f}$ then

[TABLE]

where the constant $a$ only depends on the number of overlaps of the subdomains $\Omega_{f}$. Therefore, by combining this with (2.5) we obtain that the norms $\|u\|_{L^{2}}$ and $(\sum_{f\in\Delta}\|B_{f}u\|_{L^{2}}^{2})^{1/2}$ are equivalent.

The bubble transform suggests a decomposition of the finite element spaces ${\mathcal{P}_{r}}(\mathcal{T}_{h})$ of the form

[TABLE]

In fact, this decomposition follows directly from the properties above. The spaces $Q_{f,r}^{*}$ are local spaces consisting of piecewise polynomials with support on $\Omega_{f}$. On the other hand, the sum above is in general not direct. To see this, we observe that that if $y$ is a vertex, i.e., $y\in\Delta_{0}$ then

[TABLE]

where $\lambda_{y}$ is the extended barycentric coordinate of the vertex $y$. In particular, the function $u(x)=\lambda_{y}(x)(1-\lambda_{y}(x))$ is an element of $Q_{y,r}^{\ast}$ for $r\geq 2$. Let $x_{1},x_{2},\ldots x_{k}$ be the other vertices in $\Delta_{0}(\mathcal{T}_{f,h})$ with correponding exteded barycentric coordinates $\lambda_{1},\lambda_{2},\ldots,\lambda_{k}$. Then the function $u$ can alternatively be expressed as

[TABLE]

Furthermore, $\lambda_{y}\lambda_{j}\in Q_{f_{j},r}^{\ast}$, where $f_{j}=[y,x_{j}]\in\Delta_{1}$. We conclude that the function $u$ is both in $Q_{y,r}^{\ast}$ and $\sum_{j}Q_{f_{j},r}^{\ast}$, and therefore the sum (2.7) is not direct. Similar redundancies also appear for the spaces $Q_{f,r}^{\ast}$ for simplices of higher dimensions. Therefore, if we want to utilize the decomposition (2.7), and bases for the local spaces $Q_{f,r}^{\ast}$, to represent the functions in the spaces ${\mathcal{P}_{r}}(\mathcal{T}_{h})$ we are forced to study representations by frames.

Interpreting the properties of the bubble transform for the decomposition (2.7) of ${\mathcal{P}_{r}}(\mathcal{T}_{h})$ the following result is obtained.

Theorem 2.1**.**

The decomposition (2.7) is stable in the sense that there exists $B_{f}:{\mathcal{P}_{r}}(\mathcal{T}_{h})\to{Q_{f,r}^{\ast}},~{}\forall f\in\Delta$, and a positive constant $b$ such that

[TABLE]

Furthermore, as a result of the finite overlapping property of the mesh topology, there exists a positive constant $a$ such that

[TABLE]

3. Estimates of frame condition numbers

We will utilize the decomposition (2.7) to obtain a well conditioned representation of functions in the space ${\mathcal{P}_{r}}(\mathcal{T}_{h})$. More precisely, we will combine the decomposition (2.7) with a basis for each of the spaces $Q_{f,r}^{\ast}$. Due to the redundancy of the decomposition (2.7) this will lead to a spanning set for the functions in ${\mathcal{P}_{r}}(\mathcal{T}_{h})$ where the elements are not linearly independent, i.e., we obtain representations by frames, cf. [14]. Therefore, we will give a brief review of representations by frames. In particular, we will discuss frames obtained from space decompositions.

Throughout this section we use $W$ to denote a real, finite dimensional Hilbert space with inner product $\langle\cdot,\cdot\rangle_{W}$. Roughly speaking, a frame is a set of generators which allow redundancy. In other words, if $\Phi=\{\phi_{1},\phi_{2},\cdots\}$ then each element $u\in W$ can be expressed as a linear combination $u=\sum_{k}c_{k}\phi_{k}$, but this representation is in general not unique. The condition number of the frame $\Phi$, ${\mathcal{K}}(\Phi)$, is defined as

[TABLE]

where

[TABLE]

In other words, $\alpha$ and $\beta$ are the optimal constants such that the bounds

[TABLE]

holds. Therefore, ${\mathcal{K}}(\Phi)$ is the natural concept to relate the norm of $u$ to the norm of its coefficients measured in $l^{2}$.

Remark 1**.**

If $\Phi$ is a basis then it is well known that ${\mathcal{K}}(\Phi)$ is equal to the spectral condition number of the corresponding “mass matrix,” with elements $\langle\phi_{i},\phi_{j}\rangle_{W}$. In fact, the parameters $\alpha$ and $\beta$, given in (3.1), are exactly the smallest and largest eigenvalue of the mass matrix. In the case of a frame, the mass matrix will in general be singular. In this case $\beta$ is still the largest eigenvalue of the mass matrix, while $\alpha$ is the smallest positive eigenvalue.

In fact, a similar characterization can be given in the case of frames, cf. Section 5.2 below.

3.1. Frames based on space decomposition

To give a general description of frames based on space decomposition, we assume that the space $W$ admits a decomposition of the form

[TABLE]

where $W_{j},j=1,\cdots,J$ are subspaces of $W$. The decomposition is not assumed to be direct, but we assume that there exists a positive constant $a$ such that for any $u=\{u_{j}\}\in\prod_{j=1}^{J}W_{j}$, we have

[TABLE]

Furthermore, we assume that there is a positive constant $b$ such that all $u\in W$ admits a decomposition $u=\sum_{j=1}^{J}u_{j}$, $u_{j}\in W_{j}$, where

[TABLE]

Of course, due to (2.5) and (2.6), these bounds are known to hold for the spaces ${\mathcal{P}_{r}}(\mathcal{T}_{h})$ with $L^{2}$ inner product. We will use the decomposition (3.2) to define a frame for $W$. More precisely, for each $j$ let $\{\phi_{j,k}\}_{k}$ be a basis for the space $W_{j}$. The frame $\Phi$ is then given as $\{\phi_{j,k}\}_{1\leq j\leq J,1\leq k\leq N_{j}}$, where $N_{j}$ is the dimension of $W_{j}$. For each $j$ we assume that $0<\alpha_{j}<\beta_{j}$ are the optimal constants such that

[TABLE]

In other words, $\mathcal{K}_{j}=\beta_{j}/\alpha_{j}$ is the condition number for the basis $\{\phi_{j,k}\}_{k}$ of $W_{j}$.

Theorem 3.1**.**

Assume that the decomposition (3.2) satisfies (3.3), (3.4) and (3.5). The frame $\Phi=\{\phi_{j,k}\}_{j,k}$ introduced above satisfies

[TABLE]

Proof.

Let $\alpha$ and $\beta$ be the two constants defined by (3.1). We will show that

[TABLE]

From these bounds we immediately obtain

[TABLE]

which is the desired bound. Therefore, it is enough to establish the bounds given by (3.6).

To show the first inequality, let $u$ be any element in $W$, and $u=\sum_{j}u_{j}$ a decomposition of the form (3.2) satisfying (3.4). Furthermore, let $c=\{c_{j,k}\}$ be the unique coefficients such that $u_{j}=\sum_{k}c_{j,k}\phi_{j,k}$. Then

[TABLE]

This implies

[TABLE]

and the first inequality of (3.6) follows by taking infimum over all elements of $u$ of $W$.

On the other hand, for any coefficient $c=\{c_{j,k}\}$, we define $u_{j}=\sum_{k}c_{j,k}\phi_{j,k}$ and $u=\sum_{j}u_{j}$. We then have

[TABLE]

Therefore, for any $u\in W$, we have

[TABLE]

and hence the second bound of (3.6) follows. ∎

The result above shows that the condition number of the frame $\Phi$ is bounded by the constants $a$ and $b$, derived from the decompostion (3.2) and the local condition numbers $\mathcal{K}_{j}$. In addition, the factor $\max_{j,k}\alpha_{j}/\alpha_{k}$, which we will refer to as a scaling factor, appears. This factor will be small if all the local condition numbers $\mathcal{K}_{j}$ are small, and if the local bases $\{\phi_{j,k}\}_{k}$ in addition are scaled similarly. In fact, the appearance of this factor is similar to a well known phenomenon. Consider a block diagonal matrix of the form

[TABLE]

where $I_{1}$ and $I_{2}$ are identity matrices of proper dimensions, and where $\epsilon>0$ is a real parameter. Then each block has condition number $1$, while the full matrix has condition number $\epsilon^{-1}$ due to the different scaling of the two blocks.

3.2. The bubble decomposition of ${\mathcal{P}_{r}}(\mathcal{T}_{h})$

We end this section by applying Theorem 3.1 to the decomposition (2.7) of the spaces ${\mathcal{P}_{r}}(\mathcal{T}_{h})$. For each $f\in\Delta$ let $N_{f}$ denote the dimension of the space ${Q_{f,r}^{\ast}}$. We will see below that it is possible to construct a basis for each of the spaces ${Q_{f,r}^{\ast}}$ such that all are well conditioned in $L^{2}$, and with a comparable scaling. Therefore, consider the set up when we have a basis for each of the spaces ${Q_{f,r}^{\ast}}$ of the form $\Phi_{f}=\phi_{f,1},\cdots,\phi_{f,N_{f}}$, satisfying

[TABLE]

where the positive constants $\alpha_{0}$ and $\beta_{0}$ are independent of $f\in\Delta$. By combining Theorem 3.1 with Theorem 2.1 we immediately obtain the following.

Corollary 3.1**.**

Let $\Phi=\{\Phi_{f}\}_{f\in\Delta}$ be the frame representation of the space ${\mathcal{P}_{r}}(\mathcal{T}_{h})$ just introduced, and satisfying (3.7). We have the estimate:

[TABLE]

where $a$ and $b$ are the constants appearing in (2.5) and (2.6).

Proof.

It is a consequence of (3.7) that the condition number of each local basis, $\Phi_{f}$ of ${Q_{f,r}^{\ast}}$, is bounded by $\alpha_{0}^{-1}\beta_{0}$, and that the same bound holds for the scaling factor appearing in Theorem 3.1. The result therefore follows from this theorem. ∎

Remark 2**.**

We note that in the special case when each of the local bases $\Phi_{f}=\phi_{f,1},\cdots,\phi_{f,N_{f}}$ is orthonormal then the bound above reduces to $\mathcal{K}(\Phi)\leq a^{-1}b$, i.e, the condition number of the frame is bounded entirely by the two constants given in (2.5) and (2.6).

4. Construction of bases for the local spaces

Based on the discussion above, cf. Corollary 3.1, we can conclude that to obtain a well-conditioned frame for the spaces ${\mathcal{P}_{r}}(\mathcal{T}_{h})$, it is enough to construct bases for the local spaces $Q_{f,r}^{*}$ which are uniformly well-conditioned in $L^{2}$. More precisely, it is enough to construct bases $\Phi_{f}$ for the spaces $Q_{f,r}^{*}$ such that condition (3.7) holds. In the special case when $\dim f=d$ then $Q_{f,r}^{*}=\mathaccent 23{\mathcal{P}_{r}}(f)$, and the construction of a basis for this space is well known. We return to this case in the Appendix below. When $\dim f<d$ we recall from Section 2 that $Q_{f,r}^{*}$ is defined by a pull back, with respect to the map $\lambda_{f}:\Omega\to S_{m}^{c}$, of the polynomial space $Q_{r}(S_{m}^{c})$, where $S_{m}^{c}$ is a reference simplex in $\mathbb{R}^{m+1}$ and $f\in\Delta_{m}(\mathcal{T}_{h})$. Therefore, we will construct a basis for the space $Q_{f,r}^{*}$ by utilizing a basis for $Q_{r}(S_{m}^{c})$.

4.1. Construction of local bases

Element of the space $Q_{r}(S_{m}^{c})$ is of the form $\lambda_{0}\cdot\lambda_{1}\cdots\lambda_{m}\,p\equiv(\Pi\lambda)_{m}\,p$, where $p\in\mathcal{P}_{r-m-1}(S_{m}^{c})$. Therefore, any basis for the space $\mathcal{P}_{r-m-1}(S_{m}^{c})$ leads to a corresponding basis for $Q_{r}(S_{m}^{c})$. To proceed we recall the fact from [5, formula (5.6)], that if $\phi$ is any smooth function on $S_{m}^{c}$ then

[TABLE]

where $b(\lambda)=1-\sum_{j=0}^{m}\lambda_{i}$, and $c_{f}$ is a scaling factor depending on the geometry of the macroelement $\Omega_{f}$. In fact, in the notation of [5, Section 5] the constant $c_{f}$ is given by

[TABLE]

where $f^{*}$ is a piecewise flat manifold of dimension $d-1-\dim f$ contained in $\partial\Omega_{f}$, and $J(f,q)$ is a piecewise constant function on $f^{*}$. However, for the discussion here it suffices to observe that $c_{f}=\mathcal{O}(h_{f}^{d})$, uniformly with respect to a family of shape regular meshes. Here $h_{f}$ is a local parameter representing the diameter of $T\in\mathcal{T}_{h}$, i.e. $h_{f}$ represents the size of the elements contained in $\Omega_{f}$. As a consequence, if $u$ and $v$ are orthogonal functions on $S_{m}^{c}$, with respect to the weight functions $b(\lambda)^{d-m-1}$, then $\lambda_{f}^{*}u$ and $\lambda_{f}^{*}v$ are $L^{2}$ orthogonal functions on $\Omega_{f}$. Alternatively, if $p$ and $q$ are elements of $\mathcal{P}_{r-m-1}(S_{m}^{c})$, which are orthogonal with respect to the weight function $w_{m}(\lambda):=(\Pi\lambda)_{m}^{2}b(\lambda)^{d-m-1}$, then the corresponding functions $\lambda_{f}^{*}[(\Pi\lambda)_{m}p]$ and $\lambda_{f}^{*}[(\Pi\lambda)_{m}q]$ are $L^{2}$ orthogonal functions belonging to the space $Q_{f,r}^{*}$. Furthermore, the norm of $\lambda_{f}^{*}[(\Pi\lambda)_{m}p]$ in $L^{2}(\Omega_{f})$ is equivalent to $h_{f}^{d/2}$ times the corresponding weighted $L^{2}$ norm of $p$ on $S_{m}^{c}$. Therefore, the problem of constructing $L^{2}$ orthogonal and uniformly scaled bases for the local spaces $Q_{f,r}^{*}$, is equivalent to the construction of bases for the polynomial spaces $\mathcal{P}_{r-m-1}(S_{m}^{c})$, which are orthogonal with respect to the weight function $w_{m}$, and uniformly scaled. Actually, since the scaling factor $c_{f}$ only depends on $f$ and is uniform for all the bases associated with $f$, any orthonormal bases on $\mathcal{P}_{r-m-1}(S_{m}^{c})$ with the weight $w_{m}$ will transform to bases of $Q_{r}(S_{m}^{c})$ with condition number one.

To construct a polynomial basis, which is orthogonal with respect to a polynomial weight function, corresponds to the study of Jacobi polynomials. Single variate orthogonal polynomials are of course well studied, but there are also explicit formulas for Jacobi polynomials with respect to simplexes in higher dimensions. The most popular approach to construct orthogonal polynomials in higher dimensions is to use a transform between simplexes and cubes, referred to as “the Duffy transform” or “the Koorwinder method” [9]. Furthermore, hierarchical constructions of orthogonal polynomial can be found in [4].

4.2. Numerical results

As an example, we present explicit formulas of our frames and explore the condition numbers of the matrices by numerical computation. This will verify our theoretical results and show that the frame is well-conditioned and the constants $a$ and $b$ in (3.8) are bounded on various regular or irregular meshes. We will only consider the two dimensional case, i.e., the space dimension $d$ is equal to $2$.

Let $J^{\alpha,\beta}_{s}(x)$ be the standard single variate Jacobi polynomial of degree $s$ with weight $(\alpha,\beta)$ defined on $[-1,1]$ (see also Appendix below). Define

[TABLE]

so that $\{\tilde{J}_{s}^{\alpha,\beta}(\xi)\},s=0,1,\cdots$ is a single variate orthogonal basis on $[0,1]$ with weight $\omega^{\alpha,\beta}:=(1-\xi)^{\alpha}\xi^{\beta}$. The explicit formula for our frames are derived from the corresponding basis for the spaces $Q_{r}\left(S_{m}^{c}\right)$ for $m=0,1,2$, cf. (2.3), where we recall that $Q_{r}\left(S_{2}^{c}\right)=\mathaccent 23{\mathcal{P}_{r}}(S_{2})$. Let $\lambda_{0},\lambda_{1},\cdots,\lambda_{m}$ be the barycentric coordinates of $S_{m}$ extended to $S_{m}^{c}$. The explicit bases for the spaces $Q_{r}\left(S_{m}^{c}\right)$ can be given as follows:

[TABLE]

for the three cases $m=0,1,2$. In fact, these functions are special cases of the simplicial orthogonal polynomials (6.2), cf. the Appendix below. The corresponding bases for the spaces $Q_{f,r}^{\ast}$ are given by (2.3) for $\dim f=m$. We note that the functions in (4.3) and (4.4) do not have rotational symmetry, meaning that exchanging two barycentric coordinates in the formulas will lead to different bases. In the assembling of matrices permuting the variables in (4.3) may facilitate the match of orientation.

To avoid considering effects from inaccurate numerical quadrature, in the computation we use Gauss type formulas on triangles with order 20 [23]. Since in the tests below we only consider polynomials of degree less than 10, the numerical quadrature will lead to exact integration. Due to the roundoff error, the mass matrices of the frames do not have precise zero eigenvalues. Therefore we will use a tolerance $5e{-11}$ to select nonzero eigenvalues and thus will consider all numbers near the machine precision as zero.

The bases of $Q_{f,r}^{\ast}$ defined above are $L^{2}$ orthogonal on each triangle of $\Omega_{f}$. Therefore the block of the mass matrix corresponding to contributions from a single macroelement is diagonal on each triangle. This is also observed in the numerical tests below. However, bases associated to different macroelements can still interfere with each other and the contribution of this overlapping effect to the condition number is reflected in the bound (3.8) with constants $a$ and $b$.

In the numerical tests below we consider the mass matrix in the case when the complete mesh is a single vertex macroelement. Alternatively, this matrix can be thought of as a ”local block matrix” of a mass matrix generated by a larger mesh. On this local mesh we consider all the frame functions generated by all the subsimplexes of the mesh. In all the examples below the local mesh corresponds to the macroelement of the origin, and we will consider three cases:

• •

Test 1: a macroelement consisting of three triangles with boundary vertices (1, -1), (0, 2), (-1, -1) (left of Figure 6);

• •

Test 2: distortion of Test 1, which has three triangles with boundary vertices (1.3, -0.0001), (0, 2.12), (-1, -0.0001) (right of Figure 6);

• •

Test 3: a macroelement consisting of six triangles with boundary vertices $(2,0)$, $(1,1)$, $(0,1)$, $(-1,1)$, $(-2,0)$, $(-1,-1)$, $(1,-1)$ (Figure 7).

The diagonal entries of the mass matrix are scaled to one such that the bases are scaled to have unit $L^{2}$ norms. In the tables below, dimension of the frame means the number of functions in the frame representation, while rank of the frame means the rank of the mass matrix. Alternatively, the rank of the frame can be identified as the dimension of the space of continuous piecewise polynomial of the same degree on the same mesh.

The condition numbers in Test 1 are shown in Table 1. Here $\lambda_{\min}^{\ast}$ denotes the minimal nonzero eigenvalue. Due to the roundoff error, the minimal nonzero eigenvalue slightly decreases as $r$ increases, in contrast to the monotonicity predicted by the min-max principle. Thus the condition numbers also slightly decrease.

The results on a distorted mesh, Test 2, are shown in Table 2. The condition numbers are only slightly larger than those in Test 1, showing that the condition numbers are nicely bounded even if there are very thin elements present in the mesh.

The results for Test 3 are shown in Table 3. Similar to Test 1 and 2, the condition number also remains bounded as the degree increases.

Patterns of the mass matrices, i.e., the zero–nonzero structure, are shown in Figure 9 and Figure 9 for Test 1. The pattern only depends on the ordering of the bases and the topology of the mesh. Therefore, the results for Test 1 and Test 2 are the same. In the results below the bases are in the order of the interior vertex, boundary vertices, boundary edges, interior edges and finally interior modes. The corresponding results for Test 3 are given in Figure 11 and Figure 11. Due to the locality of the frames, more sparsity appears as the mesh has more triangles.

5. Iterative methods and preconditioning

We recall the setting discussed in the introduction of this paper, where $V_{h}\subset H^{1}(\Omega)$ is a $C^{0}$ piecewise polynomial space. In the first part of this paper we have outlined how to construct representations of the spaces $V_{h}$ by frames which admit $L^{2}$ condition numbers which are bounded independently of the polynomial degree. The main purpose of this section is to present a brief discussion of the use of preconditioned iterative methods in the setting of representations by frames. In particular, we will clarify how the preconditioned iteration is effected by the boundedness of the $L^{2}$ condition number. However, before we discuss this in the setting of frames we will first review the more standard situation when the computation relies on a basis for the space $V_{h}$.

5.1. Representation by a basis

If $\Phi=\{\phi_{j}\}_{j=1}^{n}$ is a basis for the space $V_{h}$, then the two bijective maps $\tau_{h}:\mathbb{R}^{n}\mapsto V_{h}$ and $\mu_{h}=\tau_{h}^{\ast}:V_{h}^{\ast}\mapsto\mathbb{R}^{n}$, introduced in Section 1, are used to represent elements of $V_{h}$ and $V_{h}^{*}$, respectively. The stiffness matrix ${\mathbb{A}}_{h}$ admits the representation $\mathbb{A}_{h}:=\mu_{h}\mathcal{A}_{h}\tau_{h}$, where the operator $\mathcal{A}_{h}:V_{h}\to V_{h}^{*}$ is independent of the choice of basis. Discrete elliptic systems of the form (1.2), or equivalently (1.3), are most effectively solved by preconditioned iterative methods. Therefore, it seems appropriate to study the effect of the choice of basis on the complete preconditioned system.

In operator form the preconditioned system appears as

[TABLE]

where the preconditioner $\mathcal{B}_{h}$ is an operator, $\mathcal{B}_{h}:V_{h}^{*}\to V_{h}$, which is symmetric and positive definite with respect to the duality pairing between $V_{h}^{*}$ and $V_{h}$. Furthermore, its standard representation is the matrix ${{{\mathbb{B}}}_{h}}=\tau_{h}^{-1}{\mathcal{B}_{h}}\mu_{h}^{-1}=\tau_{h}^{-1}{\mathcal{B}_{h}}\tau_{h}^{-*}$. The two basic necessary conditions for the construction of an effective preconditioner are, i) the spectral condition number of ${\mathcal{B}_{h}}{\mathcal{A}_{h}}$ is well behaved, and ii) the matrix-vector products of the form ${{{\mathbb{B}}}_{h}}{{{\mathbb{A}}}_{h}}c$, for any $c\in\mathbb{R}^{n}$, can be evaluated fast. Since the stiffness matrix usually is represented by a sparse matrix, the second condition will hold if the matrix-vector products of the form ${{{\mathbb{B}}}_{h}}c$ can be evaluated fast. It is not the purpose of this paper to discuss the design of effective preconditioners. Instead we refer to [13, 20, 21] for general discussions of such constructions. Our main concern here is to discuss how the $L^{2}$ condition number of the basis influences the key properties of the preconditioned iterative method.

The convergence of a standard iterative method for the preconditioned system, such as the conjugate gradient method, will be governed by the spectral condition number of the coefficient matrix,

[TABLE]

However, this matrix is similar to the operator ${\mathcal{B}_{h}}{\mathcal{A}_{h}}$, and therefore its condition number is independent of the basis. On the other hand, the basis will effect the properties of the two matrices ${\mathbb{A}}_{h}$ and ${\mathbb{B}}_{h}$. These operators have to be evaluated in each iteration, and their conditioning will effect the numerical stability of the computations. Recall the inequality (2.2), which relates $\kappa({\mathbb{A}}_{h})$, $\kappa(\mathcal{I}_{h}\mathcal{A}_{h})$ and the $L^{2}$ condition number of the basis, $\kappa({\mathbb{M}}_{h})$. In fact, $\kappa({\mathbb{M}}_{h})$ is the only quantity on the right hand side of the inequality (2.2) which is basis dependent. Furthermore, the numerical experiments presented in Section 2 indicate that this inequality is rather sharp. Therefore, if the mass matrix is well-conditioned then $\kappa({\mathbb{A}}_{h})$ will behave approximately like the basis independent condition number $\kappa(\mathcal{I}_{h}\mathcal{A}_{h})$. This condition number reflects the properties of the elliptic operator we are approximating.

The situation for the preconditioner is similar. The matrix ${\mathbb{B}}_{h}$ admits the representation

[TABLE]

and from this we can easily derive the inequality

[TABLE]

Therefore, if the mass matrix is well conditioned, then we can conclude that $\kappa({\mathbb{B}}_{h})$ is essentially bounded by a basis independent quantity. We have therefore seen that, even if the choice of basis has no direct effect on the conditioning of the preconditioned system, an $L^{2}$ well-conditioned basis will result in matrix representations ${\mathbb{A}}_{h}$ and ${\mathbb{B}}_{h}$ with condition numbers that roughly behave like their basis independent counterparts. Below we will argue that these conclusions also hold if we allow representations by frames.

5.2. Representations by frames

Assume that we are given a frame $\Phi=\{\phi_{j}\}_{j=1}^{N}$ of $V_{h}$, where in general $N$ is larger than the dimension of $V_{h}$. The operators $\tau_{h}:\mathbb{R}^{N}\mapsto V_{h}$ and $\mu_{h}:V_{h}^{\ast}\mapsto\mathbb{R}^{N}$, are defined as above, i.e.,

[TABLE]

such that the identity $\mu_{h}f\cdot c=\langle f,\tau_{h}(c)\rangle$ holds. In this setting the operator $\tau_{h}$ is surjective, the operator $\mu_{h}:V_{h}^{\ast}\mapsto\mathbb{R}^{N}$ is injective, and as above $\tau_{h}$ and $\mu_{h}$ correspond to the dual of the other. If the $L^{2}$ condition number of the frame, $\mathcal{K}(\Phi)$, is controlled then we will argue that also in this case the matrix representations of the discrete elliptic operator $\mathcal{A}_{h}$ and suitable preconditioners $\mathcal{B}_{h}$ behaves roughly like the corresponding basis independent counterparts. In fact, this simply follows by proper block decompositions of the matrices.

The stiffness matrix, representing the coefficient operator $\mathcal{A}_{h}:V_{h}\mapsto V_{h}^{\ast}$ is still defined as $\mathbb{A}_{h}=\mu_{h}\mathcal{A}_{h}\tau_{h}$, cf. (1.4). While the operator $\mathcal{A}_{h}$ is positive definite, the stiffness matrix $\mathbb{A}_{h}$ is only positive semi-definite with $\ker(\mathbb{A}_{h})=\ker\left(\tau_{h}\right)$ and

[TABLE]

Here the orthogonal complement is with respect to the standard Euclidean inner product of $\mathbb{R}^{N}$. In fact, with respect to the orthogonal decomposition $\mathrm{Im}(\mu_{h})\oplus\ker\left(\tau_{h}\right)$ the matrix ${\mathbb{A}}_{h}$ has a block structure of the form

[TABLE]

where the matrix $\tilde{{\mathbb{A}}}_{h}$ is a positive definite matrix on $\mathrm{Im}(\mu_{h})$. The mass matrix, ${\mathbb{M}}_{h}$, with elements $\langle\phi_{i},\phi_{j}\rangle$, has a similar block structure of the form

[TABLE]

with respect to the decomposition $\mathrm{Im}(\mu_{h})\oplus\ker\left(\tau_{h}\right)$. Here the the matrix $\tilde{{\mathbb{M}}}_{h}$, mapping $\mathrm{Im}(\mu_{h})$ into itself, is positive definite. Furthermore, it follows from the observation done in Remark 1 that $\mathcal{K}(\Phi)=\kappa(\tilde{\mathbb{M}}_{h})$. The marices $\tilde{{\mathbb{A}}}_{h}$ and $\tilde{{\mathbb{M}}}_{h}$ can be related by the identity

[TABLE]

where $\tilde{\tau}_{h}$ denotes the restriction of $\tau_{h}$ to $\mathrm{Im}(\mu_{h})$. By arguing as above, this leads to

[TABLE]

which is a generalization of inequality (2.2). Therefore, we can again conclude that $\kappa(\tilde{A}_{h})$ behaves roughly like its basis independent counterpart, $\kappa(\mathcal{I}_{h}\mathcal{A}_{h})$, as long as the frame condition number $\mathcal{K}(\Phi)=\kappa(\tilde{\mathbb{M}}_{h})$ is well controlled.

In the setting of frames a preconditioner is represented by an $N\times N$ matrix $\mathbb{B}_{h}$ which is symmetric and positive definite with respect to the Euclidean inner product. Furthermore, the corresponding operator $\mathcal{B}_{h}:V_{h}^{*}\to V_{h}$ is given by the diagram

[TABLE]

i.e., $\mathcal{B}_{h}:=\tau_{h}\mathbb{B}_{h}\mu_{h}$. The operator $\mathcal{B}_{h}$ is symmetric in the sense that for $f,g\in V_{h}^{*}$,

[TABLE]

and $\mathcal{B}_{h}$ is positive definite since ${\mathbb{B}}_{h}$ has this property. Furthermore, the identity

[TABLE]

holds. We will implicitly assume that $\mathcal{B}_{h}$ has an interpretation as an operator from $V_{h}^{*}$ to $V_{h}$ which is independent of the choice of frame. In this respect, also the operator $\mathcal{B}_{h}\mathcal{A}_{h}$, and its spectral condition number, will be frame independent. In fact, we will justify this assumption in Section 5.3 below, in the case when the frame is derived from a basis of each of the local spaces ${Q_{f,r}^{\ast}}$, cf. (2.7).

Consider a linear system of the form (1.2), i.e., $\mathcal{A}_{h}u=f$, where the data $f\in V_{h}^{*}$ and $u\in V_{h}$. The preconditioned version of this system takes the form

[TABLE]

or in matrix-vector form

[TABLE]

where $c\in\mathbb{R}^{N}$ is any vector such that $\tau_{h}(c)=u$. Hence, even if the solution $u\in V_{h}$ is uniquely determined by the system, the vector $c$ is only determined up to addition of elements in $\ker\left(\tau_{h}\right)=\ker\left({\mathbb{B}}_{h}{\mathbb{A}}_{h}\right)$.

It is well known that Krylov subspace methods can be used to solve semidefinite systems, cf. for example [10]. In fact, if we consider a system of the form (5.7), with initial guess $c^{0}=0$, then the initial residual

[TABLE]

where the superscript $\angle$ indicates orthogonality with respect to the inner product generated by the positive definite matrix $\mathbb{B}_{h}^{-1}$. The matrix ${\mathbb{B}}_{h}{\mathbb{A}}_{h}$ maps $\ker(\tau_{h})^{\angle}$ to itself, and as consequence all the vectors in the associated Krylov spaces, spanned by vectors of the form $({\mathbb{B}}_{h}{\mathbb{A}}_{h})^{k}r^{0}$, will be in this space. As a consequence, we can view the iterative method as if it is restricted to the space $\ker(\tau_{h})^{\angle}$, where the system (5.7) has a unique solution. Furthermore, the convergence rate will be bounded by the spectral condition number of the coefficient matrix ${\mathbb{B}}_{h}{\mathbb{A}}_{h}$ restricted to this space. If we let $\hat{\tau}_{h}:\ker(\tau_{h})^{\angle}\mapsto V_{h}$ be the restriction of $\tau_{h}$ to $\ker(\tau_{h})^{\angle}$ then $\hat{\tau}$ is invertible and from (5.5) we obtain

[TABLE]

as a generalization of the similarity relation (5.1). In particular, this shows that the spectral condition number of the matrix ${\mathbb{B}}_{h}{\mathbb{A}}_{h}$, restricted to $\ker(\tau_{h})^{\angle}$, is equal to the spectral condition number of the operator $\mathcal{B}_{h}\mathcal{A}_{h}:V_{h}\to V_{h}$. Therefore, also in the case of representations by frames, we can conclude that the performance of a preconditioned Krylov space method is, in a proper sense, independent of the choice of representation of the spaces $V_{h}$ and $V_{h}^{*}$. However, as we have seen above, a well condition frame guarantees that the conditioning of the stiffness matrix reflects the condition number of its basis independent counterpart, cf. (5.4). A similar conclusion for the preconditioner, i.e., a generalization of inequality (5.3), is also easily established.

5.3. Representation of the preconditioner

The discussion above is based on the assumption that the matrix ${\mathbb{B}}_{h}$, representing the preconditioner $\mathcal{B}_{h}:V_{h}^{*}\to V_{h}$, is positive definite, and that the operator $\mathcal{B}_{h}$ has an interpretation which is independent of the representations of the spaces $V_{h}$ and $V_{h}^{*}$. Here we will argue that if $V_{h}$ is of the form ${\mathcal{P}_{r}}(\mathcal{T}_{h})$ and the frame $\Phi$ is generated by an overlapping decomposition of the form (2.7), then this assumption is indeed very natural. To illustrate this we consider an additive Schwartz preconditioner of the form proposed in [15]. In operator form this preconditioner has the structure

[TABLE]

where $\mathcal{B}_{h}^{0}$ is a “coarse space preconditioner”, and $\mathcal{B}_{f,h}$ are local preconditioners defined with respect to local spaces ${Q_{f,r}^{\ast}}\subset V_{h}$. More precisely, each of the operators $\mathcal{B}_{f,h}$ are preconditioner for the corresponding local operator $\mathcal{A}_{f,h}$ defined by the bilinear form $a(\cdot,\cdot)$ with respect to the local space ${Q_{f,r}^{\ast}}$, while $\mathcal{B}_{h}^{0}$ is a corresponding preconditioner defined with respect to the piecewise linear space ${\mathcal{P}_{1}}(\mathcal{T}_{h})\subset V_{h}$.

Consider the set up in Section 3.2 above, where for each $f\in\Delta=\Delta(\mathcal{T}_{h})$ the set $\Phi_{f}=\{\phi_{f,k}\}_{k=1}^{N_{f}}$ is a basis for the local space ${Q_{f,r}^{\ast}}$, and $\Phi=\{\Phi_{f}\}_{f\in\Delta}$. The natural representations of the operators $\mathcal{B}_{f,h}$ are $N_{f}\times N_{f}$ matrices of the form ${\mathbb{B}}_{f,h}=\tau_{f,h}^{-1}\mathcal{B}_{f,h}\mu_{f,h}^{-1}$, where the representations $\tau_{f,h}$ and $\mu_{f,h}$ are defined as above, but now with respect to the local spaces ${Q_{f,r}^{\ast}}$. Furthermore, since $\Phi_{f}$ is a basis for this space, each of the maps $\tau_{f,h}$ and $\mu_{f,h}$ is invertible, and as a consequence, the matrix ${\mathbb{B}}_{f,h}$ is positive definite. The matrix ${\mathbb{B}}_{h}$, representing the preconditioner $\mathcal{B}_{h}$, will now be of the form

[TABLE]

where $N=\sum_{f}N_{f}$. Here ${\mathbb{B}}_{h}^{0}$ is a symmetric and positive semidefinite matrix representing the operator $\mathcal{B}_{h}^{0}$, while the block diagonal matrix, ${\operatorname{diag}}\{{\mathbb{B}}_{f,h}\}$, is symmetric and positive definite, since each block has this property. We refer to [15] for more details and tests of numerical performance.

Remark 3**.**

The additional global operator $\mathcal{B}_{h}^{0}$ appears in the preconditioner proposed in [15] in order to obtain condition numbers which are independent of the mesh size $h$. In the present paper, where we consider the mesh $\mathcal{T}_{h}$ to be fixed, we could have dropped this term and still obtain condition numbers which are bounded uniformly with respect to the polynomial degree $r$.

6. Concluding remarks

The purpose of this paper is to discuss how to represent $H^{1}$ finite element spaces of high degree. More precisely, we study the spaces $\mathcal{P}_{r}(\mathcal{T}_{h})$, consisting of $C^{0}$ piecewise polynomial spaces of degree $r$ with respect to a simplicial mesh $\mathcal{T}_{h}$. When the degree $r$ grows, the choice of basis for the spaces will intuitively affect the properties of the corresponding linear systems derived from a finite element discretization. The construction outlined in this paper, based on properties of the bubble transform and of Jacobi polynomials with respect to simplexes, leads to frames with $L^{2}$ condition numbers which are independent of the polynomial degree. In this respect, we have been able to present a procedure to construct well-conditioned frames. Furthermore, we have shown that such frames leads to stiffness matrices, and matrix representations of corresponding preconditioners, with condition numbers that behave roughly like their operator counterparts.

A possible disadvantage of a representation by a frame instead of a basis, is that the number of unknowns is increasing. For the frames proposed above the following table compares the dimension of the frame with a corresponding standard basis.

Basis Frame

1D $V+(r-1)E=rV-(r-1)$ $V+rV+(r-1)E=(2r-1)V-(r-1)$

2D $V+(r-1)E+\left(\frac{1}{2}(r+2)(r+1)-3\right)F$ $(r+1)V+\frac{1}{2}r(r-1)E+\left(\frac{1}{2}(r+2)(r+1)-3\right)F$

3D $V+(r-1)E+\left(\frac{1}{2}(r+2)(r+1)-3\right)F$ $(r+1)V+\frac{1}{2}r(r-1)E+\frac{1}{6}r(r-1)(r-2)F$

$\quad+\left(\frac{1}{6}(r+3)(r+2)(r+1)-4\right)T$ $\quad+\left(\frac{1}{6}(r+3)(r+2)(r+1)-4\right)T$

Here $V,E,F$ and $T$ are the number of vertices, edges, faces and 3D cells in the triangulation. We observe that even if the frame representations have redundancy, the asymptotic order of the total dimensions remains the same. Of course, in addition to condition numbers there are potentially a number of other criteria that could have been used to choose a basis or a frame, for example sparsity and fast evaluation of the stiffness matrix and the preconditioner. However, such discussions are outside the scope of the present paper.

Appendix. Jacobi polynomials on simplexes

The purpose of this appendix is to present precise formulas for basis functions of the local spaces $Q_{f,r}^{*}$. We will first consider the case when $\dim f<d$. We recall from Section 4 above that if we can construct bases $\{p_{\bm{s}}\}$ for the polynomial spaces $\mathcal{P}_{r-m-1}(S_{m}^{c})$, which are orthonormal with respect to the weighted $L^{2}$ inner product with weight $w_{m}(\lambda)=(\Pi\lambda)_{m}^{2}b(\lambda)^{d-m-1}$, then the corresponding set $\{\lambda_{f}^{*}[(\Pi\lambda)_{m}p_{\bm{s}}]\}$ will be uniformly scaled bases for the spaces $Q^{\ast}_{f,r}$. Therefore, we will briefly outline how such bases for the spaces $\mathcal{P}_{r-m-1}(S_{m}^{c})$ can be constructed. This discussion is based on multi-variate Jacobi polynomials on simplexes, and is mostly taken from [4, Section 5.3], cf. also [18] where the scaling of the polynomials is corrected.

We use $J_{s}^{\alpha,\beta}(x)$ to denote the orthonormal Jacobi polynomial on interval $[-1,1]$ with weight $w^{\alpha,\beta}:=(1-x)^{\alpha}(1+x)^{\beta}$, i.e.

[TABLE]

Using a linear transform $\xi=(x+1)/2$, we get Jacobi polynomials on the reference interval $[0,1]$:

[TABLE]

where

[TABLE]

We introduce the notation

[TABLE]

so that $\{\tilde{J}_{s}^{\alpha,\beta}(\xi)\},s=0,1,\cdots$ is a single variate orthogonal basis on $[0,1]$ with weight $w^{\alpha,\beta}$.

We define

[TABLE]

Given integer $s\geq 0$ and $0\leq m<d$, we will define a basis for the space $\mathcal{P}_{\bm{s}}\left(S_{m}^{c}\right)$. We consider the polynomials $J_{\bm{s}}(\lambda)=J_{\bm{s}}(\lambda_{0},\lambda_{1},\cdots,\lambda_{m})$ given by

[TABLE]

where $\bm{s}=(s_{0},\dots,s_{m})$ is a multi-index. Here the constants $a_{j}$ and $c_{\bm{s}}$ are given by $a_{j}=2\sum_{i=j+1}^{m}s_{i}+d+2m-3j-1$ and

[TABLE]

The polynomials (6.2) form a mutually orthonormal basis for $\mathcal{P}_{s}(S_{m}^{c})$ with the weight $w_{m}$.

For $f\in\Delta_{d}$, the construction is similar and more standard. To construct a basis for $Q_{f,r}^{\ast}=\mathaccent 23{\mathcal{P}}_{r}(f)$, we can follow [4] to construct $L^{2}$ orthonormal bases $\left\{q_{\bm{s}}\right\}$ for $\mathcal{P}_{r-d-1}(S_{d})$ with weight $(\Pi\lambda)_{d}^{2}$. The set $\{\lambda_{f}^{*}[(\Pi\lambda)_{d}q_{\bm{s}}]\}$ will be uniformly scaled bases for the spaces $Q^{\ast}_{f,r}$, where $\lambda_{f}:f\mapsto S_{d}$ is defined by the barycentric coordinates.

Acknowledgments.

The authors are grateful to Douglas N. Arnold, Richard S. Falk and Jinchao Xu for several discussions about the results of this paper.

The research leading to these results has received funding from the European Research Council under the European Union’s Seventh Framework Programme (FP7/2007-2013) / ERC grant agreement 339643. The research of the first author leading to the results of this paper was partly carried out during his affiliation with the University of Oslo, partly supported by China Scholarship Council (CSC), project 201506010013.