TY - JOUR
AU - Jüttler, Bert
AB - Abstract In this article, we prove the inf–sup stability of an adaptive isogeometric discretization of the Stokes problem. The discretization is based on the hierarchical generalization of the isogeometric Taylor–Hood and Sub-Grid elements, which were described by Bressan & Sangalli (2013, Isogeometric discretizations of the Stokes problem: stability analysis by the macroelement technique. IMA J. Numer. Anal., 33, 629–651) for tensor-product splines. In order to extend the existing proof to the hierarchical setting, we need to adapt some of the steps considerably. In particular, the required local approximation estimate is obtained by analysing the properties of the quasi-interpolant of Speleers & Manni (2016, Effortless quasi-interpolation in hierarchical spaces. Numer. Math., 132, 155–184) with respect to certain Sobolev norms. In addition to the theoretical results, we also perform numerical tests in order to analyse the dependency of the inf–sup constant on the mesh regularity assumptions. Finally, the article also presents a numerical convergence test of the resulting adaptive method on a T-shaped domain. 1. Introduction The Taylor–Hood (TH) element is a well known mixed discretization, both of incompressible elasticity and incompressible fluid dynamic. The generalization of the TH method to isogeometric analysis (IGA) was initially studied by Bazilevs et al. (2006) for $$C^0$$ discretizations (i.e., isoparametric FEM) and later extended to smoother function spaces: $$C^1$$ by Bressan (2011) and $$C^k$$ discretizations, $$k\geq0$$ by Bressan & Sangalli (2013). The latter contains general conditions for stability and two families of inf–sup-stable pairs: the isogeometric TH and the isogeometric Sub-Grid (SG) element. In the SG method, stability is obtained by refining the velocity grid, and the same idea was independently proposed by Rüberg & Cirak (2012). The cited articles deal only with tensor-product spline spaces. In this article, we recall the theory and extend the stability result to adaptive methods based on hierarchical splines. To our knowledge, this is the first extension of TH and SG methods to nontensor-product splines. Note that there are also div-conforming elements for the Stokes problem, which we do not consider in this article. Div-conforming methods with tensor-product splines were studied by Buffa et al. (2011) and Evans & Hughes (2013a,b). Then they were extended to T-splines (Buffa et al., 2014), LR-splines (Johannessen et al., 2015) and hierarchical splines (Evans et al., 2017). In order to make the article self-contained, we recall the structure of the proof for tensor-product splines in Section 2. In Section 3, we quickly introduce hierarchical splines, we prove the properties required by the results recalled in Section 2 and conclude that a discretization based on hierarchical splines is stable under mild assumptions on the mesh structure (Theorem 3.1). In Section 4, we report numerically computed inf–sup constants for meshes that conform and do not conform to the hypothesis of the stability result. In Section 5, we report the results of the application of a simple adaptive method based on hierarchical splines to a problem having corner singularities. 2. Preliminaries This section is composed of three independent parts: the description of the inf–sup theory in the context of the Stokes problem, a quick description of isogeometric analysis, and a summary of the techniques used for assessing stability of the TH element in FEM and the TH and SG method in IGA. 2.1 Stokes problem and inf–sup stability The Stokes system of equations is a linearized model of the stationary flow of a nearly incompressible fluid under the assumption of negligible inertia. The unknowns are the velocity of the fluid $${\mathbf{u}}$$ and the pressure $$p$$ in the domain $${\it{\Omega}}\subset {\mathbb{R}}^n$$. The data are the fluid viscosity $$\nu$$, the volume forces $${\mathbf{f}}$$, the boundary forces $${\mathbf{u}}_N$$ on $${\it{\Gamma}}_N\subseteq{\partial}{{\it{\Omega}}}$$ and the velocity values $${\mathbf{u}}_D$$ on $${\it{\Gamma}}_D={\partial}{{\it{\Omega}}}\setminus {\it{\Gamma}}_N$$. We also include a source term $$g$$ that is usually absent in physical models. The Stokes system is   \begin{equation} \left \{\!\begin{aligned} -\nu{{\it{\Delta}}}{{\mathbf{u}}}+ {\mathrm{grad}\,{p}} &= f\\ -{\text{div}}{{\mathbf{u}}}&=g\\ {\mathbf{u}}&={\mathbf{u}}_D \text{ on }{\it{\Gamma}}_D\\ {\mathrm{grad}\,{u}}\cdot {\mathbf{n}}-{\mathbf{n}} p&={\mathbf{u}}_N \text{ on }{\it{\Gamma}}_N. \end{aligned}\right . \end{equation} (2.1) Note that if $${\it{\Gamma}}_N=\emptyset$$, then $$p$$ is determined solely up to a constant because only its gradient $${\mathrm{grad}\,{p}}$$ appears in (2.1). Similarly if $${\it{\Gamma}}_D=\emptyset$$, then $${\mathbf{u}}$$ is determined only up to a constant vector. These degenerate cases can be handled by constraining the mean value of $$p$$ or $${\mathbf{u}}$$, respectively. In order to simplify the notation, we fix $$\nu=1$$ and write the formulas for the case $${\mathbf{u}}_D=0$$ and $${\mathbf{u}}_N=0$$. The general case can be reduced to this by the following strategies: If $${\mathbf{u}}_D \neq 0$$, then it is possible to choose $${\mathbf{z}}:{\it{\Omega}}\to{\mathbb{R}}^n$$ such that $${\mathbf{z}}|_{{\it{\Gamma}}_D}={\mathbf{u}}_D$$ and solve for $$\tilde{{\mathbf{u}}}={\mathbf{u}}-{\mathbf{z}}$$ with data $$\tilde{{\mathbf{f}}}= {\mathbf{f}}-{{\it{\Delta}}}{{\mathbf{z}}}$$ and $$\tilde{g}= g-{\text{div}}{{\mathbf{z}}}$$. If $${\mathbf{u}}_N \neq 0$$, then the boundary term can be included in the force term in the weak formulation of the problem. The Stokes problem can be written as a mixed formulation by setting   \begin{equation*} {a}({\mathbf{w}},{\mathbf{u}})=\int_{\it{\Omega}} {\mathrm{grad}\,{{{\mathbf{w}}}}}:{\mathrm{grad}\,{{{\mathbf{u}}}}}, \qquad {b}({\mathbf{w}},q)=-\int_{\it{\Omega}} {\text{div}}{{\mathbf{w}}}\ q. \end{equation*} Then (2.1) becomes: find $$ ({\mathbf{u}},p) \in V \times P $$ such that $$\forall ({\mathbf{w}},q)\in V \times P $$  \begin{equation} \left \{\!\begin{aligned} {a}({\mathbf{u}}, {\mathbf{w}}) + {b}({\mathbf{w}},p) &= f({\mathbf{w}})\\ {b}({\mathbf{u}},q) &= g(q), \end{aligned}\right . \end{equation} (2.2) where $$V$$ and $$P$$ are function spaces, $${a}:V \times V\to{\mathbb{R}}$$ and $${b}:V \times P\to{\mathbb{R}}$$ are continuous bilinear forms and $$f\in V^*$$, $$g\in P^*$$ represent both the source and the boundary data. The kernel of $${B}:V\to P^*$$ defined by $${B}( {\mathbf{w}})(q)={b}({\mathbf{w}},q)$$ plays an important role in the analysis of problem (2.2) as described in Theorem 1.1 of Brezzi & Fortin (1991). A specialization for symmetric $${a}$$ and $$\ker {B}=\{0\}$$ is reported below: Theorem 2.1 There exists a unique solution $$({\mathbf{u}},p)$$ of (2.2) if (i) $$\exists {\alpha} > 0:\, \forall {\mathbf{w}} \in \ker {B},\, {a}({\mathbf{w}}, {\mathbf{w}}) \geq {\alpha} {\Vert{{\mathbf{w}}}\Vert}_{V}^2$$ (ii) $$\exists {\beta} > 0:\, \forall q \in P,\, \exists {\mathbf{w}} \in V:\, {b}({\mathbf{w}},q) \geq {\beta} {\Vert{{\mathbf{w}}}\Vert}_{V} {\Vert{q}\Vert}_{P}$$. Moreover the solution satisfies:   \begin{align} {\Vert{ {\mathbf{u}}}\Vert}_{V} &\leq {\alpha}^{-1}{\Vert{f}\Vert}_{V^*}+{\beta}^{-1}(1+{\alpha}^{-1} {\Vert{{a}}\Vert}){\Vert{g}\Vert}_{P^*}\\ \end{align} (2.3)  \begin{align} {\Vert{p}\Vert}_{P} &\leq {\beta}^{-1}(1+{\alpha}^{-1} {\Vert{{a}}\Vert}){\Vert{f}\Vert}_{V^*}+{\beta}^{-2} {\Vert{{a}}\Vert}(1+{\alpha}^{-1} {\Vert{{a}}\Vert}){\Vert{g}\Vert}_{P^*}. \end{align} (2.4) Condition (2.1) is called coercivity on the kernel and condition (2.1) is called inf–sup condition. A conformal discretization of (2.2) is a restriction to a finite dimensional subspace: find $$ ({\mathbf{u}}_h,p_h) \in V_h \times P_h \subseteq V \times P $$ such that $$\forall ({\mathbf{w}},q)\in V_h \times P_h $$  \begin{equation} \left\{\!\begin{aligned} {a}({\mathbf{u}}_h, {\mathbf{w}}) + {b}({\mathbf{w}},p_h) &= f({\mathbf{w}})\\ {b}({\mathbf{u}}_h,q) &= g(q). \end{aligned} \right . \end{equation} (2.5) Theorem 2.1 applies to (2.5), but the validity of conditions (2.1) and (2.1) for the pair $$V,P$$ does not imply their validity for the pair $$V_h,P_h$$. The main result concerning the discretization error is Theorem 2.1 in Brezzi & Fortin (1991), which in our case can be written as: Theorem 2.2 The difference between the solutions $$({\mathbf{u}},p)$$ of (2.2) and $$({\mathbf{u}}_h,p_h)$$ of (2.5) is bounded by   \begin{align} {\Vert{{\mathbf{u}}-{\mathbf{u}}_h}\Vert}_V &\leq \left( 1+\frac{{\Vert{{a}}\Vert}}{{\alpha}_h}\right)\left( 1+\frac{{{\Vert{{b}}\Vert}}}{{\beta}_h}\right) \inf_{{\mathbf{w}}\in V_h}{\Vert{{\mathbf{u}}-{\mathbf{w}}}\Vert}_V +\frac{{{\Vert{{b}}\Vert}}}{{\alpha}_h}\inf_{q\in P_h}{\Vert{p-q}\Vert}_P;\\ \end{align} (2.6)  \begin{align} {\Vert{p-p_h}\Vert}_{P} &\leq \left( 1+\frac{{{\Vert{{b}}\Vert}}}{{\beta}_h}\right) \inf_{q\in P_h}{\Vert{p-q}\Vert}_P+\frac{{\Vert{{a}}\Vert}}{{\beta}_h}{\Vert{{\mathbf{u}}-{\mathbf{u}}_h}\Vert}_V \end{align} (2.7) provided that conditions (2.1) and (2.1) hold with constants $${\alpha}_h, {\beta}_h$$ for the pair $$V_h,P_h$$. Appropriate spaces for the Stokes problem are   \begin{equation} V=\begin{cases} \{{\mathbf{w}} \in H^1({\it{\Omega}})^n: {\mathbf{w}} |_{{\it{\Gamma}}_D}=0 \} & {\it{\Gamma}}_D \text{ has positive}\,n-1\text{-dimensional volume}\\ \{{\mathbf{w}} \in H^1({\it{\Omega}})^n: \int_{\it{\Omega}}{\mathbf{w}} =0 \}& \text{otherwise} \end{cases} \end{equation} (2.8) with norm $${\Vert{{\mathbf{w}}}\Vert}_V={|{{\mathbf{w}}}|}_{H^1({\it{\Omega}})^n}={a}({\mathbf{w}},{\mathbf{w}})^{\frac{1}{2}}$$ and   \begin{equation} P=\begin{cases}L^2({\it{\Omega}}) & {\it{\Gamma}}_N \text{ has positive}\,n-1\text{-dimensional volume}\\ \{q\in L^2({\it{\Omega}}):\int_{\it{\Omega}} q=0\} & \text{otherwise} \end{cases} \end{equation} (2.9) with norm $${\Vert{q}\Vert}_P= {\Vert{q}\Vert}_{L^2}$$. Then $${\alpha}={\alpha}_h={\Vert{{a}}\Vert}=1$$, $${{\Vert{{b}}\Vert}}= 1$$ and (2.1) is fulfilled for all pairs $$V_h,P_h$$. Condition (2.1) for $$V,P$$ depends on the shape of the domain. A positive $${\beta}$$ as in (2.1) exists if $${\it{\Omega}}$$ is a John domain (Acosta et al., 2006; Costabel et al., 2015) and thus, in particular, if $${\it{\Omega}}$$ is a Lipschitz domain. The discrete inf–sup condition depends on the choice of $$V_h,P_h$$. When it holds, Theorem 2.2 provides the following error estimates:   \begin{align} {\Vert{{\mathbf{u}}-{\mathbf{u}}_h}\Vert}_V &\leq 2\frac{1+{\beta}_h}{{\beta}_h} \inf_{{\mathbf{w}}\in V_h}{\Vert{{\mathbf{u}}-{\mathbf{w}}}\Vert}_V +\inf_{q\in P_h}{\Vert{p-q}\Vert}_P;\\ \end{align} (2.10)  \begin{align} {\Vert{p-p_h}\Vert}_P &\leq 2 \frac{1+{\beta}_h}{{\beta}_h^2}\inf_{{\mathbf{w}}\in V_h}{\Vert{{\mathbf{u}}-{\mathbf{w}}}\Vert}_V + \frac{2+{\beta}_h}{{\beta}_h} \inf_{q\in P_h}{\Vert{p-q}\Vert}_P. \end{align} (2.11) The above inequalities show the competing behaviour of the approximation properties of $$P_h$$ and the inf–sup condition. On one hand, choosing a larger space $$P_h$$ decreases the approximation error $$\inf_{q\in P_h}{\Vert{p-q}\Vert}_P$$; on the other hand, it can simultaneously decrease the inf–sup constant $${\beta}_h$$. Vice-versa choosing a smaller space $$P_h$$ provides weaker approximation properties and can increase the inf–sup constant. These inequalities also show that the existence of a common lower bound $${\beta}_0\leq {\beta}_h$$ for all the pairs $$V_h,P_h$$ defined by a method implies the quasi-optimality of the error for the method. 2.2 Isogeometric analysis Isogeometric analysis (IGA) is a framework that streamlines the design process of technological products. It was introduced by Hughes et al. (2005), it is described in the monograph of Cottrell et al. (2009) and analysed mathematically in Beirão da Veiga et al. (2014). It aims at integrating computer aided design (CAD) and computer aided engineering (CAE) by sharing the same geometry representation techniques. In the traditional approach, the geometry provided by CAD must be converted to a mesh on which the simulation can be performed. Then, in case of geometry optimization, the result is re-imported into the CAD system. Experience suggests that the conversions of the geometry representation is responsible for a significant percentage of the product development time. IGA tries to optimize the complexity of the whole process by reducing the conversion cost. IGA uses common industrial technology: in particular, NURBS functions and parametrizations. The physical domain $${\it{\Omega}}$$ is assumed to be the image of an NURBS map $$G:\widehat{{\it{\Omega}}}=[0,1]^n\to {\mathbb{R}}^k$$. This is the most common method to represent a domain in CAD. Then the isoparametric approach is applied: the discretization is constructed on the parametric domain $$\widehat{{\it{\Omega}}}$$ and pushed to the physical domain $${\it{\Omega}}$$. Following the IGA convention, quantities related to the parametric domain are marked with a hat: $$\hat\ $$. The discretization spaces in IGA are the push-forward of spline (or NURBS) spaces defined in the parametric domain:   \begin{align} V_h &=\left\{{\mathbf{w}}=\hat{{\mathbf{w}}}{\circ} G^{-1}:\, \hat{{\mathbf{w}}}\in \widehat V_h \right\}\cap V\\ \end{align} (2.12)  \begin{align} P_h &=\left\{q=\hat{q} {\circ} G^{-1}:\, \hat{q}\in \widehat P_h\right \} \cap P. \end{align} (2.13) The approximation error is studied in the parametric domain and then the result is transferred on the physical domain (Beirão da Veiga et al., 2014). This requires the following norm equivalences: for each $$\phi:{\it{\Omega}}\to{\mathbb{R}}$$  \begin{align} {\Vert{\hat{\phi}}\Vert}_{H^1}&\cong{\Vert{\hat{\phi}{\circ} G^{-1}}\Vert}_{H^1}={\Vert{\phi}\Vert}_{H^1}\\ \end{align} (2.14)  \begin{align} {\Vert{\hat{\phi}}\Vert}_{L^2}&\cong{\Vert{\hat{\phi} {\circ} G^{-1}}\Vert}_{L^2}={\Vert{\phi}\Vert}_{L^2}. \end{align} (2.15) The equivalences above depend on the regularity of $$G$$. In particular, they are implied by the following assumptions: $$G:\widehat{{\it{\Omega}}}\to{\it{\Omega}}$$ is a homeomorphism: a continuous bijective map with continuous inverse. Both $${\mathrm{grad}\,{G}}$$ and $${\mathrm{grad}\,{G}}^{-1}$$ are continuous and bounded on their respective domains. The discrete spaces have the same (or lower) local smoothness than the parametrization $$G$$. This means that the spline spaces must allow for discontinuities of the $$k$$ derivatives at $${\mathbf{x}}$$ if $$G$$ is not $$C^k$$ at $${\mathbf{x}}$$. The second condition can be weakened to piecewise continuity of $${\mathrm{grad}\,{G}}$$ and $${\mathrm{grad}\,{G}}^{-1}$$, but then the inf–sup-stability must be proved on each domain on which these are continuous. This allows to consider multipatch domains. 2.3 Discrete inf–sup The inf–sup stability of the TH discretizations in FEM was studied by Engelman & Bercovier (1979), Verfürth (1984) and Stenberg (1990). Stenberg combined the previous results and proved the stability of the generalized TH elements. The proof is based on the inf–sup stability of the pair $$V,P$$, an approximation operator $$\wp:V\to V_h$$, the trick from Verfürth (1984) and the macro-element technique. Lemma 2.1 (Verfürth’s trick) The inf–sup condition (2.1) for the pair $$V_h,P_h$$ is implied by the following two auxiliary conditions:   \begin{align} \forall q \in P_h,\, \exists {\mathbf{w}}_1 \in V_h:\ {b}({\mathbf{w}}_1,q) &\geq \big(c_1-\eta(q) \big) {\Vert{{\mathbf{w}}_1}\Vert}_{V} {\Vert{q}\Vert}_{P}\\ \end{align} (2.16)  \begin{align} \forall q \in P_h,\, \exists {\mathbf{w}}_2 \in V_h:\ {b}({\mathbf{w}}_2,q) &\geq c_2 \eta(q){\Vert{{\mathbf{w}}_2}\Vert}_{V} {\Vert{q}\Vert}_{P}, \end{align} (2.17) where $$c_1,c_2>0$$ and $$\eta(q)\geq 0$$. Proof. By choosing either $${\mathbf{w}}={\mathbf{w}}_1$$ or $${\mathbf{w}}={\mathbf{w}}_2$$ it follows that (2.1) holds for the pair $$V_h,P_h$$ with   \begin{equation} {\beta}_h\geq \min_{q\in P_h}\max\{c_1 - \eta(q), c_2 \eta(q)\}\geq \min_{t \in {\mathbb{R}}_+}\max\{c_1 - t, c_2 t\} = \frac{c_1c_2}{1+c_2}. \end{equation} (2.18) □ Lemma 2.2 (approximation) The condition in (2.16) holds with   \begin{equation} c_1=\frac{{\beta}}{{\Vert{\wp}\Vert}}\quad \text{and}\quad \eta(q)=\frac{{\Vert{\rho {\mathrm{grad}\,{q}}}\Vert}_{L^2({\it{\Omega}})}}{{\Vert{q}\Vert}_{L^2({\it{\Omega}})}} \end{equation} (2.19) provided that there exists $$\wp:V\to V_h$$ such that   \begin{equation} {\Vert{\rho^{-1}( {\mathbf{w}}- \wp {\mathbf{w}})}\Vert}_{L^2({\it{\Omega}})} \leq {|{{\mathbf{w}}}|}_{H^1({\it{\Omega}})}, \end{equation} (2.20) where $$\rho:{\it{\Omega}}\to{\mathbb{R}}$$ is a multiple of the mesh size: on each element $$E$$, $$\rho(E)=c R(E)$$ for some $$c>0$$. Lemma 2.3 (macro-elements). The condition in (2.17) holds for $$\eta(q)$$ as in (2.19) with constant   \begin{equation} c_2\geq \frac{c_{\mathrm{SD}}}{m} \end{equation} (2.21) if the domain can be covered by subdomains $$\{M_1,\dots,M_N\}$$ such that at most $$m$$ subdomains overlap (each point $$p\in{\it{\Omega}}$$ is contained in at most $$m$$ sets $${\text{int}}{M_i}$$); on each $$M_i$$ the following local inf–sup condition holds:   \begin{equation} \forall q\in P_h,\, \exists {\mathbf{w}}_i \in V_h\cap H^1_0(M_i)^n:\ {b}({\mathbf{w}}_i,q) \geq c_{\mathrm{SD}} {\Vert{\rho {\mathrm{grad}\,{q}}}\Vert}_{L^2(M_i)}{|{{\mathbf{w}}_i}|}_{V}. \end{equation} (2.22) Theorem 2.3 If there exists a continuous approximation operator $$\wp:V\to V_h$$ as in Lemma 2.1 and a subdomain covering of $${\it{\Omega}}$$ as in Lemma 2.3, then the inf–sup condition holds with   \begin{equation} {\beta}_h\geq \frac{{\beta} c_{\mathrm{SD}}}{(m+c_{\mathrm{SD}}){\Vert{\wp}\Vert}}>0. \end{equation} (2.23) Proof. Combine Lemmas 2.1, 2.2 and 2.3. □ The strategy to use Theorem 2.3 effectively is to consider subdomains that are images through a similarity of macro-elements in some catalogue $$\mathfrak{C}$$ of stable macro-elements. Each macro-element in $$\mathfrak{C}$$ is a triplet $$(M_\sigma,V_\sigma,P_\sigma)$$, where $$M_\sigma\subset {\mathbb{R}}^n$$, $$V_\sigma\subset H^1_0(M_\sigma)^n $$ is a velocity space and $$P_\sigma\subset L^2(M_\sigma)$$ is a pressure space. A macro-element is stable if $$C_\sigma>0$$ exists such that:   \[ \forall q\in P_\sigma,\,\exists {\mathbf{w}}\in V_\sigma:\ {b}({\mathbf{w}},q) \geq C_{\sigma} {\Vert{\rho {\mathrm{grad}\,{q}}}\Vert}_{L^2({\mathbb{R}}^n)}{|{{\mathbf{w}}}|}_{H^1({\mathbb{R}}^n)^n}. \] The macro-elements in the catalogue $$\mathfrak{C}$$ are indexed by a shape parameter $$\sigma\in{\mathbb{R}}^s$$. The constant $$C_\sigma$$ is preserved under similarities and depends continuously on $$\sigma$$. If the shape parameter $$\sigma$$ varies in a compact set $$\mathfrak{S}$$, then Theorem 2.3 applies and condition (2.1) holds with   \begin{equation} {\beta}_h\geq \frac{{\beta} C_{\mathfrak{S}}}{(m+C_{\mathfrak{S}}){\Vert{ \wp }}\Vert}, \end{equation} (2.24) where $$\displaystyle C_\mathfrak{S}=\min_{\sigma\in \mathfrak{S}} \{C_{\sigma}\}$$. The macro-element collection approach does not apply directly to IGA because of the distortion caused by the parametrization of the domain. The solution is to use a macro-element collection $${\mathfrak{C}}$$ on the parametric domain $$\widehat{{\it{\Omega}}}$$ and to transfer the result to the physical macro-elements. The effect of the mapping on the inf–sup was studied by Bressan (2011) by decomposing $$G$$ in an affine part $$F$$ and a reminder term $$T$$:   \begin{equation} G(\hat{{\mathbf{x}}}) = F(\hat{{\mathbf{x}}})+T(\hat{{\mathbf{x}}}). \end{equation} (2.25) $$F(\hat{{\mathbf{x}}})$$ is given by the averages of $$G$$ and $${\mathrm{grad}\,{G}}$$ on the parametric macro-element $$\widehat M$$. It preserves the inf–sup up to a factor that depends on the regularity of $$G$$. The reminder $$T$$ can be bounded in terms of the size of the parametric macro-element $$\widehat M$$. The compactness of $$\mathfrak{S}$$ is obtained by requiring shape-regularity of the elements of the parametric macro-elements (and thus of the mesh). The shape regularity of a set $$S$$ is the ratio between its radius $$R(S)$$ and that of the biggest in-circle $$r(S)$$. A mesh $${\mathscr{M}}$$ is shape regular if there is $$C_{\mathrm{shape}}$$ such that for all elements $$E \in {\mathscr{M}}$$  \begin{equation} \frac{R(E)}{r(E)}\leq C_{\mathrm{shape}}. \end{equation} (2.26) A consequence of shape-regularity is local quasi-uniformity of the macro-elements. A mesh $${\mathscr{M}}$$ is locally quasi-uniform if there is $$C_{\mathrm{uniform}}$$ such that for each pair of neighbouring elements $$E_1, E_2\in {\mathscr{M}}$$  \begin{equation} \frac{R(E_1)}{R(E_2)}\leq C_{\mathrm{uniform}}. \end{equation} (2.27) As seen later, local quasi-uniformity is required for the construction of the approximation operator presented here. Lemma 2.4 (Bressan, 2011) The local inf–sup condition of Lemma 2.3 for the pair $$V_h,P_h$$ holds with a constant $$C_{\mathfrak{S},G}$$ provided that the parametric mesh can be covered with parametric macro-elements of a stable collection $$\mathfrak{C}$$ whose shape parameter is in the compact set $$\mathfrak{S}$$; the parametric macro-element diameters are smaller than a constant $$h_{\mathfrak{S},G}$$ depending on $$\mathfrak{S}$$ and $$G$$. Bressan & Sangalli (2013) proposed two types of stable parametric macro-elements: the TH type and the Sub-Grid type. Both depend on four parameters: the degrees and smoothnesses of the velocity and pressure. These parameters will be denoted according to the following table:   \[\begin{array}{r|cc} & V_\sigma & P_\sigma\\\hline \textrm{degree} & d_V & d_P\\ \textrm{smoothness} & s_V & s_P\end{array}\] In both types the equality $$d_V=d_p+1$$ ensures the same order of the approximation error for the velocity and for the pressure space, and so this is the preferred choice. In both cases, the shape parameter $$\sigma$$ contains the normalized knot-values of the knot vectors associated to the tensor-product spline space used to discretize the components of the velocity space. TH-type macro-elements are stable provided that   \[d_V-s_V> d_P-s_P.\] A parametric macro-element is described by $$n$$$$(L+1)$$-tuples of break points   \[ {\it{\Theta}}_i=[\theta_{i,1},\dots, \theta_{i,L+1}],\quad i=1,\dots,n \] with $$\theta_{1,1}=0,\theta_{1,L+1}=1$$, where   \begin{equation} L={\left\lceil{\frac{s_V+s_P+2}{(d_V-s_V)-(d_P-s_P)}}\right\rceil}. \end{equation} (2.28) The knot vectors of the tensor-product spline spaces for the velocity and pressure are obtained by repeating the break points according to $$d_V,s_V,d_P,s_P$$. The multiplicity of each break point in the velocity knot vectors is $$d_V-s_V$$. The multiplicity of the internal break point in the pressure knot vectors is $$d_P-s_P$$ and that of end knots is $$d_P+1$$. The shape parameter is   \[\sigma=(\theta_{1,2},\dots,\theta_{1,L},\theta_{2,1},\dots,\theta_{n,L+1}).\] SG-type macro-elements are stable provided that   \[2(d_V-s_V)> d_P-s_P.\] Again a parametric macro-element is described by $$n$$$$(L+1)$$-tuples of break points with   \begin{equation} L=2 {\left\lceil{\frac{s_V+s_P+2}{2(d_V-s_V)-(d_P-s_P)}}\right\rceil}. \end{equation} (2.29) The knot vectors of the velocity space are obtained as in the TH-case. Those of the pressure space by using every other break point so that every element of the pressure contains $$2^n$$ elements of the velocity. The shape parameter is as in the TH-type. The SG-type macro-elements can achieve maximal smoothness for both the velocity and the pressure spaces. 3. Application to hierarchical splines Hierarchical spline spaces were introduced by Forsey & Bartels (1988). They have been studied extensively by Giannelli et al. (2012), Giannelli & Jüttler (2013a), Giannelli & Jüttler (2013b), Mokriš et al. (2014), Giannelli et al. (2014) and Berdinsky et al. (2014). The application of hierarchical splines to IGA is under study by many authors including Vuong et al. (2011), Kuru et al. (2014), Bornemann & Cirak (2013), Buffa & Giannelli (2016) and Buffa et al. (2016). The main idea of the hierarchical spline basis is to enrich a tensor-product spline basis by adding functions from finer spaces while, simultaneously, removing redundant functions. Let $$\widehat S_0,\dots,\widehat S_\ell$$ be a sequence of tensor-product spline spaces (of the same degree) such that $$\widehat S_i \subseteq\widehat S_{i+1}$$ and $$\widehat{\mathfrak{B}}_0,\dots,\widehat{\mathfrak{B}}_\ell$$ be the corresponding bases. Let   \[\widehat{\it{\Omega}}_0=\widehat{\it{\Omega}}\supseteq\dots\supseteq\widehat{\it{\Omega}}_\ell\supseteq\widehat{\it{\Omega}}_{\ell+1}=\emptyset\] be a sequence of domains such that for $$i=1,\dots,\ell$$ the set $$\widehat{\it{\Omega}}_i$$ is a union of elements of the mesh $$\widehat{\mathscr{M}}_{i-1}$$ associated to $$\widehat S_{i-1}$$. The hierarchical basis $$\widehat{\mathscr{H}}$$ is defined by   \begin{equation} \widehat{\mathscr{H}} = \{ \hat\psi : \hat\psi \in \widehat{\mathfrak{B}}_i,\ {\text{supp}}\psi \subseteq \widehat{\it{\Omega}}_{i} \text{ and } {\text{supp}} \hat\psi\nsubseteq\widehat{\it{\Omega}}_{i+1}, i=0,\dots,\ell \}, \end{equation} (3.1) and the hierarchical space is $$\widehat{\mathbb{H}} = {\text{span}} \widehat{\mathscr{H}}$$. The truncated basis $$\widehat{\mathscr{H}}_T$$ is a modification of the standard basis. Each function $$\hat\psi {\,\in\,} \widehat{\mathscr{H}}$$ is replaced by its truncated version $$\hat\psi_T$$. The truncated basis has two important properties: it is a partition of unity and it preserves polynomial coefficients. This means that the coefficients of the expansion of a polynomial $$\hat g$$ in the truncated basis are the same as that of the tensor bases:   \begin{equation} \forall i=0,\dots,\ell,\quad \hat g=\sum_{\hat \psi \in \widehat{\mathfrak{B}}_i} c_{\hat\psi} \hat\psi \qquad{\Rightarrow}\qquad g=\sum_{\hat\psi_T \in \widehat{\mathscr{H}}_T} c_{\hat\psi} \hat\psi_T. \end{equation} (3.2) This property is important in the construction of the approximation operators. In the following, $$\widehat{\mathscr{M}}$$ and $${\mathscr{M}}$$ will refer to the parametric and physical meshes, i.e., to the collection of the elements of the space. 3.1 Approximation operator for hierarchical splines We describe an approximation operator $$\wp$$ that fulfils the requirements of Lemma 2.2. The operator is defined on the parametric domain and then it is pushed to the physical domain:   \[ \wp(\phi)=\hat\wp( \hat \phi)\circ G^{-1}=\hat\wp(\phi\circ G) \circ G^{-1}. \] The approximation operators on spline spaces are often described by providing the coefficients with respect to a basis. To avoid heavy notation, we set the analysis directly in the physical domain and we detail the construction for scalar valued functions; the extension to vector valued functions is obtained by applying the operator component-wise. Let $${\mathscr{H}}_T$$ be the push-forward of the truncated basis, then $$\wp:H^1({\it{\Omega}})\to {\mathbb{H}}$$ can be written as   \begin{equation} \wp \phi=\sum_{\psi\in{\mathscr{H}}_T} \lambda_\psi(\phi) \psi, \end{equation} (3.3) where the $$\lambda_\psi$$ are in $$H^1({\it{\Omega}})^*$$. Speleers & Manni (2016) showed that a set of $$\hat \lambda_\psi$$ for hierarchical bases can be constructed from those for the spaces $$ \widehat S_0,\dots, \widehat S_\ell$$. Moreover, following their construction, the following properties hold (rewritten in the physical domain): $$\wp$$ reproduces constant functions: $$\forall g$$ constant function, $$\displaystyle \sum_{\psi\in{\mathscr{H}}_T} \lambda_{\psi}(g) \psi= g$$. The operators $$\lambda_{\psi}$$ can be represented as functions in $$L^2({\it{\Omega}})$$ with local support and bounded norms. More precisely with support included in $${\text{supp}} \psi$$ and norm   \[ {\Vert{\lambda_{\psi}}\Vert}_{L^2({\it{\Omega}})} \leq C {|{{\text{supp}} \psi}|}^{-1/2}. \] The above imply the quasi-optimal approximation estimates with respect to $$L^\infty$$ norm (Speleers & Manni, 2016, Theorem 6). Here the result is rewritten for the Sobolev norms of interest, and a global estimate is derived as a corollary. Lemma 3.1 If $${\mathscr{M}}$$ is shape-regular, the approximation error of $$\wp$$ on a mesh element $$K$$ is bounded by   \begin{equation} {\Vert{\phi - \wp \phi }\Vert}_{L^2(K)}\leq C_{\mathrm{app}} h_{K} {|{\phi }|}_{H^1(\tilde K)}, \end{equation} (3.4) where $$\tilde K$$ is the union of the supports of the basis functions in $${\mathscr{H}}_T$$ that are active on $$K$$; $$C_{\mathrm{app}}$$ is a constant that depends on the regularity of the parametric mesh (2.26), (2.27) and of the parametrization $$G$$. Proof. From the reproduction of constants, we deduce   \begin{equation} \begin{aligned} {\Vert{\wp\phi -\phi }\Vert}_{L^2(K)} &= {\Vert{\wp (\phi -c) -(\phi -c)}\Vert}_{L^2(K)}\\ &\leq\inf_{c\in{\mathbb{R}}} \left ({\Vert{\wp (\phi -c)}\Vert}_{L^2(K)} + {\Vert{\phi -c}\Vert}_{L^2(K)} \right ) \end{aligned} \end{equation} (3.5) then, since $$\wp (\phi -c)|_{K}$$ depends only on $$(\phi -c)|_{\tilde K}$$,   \begin{equation} \begin{aligned} {\Vert{\wp\phi -\phi }\Vert}_{L^2(K)} &|\leq ({\Vert{\wp}\Vert} +1) \inf_c {\Vert{\phi -c}\Vert}_{L^2(\tilde K)}\\ &\leq ({\Vert{\wp}\Vert} +1) C_P(\tilde K) h_{\tilde K} {|{\phi }|}_{H^1(\tilde K)}, \end{aligned} \end{equation} (3.6) where $$C_P(\tilde K)$$ is the Poincaré constant of $$\tilde K$$ and depends on the shape of $$\tilde K$$; $$h_{\tilde K}$$ is the diameter of $$\tilde K$$; $${\Vert{\wp}\Vert}$$ is the norm of $$\wp$$ as an operator $$L^2(\tilde K)\to L^2(K)$$. Equation (3.4) follows from equation (3.6) if $$C_P(\tilde K) h_{\tilde K} \cong h_K$$ and $${\Vert{\wp}\Vert}$$ is uniformly bounded. The former is implied by the regularity assumptions on $$G$$ and the mesh. The latter is a consequence of the assumptions on the $$\lambda_\psi$$ and of the partition of the unity property of the truncated basis. □ Corollary 3.1 The approximation error of $$\wp$$ is bounded by:   \begin{equation} {\Vert{\rho( \phi - \wp \phi )}\Vert}_{L^2({\it{\Omega}})}^2 \leq {|{\phi }|}_{H^1({\it{\Omega}})}, \end{equation} (3.7) where $$\rho = C_{\mathrm{app}}^{-1} h_K^{-1} C_{\mathrm{over}}^{-1/2}$$ and $$\displaystyle C_{\mathrm{over}}=\max_{K\in {\mathscr{M}}}\#\{\psi\in{\mathscr{H}}: {\text{supp}}\psi \supseteq K\} \max_{\psi\in{\mathscr{H}}_T}\#\{K \subseteq{\text{supp}} \psi\}. $$ Proof. Rewrite (3.4) as $$ {\Vert{ C_{\mathrm{app}}^{-1} h_K^{-1}( \phi - \wp \phi ) }\Vert}_{L^2(K)}^2\leq {|{\phi }|}_{H^1(\tilde K)}^2$$ and sum over all elements in the mesh   \begin{eqnarray} {\Vert{C_{\mathrm{app}}^{-1} h_K^{-1}( \phi - \wp \phi )}\Vert}_{L^2({\it{\Omega}})}^2 &=&\sum_{K\in{\mathscr{M}}}{\Vert{ C_{\mathrm{app}}^{-1} h_K^{-1}( \phi - \wp \phi )}\Vert}_{L^2(K)}^2\nonumber\\ &\leq&\sum_{K\in{\mathscr{M}}}{|{\phi }|}_{H^1(\tilde K)}^2 = \sum_{K\in{\mathscr{M}}}\sum_{\substack{Q\in{\mathscr{M}}:\\Q \subseteq\tilde K}}{|{\phi }|}_{H^1(Q)}^2\nonumber\\ &=&\sum_{Q\in{\mathscr{M}}} \#\{K\in{\mathscr{M}}: Q \subseteq\tilde K\} {|{\phi }|}_{H^1(Q)}^2\\ &\leq&\sum_{Q\in{\mathscr{M}}} \#\{i: Q \subseteq{\text{supp}} \psi_i\}\#\{K \subseteq{\text{supp}} \psi_i\} {|{\phi }|}_{H^1(Q)}^2\nonumber\\ &\leq& C_{\mathrm{over}}{|{\phi }|}_{H^1({\it{\Omega}})}^2. \nonumber \end{eqnarray} (3.8) In the common case of levels corresponding to dyadic refinement in the parametric domain, the mesh regularity assumptions are equivalent to the fact that the boundaries of $${\it{\Omega}}_i$$ and $${\it{\Omega}}_{i+j}$$ are disjoint for $$j > \log_2 (C_{\mathrm{uniform}})$$. Unfortunately, this is not sufficient for a common upper bound of $$C_{\mathrm{over}}$$. All the required conditions are implied by admissibility of the mesh as defined by Buffa & Giannelli (2016). 3.2 TH- and SG-type methods with hierarchical splines Let $$\widehat V_0\subset \dots\subset\widehat V_\ell$$ and $$\widehat P_0\subset\dots\subset\widehat P_\ell$$ be two hierarchies of tensor-product spline spaces on $$\widehat{\it{\Omega}}$$ such that the pairs $$V_i,P_i$$ are either all of type TH or all of type SG. Then let $$\widehat{\mathbb{H}}_V$$ and $$\widehat{\mathbb{H}}_P$$ be, respectively, associated spaces $$\{\widehat V_i\}$$ and $$\{\widehat P_i\}$$ of the hierarchical spline spaces and a common hierarchy of nested domains $$\widehat {\it{\Omega}}_0\supseteq\dots\supseteq\widehat{\it{\Omega}}_\ell\subseteq\widehat{\it{\Omega}}_{\ell+1}=\emptyset$$. The discretization spaces are the push-forward of $$\widehat{\mathbb{H}}_V$$ and $$\widehat{\mathbb{H}}_P$$:   \begin{equation} \begin{aligned} V_h={\mathbb{H}}_V^n\cap V=\{\hat{{\mathbf{w}}} \circ G^{-1}:\hat{{\mathbf{w}}}\in \widehat{\mathbb{H}}_V^n \}\cap V,\\ P_h={\mathbb{H}}_P\cap P=\{\hat q \circ G^{-1}:\hat{q}\in \widehat{\mathbb{H}}_P \}\cap P. \end{aligned} \end{equation} (3.9) For these spaces, stability follows from Theorem 2.3 if $${\it{\Omega}}$$ can be covered by macro-elements. Theorem 3.1 If the parametric mesh can be covered by macro-elements of type SG or TH in a collection $$\mathfrak{S}$$ defined by shape regularity and local quasi-uniformity, and such that the diameter of the covering macro-elements is less then $$h_{\mathfrak{S}, G}$$ (see Lemma 2.4), then the method is inf–sup stable and the inf–sup constant $${\beta}_h$$ depends only on $$\mathfrak{S}$$, $$G$$ and $$C_{\mathrm{over}}$$. Consequently, the a priori estimates for the error (2.10) and (2.11) hold. The request of a macro-element covering constrains the class of meshes to which the result can be applied. Indeed, TH and SG macro-elements have a tensor-product structure, and thus it is necessary that each ring   \begin{equation} {\it{\Delta}}_i= {\it{\Omega}}_i\setminus {\it{\Omega}}_{i+1} \end{equation} (3.10) is covered separately. In practice, given a sequence of domains $$ {\it{\Omega}}_0'\supseteq\dots\supseteq {\it{\Omega}}_\ell'\supseteq{\it{\Omega}}_{\ell+1}'=\emptyset$$ that does not satisfy the condition above, for instance, the result of a refinement step, it is possible to construct a corresponding sequence $$ {\it{\Omega}}_0\supseteq\dots\supseteq {\it{\Omega}}_\ell\supseteq{\it{\Omega}}_{\ell+1}=\emptyset$$ that satisfies Theorem 3.1. A possible solution for dyadic refinement is to associate to each level $$i$$ a partition $$\mathscr{P}_i$$ of $${\it{\Omega}}$$ into nonoverlapping closed macro-elements. The different partitions must be compatible, meaning that each macro-element in $$\mathscr{P}_i$$ is the union of $$2^n$$ macro-elements in $$\mathscr{P}_{i+1}$$. Lemma 3.2 For any $$ {\it{\Omega}}_0'\supseteq\dots\supseteq {\it{\Omega}}_\ell'$$, the sequence defined is as follows. For $$i=1,\dots,\ell$$  \begin{equation} {\it{\Omega}}_i= \mathscr{E}_{i-1} ({\it{\Omega}}_{i+1}\cup {\it{\Omega}}_i'), \end{equation} (3.11) where $$\mathscr{E}_i (Q) = \bigcup \{M \in \mathscr{P}_i : M\cap \bar{Q}\neq\emptyset\} $$ has rings $${\it{\Delta}}_i= {\it{\Omega}}_i\setminus {\it{\Omega}}_{i+1}$$ that can be covered by macro-elements from $$\mathscr{P}_{i}$$. Moreover, each element of $${\mathscr{M}}$$ is contained in exactly one macro-element ($$m=1$$ in Lemma 2.3). Proof. By construction both $${\it{\Omega}}_i$$ and $${\it{\Omega}}_{i+1}$$ are unions of macro-elements of $$\mathscr{P}_{i}$$. Since the macro-elements are nonoverlapping, then $${\it{\Delta}}_i$$ is covered by the   \[ \bigcup\{ M\in \mathscr{P}_{i}: M \subseteq {\it{\Omega}}_i \text{ and } M \not\subseteq {\it{\Omega}}_{i+1}\}. \] □ The construction in Lemma 3.2 ensures that on each element $$K$$ of level $$i$$ are active only functions $$\psi\in \mathfrak{B}_j$$ with $$j=i,i-1$$. This means that the resulting mesh is 2-admissible according to Buffa et al. (2016). Moreover, in the common case of dyadic refinement the construction implies both shape-regularity and local-quasi-uniformity. In conclusion, we showed that it is possible to construct and refine a pair of discretization spaces in such a way that the hypotheses of Theorem 3.1 are fulfilled with a common lower bound on $$\beta_h$$. Thus the convergence discretization error is controlled by the best approximation error up to a constant factor. 4. Numerical stability tests The inf–sup constant $${\beta}_h$$ can be computed numerically by solving the following eigen problem:   \begin{equation} {B}_h N_{V}^{-1} {B}_h^T q = \lambda N_{P} q, \end{equation} (4.1) where $${B}_h: V_h \to P_h^*$$, $$N_{V}: V_h \to V_h^* $$ and $$N_{P}: P_h \to P_h^*$$ are, respectively, the matrices corresponding to $$ {b}_h$$, the scalar product on $$V$$ and the scalar product on $$P$$. Let $$\lambda_{\min}$$ be the smallest eigenvalue of problem (4.1) then, as described in Brezzi & Fortin (1991),   \[{\beta}_h=\sqrt{\lambda_{\min}}.\] The tests focus on the dependence of $${\beta}_h$$ on the polynomial degree ($$d_P$$), the number of levels of the hierarchical space $$(\ell)$$, the type of discretization (TH or SG) and the hypothesis of Theorem 3.1. In particular, $${\beta}_h$$ is computed on meshes both fulfilling and not fulfilling the hypothesis of Theorem 3.1. In all tests, the domain is $${\it{\Omega}}=\widehat{\it{\Omega}}=[0,1]^2$$ parametrized by the identity, the spaces are balanced and of maximum smoothness ($$d_V=d_P+1,\,s_P=d_P-1$$) and the refinement of the spaces $$S_i$$ is dyadic. 4.1 Macro-element inf–sup As a first test, we computed the inf–sup on a single 2D macro-element. The results can be read in Table 1. We observe that the inf–sup constant $${\beta}_h$$ decreases fast as a function of the degree (and the smoothness). The result can be interpreted by observing that on the macro-element the velocities are constrained to decrease to zero with $$C^{s_V}$$ smoothness toward the boundary. This condition is unnatural in $$H^1({\it{\Omega}})^n$$ and can be enforced only on the discrete space. The effect of the degree and of the smoothness can be made clearer by the following 1D example. Table 1 Numerically computed inf–sup on 2D macro-elements for various pressure degrees $$d_P$$, $$d_V=d_P+1$$. For macro-elements of type TH both pressure and velocity are $$C^{d_P-1}$$ while for SG-type the pressure is $$C^{d_P-1}$$ and velocity $$C^{d_P}$$ $$d_P$$  2  3  4  5  6  TH  $$4.68 \ 10^{-1}$$  $$ 4.86 \ 10^{-2} $$  $$4.96 \ 10^{-4}$$  $$ 1.72 \ 10^{-6}$$  $$3.23 \ 10^{-9}$$  SG  $$1.73 \ 10^{-1}$$  $$ 2.98\ 10^{-3} $$  $$8.90 \ 10^{-6}$$  $$ 1.06 \ 10^{-8}$$  $$7.00 \ 10^{-12}$$  $$d_P$$  2  3  4  5  6  TH  $$4.68 \ 10^{-1}$$  $$ 4.86 \ 10^{-2} $$  $$4.96 \ 10^{-4}$$  $$ 1.72 \ 10^{-6}$$  $$3.23 \ 10^{-9}$$  SG  $$1.73 \ 10^{-1}$$  $$ 2.98\ 10^{-3} $$  $$8.90 \ 10^{-6}$$  $$ 1.06 \ 10^{-8}$$  $$7.00 \ 10^{-12}$$  Table 1 Numerically computed inf–sup on 2D macro-elements for various pressure degrees $$d_P$$, $$d_V=d_P+1$$. For macro-elements of type TH both pressure and velocity are $$C^{d_P-1}$$ while for SG-type the pressure is $$C^{d_P-1}$$ and velocity $$C^{d_P}$$ $$d_P$$  2  3  4  5  6  TH  $$4.68 \ 10^{-1}$$  $$ 4.86 \ 10^{-2} $$  $$4.96 \ 10^{-4}$$  $$ 1.72 \ 10^{-6}$$  $$3.23 \ 10^{-9}$$  SG  $$1.73 \ 10^{-1}$$  $$ 2.98\ 10^{-3} $$  $$8.90 \ 10^{-6}$$  $$ 1.06 \ 10^{-8}$$  $$7.00 \ 10^{-12}$$  $$d_P$$  2  3  4  5  6  TH  $$4.68 \ 10^{-1}$$  $$ 4.86 \ 10^{-2} $$  $$4.96 \ 10^{-4}$$  $$ 1.72 \ 10^{-6}$$  $$3.23 \ 10^{-9}$$  SG  $$1.73 \ 10^{-1}$$  $$ 2.98\ 10^{-3} $$  $$8.90 \ 10^{-6}$$  $$ 1.06 \ 10^{-8}$$  $$7.00 \ 10^{-12}$$  Example 4.1 Consider a 1D element of length $$1$$, $$K=[0,1]$$. Now let $$M=[0,L]$$ be a macro-element and let   \[q=\begin{cases} (1-x)^{d_P} & x\leq 1\\ 0 & x>1 \end{cases} \] be the pressure. Then   \[ {b}(w,q) = c d_P \binom{d_P+r_V+1}{d_P}^{-1}, \] where $$c$$ is the coefficient of $$x^{r_V}$$ in $$w|_{[0,1]}$$. This is not compensated by the behaviour of the norms:   \[ {\Vert{\rho {\mathrm{grad}\,{q}}}\Vert}_{L^2(M)}=d_P (2d_P-1)^{-1/2}, \qquad {\Vert{w}\Vert}_{L^2(M)}\geq {\Vert{w}\Vert}_{L^2([0,1])}=c (2r_V+1)^{-1/2}. \] 4.2 Graded meshes The meshes of this test are refined at the centre in such a way that the $${\it{\Omega}}_i$$ are almost concentric squares and the $${\it{\Delta}}_i={\it{\Omega}}_i\setminus{\it{\Omega}}_{i+1}$$ are rings with a thickness of at least $$2 d_P+1$$ pressure elements. In this way, each $${\it{\Delta}}_i$$ can be covered with macro-elements. See the meshes in Fig. 1. The results are reported in Table 2 for the TH element and in Table 3 for the SG element. Observe that the discrete inf–sup constant is independent of the number of level as predicted by the theory. An interesting fact is that the computed discrete inf–sup constants exceed the continuous inf–sup constant that according to Costabel & Dauge (2015) and Costabel et al. (2015) is contained in   \[ \left [\sin\left(\frac{\pi}{8}\right), \sqrt{\frac{1}{2}-\frac{1}{\pi}}\right]. \] Fig. 1. View largeDownload slide Graded meshes. First row: the meshes of velocity (left) and pressure (right) for the SG method with parameters $$d_v=3, d_P=2,s_V=2,s_P=1$$ and 5 levels. Second row: the meshes of velocity (left) and pressure (right) for the TH method with parameters $$d_v=4,d_P=3,s_V=s_P=2$$ and 4 levels. Fig. 1. View largeDownload slide Graded meshes. First row: the meshes of velocity (left) and pressure (right) for the SG method with parameters $$d_v=3, d_P=2,s_V=2,s_P=1$$ and 5 levels. Second row: the meshes of velocity (left) and pressure (right) for the TH method with parameters $$d_v=4,d_P=3,s_V=s_P=2$$ and 4 levels. Table 2 Computed inf–sup constant on graded mesh for the TH method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=s_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  2  0.4600  0.4574  0.4568  0.4562  0.4562  0.4562  0.4562  3  0.4537  0.4521  0.4514  0.4511  0.4508  0.4508  0.4508  4  0.4500  0.4485  0.4479  0.4477  0.4475  0.4475  0.4475  $$d_P$$  2  3  4  5  6  7  8  2  0.4600  0.4574  0.4568  0.4562  0.4562  0.4562  0.4562  3  0.4537  0.4521  0.4514  0.4511  0.4508  0.4508  0.4508  4  0.4500  0.4485  0.4479  0.4477  0.4475  0.4475  0.4475  Table 2 Computed inf–sup constant on graded mesh for the TH method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=s_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  2  0.4600  0.4574  0.4568  0.4562  0.4562  0.4562  0.4562  3  0.4537  0.4521  0.4514  0.4511  0.4508  0.4508  0.4508  4  0.4500  0.4485  0.4479  0.4477  0.4475  0.4475  0.4475  $$d_P$$  2  3  4  5  6  7  8  2  0.4600  0.4574  0.4568  0.4562  0.4562  0.4562  0.4562  3  0.4537  0.4521  0.4514  0.4511  0.4508  0.4508  0.4508  4  0.4500  0.4485  0.4479  0.4477  0.4475  0.4475  0.4475  Table 3 Computed inf–sup constant on graded mesh for the SG method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=d_P,s_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  2  0.4664  0.4632  0.4623  0.4616  0.4616  0.4616  0.4616  3  0.4578  0.4558  0.4549  0.4546  0.4542  0.4542  0.4542  4  0.4531  0.4513  0.4506  0.4504  0.4502  0.4502  0.4502  $$d_P$$  2  3  4  5  6  7  8  2  0.4664  0.4632  0.4623  0.4616  0.4616  0.4616  0.4616  3  0.4578  0.4558  0.4549  0.4546  0.4542  0.4542  0.4542  4  0.4531  0.4513  0.4506  0.4504  0.4502  0.4502  0.4502  Table 3 Computed inf–sup constant on graded mesh for the SG method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=d_P,s_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  2  0.4664  0.4632  0.4623  0.4616  0.4616  0.4616  0.4616  3  0.4578  0.4558  0.4549  0.4546  0.4542  0.4542  0.4542  4  0.4531  0.4513  0.4506  0.4504  0.4502  0.4502  0.4502  $$d_P$$  2  3  4  5  6  7  8  2  0.4664  0.4632  0.4623  0.4616  0.4616  0.4616  0.4616  3  0.4578  0.4558  0.4549  0.4546  0.4542  0.4542  0.4542  4  0.4531  0.4513  0.4506  0.4504  0.4502  0.4502  0.4502  4.3 Nonmacro-element-covered meshes The meshes of this test are refined similarly to those in the previous test, but the thickness of the $${\it{\Delta}}_i$$ is not enough to allow for a macro-element covering (see Fig. 2). As a consequence, Theorem 3.1 does not apply. The results are reported in Table 4 for the TH element and in Table 5 for the SG element. Surprisingly, the inf–sup constant in this case is greater (and thus better) than in the previous case. Our interpretation is that this refinement increases the dimension of the velocity space compared to that of the pressure space even more, and hence the inf–sup constant increases. Fig. 2. View largeDownload slide Nonmacro-elements-covered meshes. On the left, the velocity mesh for the SG method with parameters $$d_v=5,d_P=4,s_V=4,s_P=3$$ and 7 levels; on the right, the velocity mesh for TH method with the same parameters. Fig. 2. View largeDownload slide Nonmacro-elements-covered meshes. On the left, the velocity mesh for the SG method with parameters $$d_v=5,d_P=4,s_V=4,s_P=3$$ and 7 levels; on the right, the velocity mesh for TH method with the same parameters. Table 4 Computed inf–sup constant on nonmacro-elements-covered meshes for the TH method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=s_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  2  0.4705  0.4705  0.4705  0.4705  0.4705  0.4705  0.4705  3  0.4604  0.4593  0.4593  0.4593  0.4593  0.4593  0.4593  4  0.4564  0.4564  0.4564  0.4564  0.4564  0.4564  0.4564  $$d_P$$  2  3  4  5  6  7  8  2  0.4705  0.4705  0.4705  0.4705  0.4705  0.4705  0.4705  3  0.4604  0.4593  0.4593  0.4593  0.4593  0.4593  0.4593  4  0.4564  0.4564  0.4564  0.4564  0.4564  0.4564  0.4564  Table 4 Computed inf–sup constant on nonmacro-elements-covered meshes for the TH method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=s_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  2  0.4705  0.4705  0.4705  0.4705  0.4705  0.4705  0.4705  3  0.4604  0.4593  0.4593  0.4593  0.4593  0.4593  0.4593  4  0.4564  0.4564  0.4564  0.4564  0.4564  0.4564  0.4564  $$d_P$$  2  3  4  5  6  7  8  2  0.4705  0.4705  0.4705  0.4705  0.4705  0.4705  0.4705  3  0.4604  0.4593  0.4593  0.4593  0.4593  0.4593  0.4593  4  0.4564  0.4564  0.4564  0.4564  0.4564  0.4564  0.4564  Table 5 Computed inf–sup constant on nonmacro-elements-covered graded meshes for the SG method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=d_P,s_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  2  0.4817  0.4817  0.4817  0.4817  0.4817  0.4817  0.4817  3  0.4661  0.4645  0.4645  0.4645  0.4645  0.4645  0.4645  4  0.4610  0.4610  0.4610  0.4610  0.4610  0.4610  0.4610  $$d_P$$  2  3  4  5  6  7  8  2  0.4817  0.4817  0.4817  0.4817  0.4817  0.4817  0.4817  3  0.4661  0.4645  0.4645  0.4645  0.4645  0.4645  0.4645  4  0.4610  0.4610  0.4610  0.4610  0.4610  0.4610  0.4610  Table 5 Computed inf–sup constant on nonmacro-elements-covered graded meshes for the SG method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=d_P,s_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  2  0.4817  0.4817  0.4817  0.4817  0.4817  0.4817  0.4817  3  0.4661  0.4645  0.4645  0.4645  0.4645  0.4645  0.4645  4  0.4610  0.4610  0.4610  0.4610  0.4610  0.4610  0.4610  $$d_P$$  2  3  4  5  6  7  8  2  0.4817  0.4817  0.4817  0.4817  0.4817  0.4817  0.4817  3  0.4661  0.4645  0.4645  0.4645  0.4645  0.4645  0.4645  4  0.4610  0.4610  0.4610  0.4610  0.4610  0.4610  0.4610  4.4 Nongraded meshes In this test, the meshes are refined only in half of the domain so that the ratio between the size of neighbouring elements increases geometrically on the number of levels (Fig. 3). In this case, not only the macro-element covering is missing, but also the approximation properties of $$\wp$$ deteriorates as the number of levels increases. The results are presented in the Tables 6 and 7. Again we observe that the inf–sup constant is stable on both degree and the number of levels. Fig. 3. View largeDownload slide Nongraded meshes. Pressure meshes for $$d_P=2,s_P=1$$ with 2 levels (left), and 8 levels (right). Fig. 3. View largeDownload slide Nongraded meshes. Pressure meshes for $$d_P=2,s_P=1$$ with 2 levels (left), and 8 levels (right). Table 6 Computed inf–sup constant on nongraded meshes for the TH method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=s_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  1  0.4683  0.4716  0.4690  0.4631  0.4627  0.4623  0.4621  2  0.4685  0.4728  0.4713  0.4629  0.4628  0.4559  0.4559  3  0.4757  0.4669  0.4599  0.4599  0.4536  0.4536  0.4536  4  0.4557  0.4598  0.4570  0.4567  0.4516  0.4516  0.4516  $$d_P$$  2  3  4  5  6  7  8  1  0.4683  0.4716  0.4690  0.4631  0.4627  0.4623  0.4621  2  0.4685  0.4728  0.4713  0.4629  0.4628  0.4559  0.4559  3  0.4757  0.4669  0.4599  0.4599  0.4536  0.4536  0.4536  4  0.4557  0.4598  0.4570  0.4567  0.4516  0.4516  0.4516  Table 6 Computed inf–sup constant on nongraded meshes for the TH method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=s_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  1  0.4683  0.4716  0.4690  0.4631  0.4627  0.4623  0.4621  2  0.4685  0.4728  0.4713  0.4629  0.4628  0.4559  0.4559  3  0.4757  0.4669  0.4599  0.4599  0.4536  0.4536  0.4536  4  0.4557  0.4598  0.4570  0.4567  0.4516  0.4516  0.4516  $$d_P$$  2  3  4  5  6  7  8  1  0.4683  0.4716  0.4690  0.4631  0.4627  0.4623  0.4621  2  0.4685  0.4728  0.4713  0.4629  0.4628  0.4559  0.4559  3  0.4757  0.4669  0.4599  0.4599  0.4536  0.4536  0.4536  4  0.4557  0.4598  0.4570  0.4567  0.4516  0.4516  0.4516  Table 7 Computed inf–sup constant on nongraded meshes for the SG method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=d_P,_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  1  0.5256  0.4981  0.4985  0.4804  0.4804  0.4804  0.4804  2  0.4952  0.4827  0.4827  0.4703  0.4703  0.4614  0.4614  3  0.4901  0.4759  0.4648  0.4651  0.4577  0.4577  0.4577  4  0.4764  0.4707  0.4614  0.4615  0.4551  0.4551  0.4551  $$d_P$$  2  3  4  5  6  7  8  1  0.5256  0.4981  0.4985  0.4804  0.4804  0.4804  0.4804  2  0.4952  0.4827  0.4827  0.4703  0.4703  0.4614  0.4614  3  0.4901  0.4759  0.4648  0.4651  0.4577  0.4577  0.4577  4  0.4764  0.4707  0.4614  0.4615  0.4551  0.4551  0.4551  Table 7 Computed inf–sup constant on nongraded meshes for the SG method depending on degree and number of levels. The parameters are $$d_V=d_P+1,s_V=d_P,_P=d_P-1$$ $$d_P$$  2  3  4  5  6  7  8  1  0.5256  0.4981  0.4985  0.4804  0.4804  0.4804  0.4804  2  0.4952  0.4827  0.4827  0.4703  0.4703  0.4614  0.4614  3  0.4901  0.4759  0.4648  0.4651  0.4577  0.4577  0.4577  4  0.4764  0.4707  0.4614  0.4615  0.4551  0.4551  0.4551  $$d_P$$  2  3  4  5  6  7  8  1  0.5256  0.4981  0.4985  0.4804  0.4804  0.4804  0.4804  2  0.4952  0.4827  0.4827  0.4703  0.4703  0.4614  0.4614  3  0.4901  0.4759  0.4648  0.4651  0.4577  0.4577  0.4577  4  0.4764  0.4707  0.4614  0.4615  0.4551  0.4551  0.4551  4.5 Observation on the numerically computed inf–sup constants The computed inf–sup constant on macro-elements shows that the lower bound obtained by the theory is very pessimistic. The strong dependence on the smoothness is not observed and the lack of a macro-element cover does not influence the inf–sup constant. Instead, the computed inf–sup constant of the TH- and SG-type methods are close to the continuous inf–sup constant, sometimes even greater, regardless of the degree, the mesh grading and the macro-element covering. 5. Adaptive implementation The main motivation for the extension of the inf–sup stability to hierarchical splines is adaptivity. In this section, we report a numerical experiment based on the TH adaptive IGA scheme. We solved a 2D flow inside a T-shaped pipe (Fig. 4). The forcing term is constant and directed along the $$x$$ axis. The homogeneous Dirichlet condition is imposed on the boundary except for the inlet $${\it{\Gamma}}_{\mathrm{in}}$$ and the outlet $${\it{\Gamma}}_{\mathrm{out}}$$ where homogeneous Neumann conditions are enforced. Fig. 4. View largeDownload slide The domain and its subdivision in four patches. Fig. 4. View largeDownload slide The domain and its subdivision in four patches. The regularity of the continuous solution is ruined by the presence of two concave corners in the domain. More precisely, the solution can be decomposed in a regular part $${\mathbf{u}}_{\mathrm{reg}}$$, whose regularity depends only on ($${\mathbf{f}}$$ and $$g$$), and a singular part $${\mathbf{u}}_{\mathrm{sing}}$$, whose regularity depends on the geometry of the domain. The behaviour of the singular part in the neighbourhood of a concave corner is described by   \begin{equation} {\mathbf{u}}_{\mathrm{sing}}=\sum_{0 = q}^Q r^{\lambda}\psi_q(\theta)\log(r)^q, \end{equation} (5.1) where $$\lambda$$ and $$Q$$ can be computed from the geometry (Dauge, 1989). In the test case (a corner singularity in 2D with an angle of $$\frac{3}{2}\pi$$), the exponent $$\lambda $$ is the smallest positive solution of   \begin{equation} \sin\Big(\frac{3}{2}\pi \lambda^2\Big) - \lambda^2 =0 \end{equation} (5.2) which approximately equals $$0.54448...$$. This means that the solution $$({\mathbf{u}},p)$$ is at most in $$H^{\lambda+1}({\it{\Omega}})^2 \times H^{\lambda}({\it{\Omega}})$$ because higher order derivatives are not square integrable. Consequently, the highest order of convergence that can be achieved with uniform refinement is   \[\text{dofs}^{-\frac{\lambda}{2}} \approx\text{dofs}^{-0.272...}.\] The domain is represented by a 4-patch geometry as shown in Fig. 4. Each patch corresponds to a different parametrization; in this case, the parametrizations are four different translations. The discretization spaces on each patch are of TH-type with parameters $$d_V=3,d_P=2,s_V=s_P=1$$. The spaces on the multipatch domain are obtained by gluing the patch-local spaces using the technique proposed by Buchegger et al. (2016). In this way, the basis functions are $$C^1$$ on the whole domain. The error is estimated using the standard residual-based estimator, but due to the $$C^1$$ continuity of the basis function the jump terms are 0 and can be ignored. Elements are marked according to Dörfler’s technique with parameter $$1/5$$, i.e., the error estimates of the refined elements sum up to a $$1/5$$ of the total estimated error. The mesh was correctly refined near the singularities as can be observed in Fig. 5. The convergence of the error estimator for both uniform and adaptive refinement is reported in Fig. 6. Uniform refinement achieves the predicted order of convergence ($$0.27...)$$ while adaptive refinement achieves an order of $$0.94...$$ on average and in line with a theoretical convergence rate of $$\frac{3}{2}$$ in the last two iterations. Figure 7 shows the magnitude of the velocity. Fig. 5. View largeDownload slide Mesh after 10 refinements steps. Fig. 5. View largeDownload slide Mesh after 10 refinements steps. Fig. 6. View largeDownload slide Convergence plot of the adaptive and uniform refining strategy for the problem on the T-shaped domain. Fig. 6. View largeDownload slide Convergence plot of the adaptive and uniform refining strategy for the problem on the T-shaped domain. Fig. 7. View largeDownload slide Magnitude of the computed velocity field. Fig. 7. View largeDownload slide Magnitude of the computed velocity field. 6. Summary This article extends the stability result of the isogeometric TH- and SG-type methods to hierarchical splines. This is achieved by defining an appropriate approximation operator and then following the same techniques used by Bressan & Sangalli (2013). The inf–sup stability can be achieved for domains of any dimension, and the method is thus applicable to both 2D and 3D problems. The convergence of a corresponding adaptive TH scheme was tested numerically on a 2D problem with corner singularities giving good results. Numerical experiments hint that the lower bound for the discrete inf–sup constant $${\beta}_h$$ is independent of the degree, mesh-quasi-uniformity and other parameters appearing in the theoretical analysis. From this point of view, the theory is not completely satisfactory. Funding European Research Council (GA 324340 (EXAMPLE)) and (GA 678727 (MOTOR)). References Acosta G., Durán R. G. & Muschietti M. A. ( 2006) Solutions of the divergence operator on John domains. Adv. Math. , 206, 373– 401. Google Scholar CrossRef Search ADS   Bazilevs Y., Beirão da Veiga L., Cottrell J. A., Hughes T. J. R. & Sangalli G. ( 2006) Isogeometric analysis: approximation, stability and error estimates for h-refined meshes. Math. Models Methods Appl. Sci. , 16, 1031– 1090. Google Scholar CrossRef Search ADS   Beirão da Veiga L., Buffa A., Sangalli G. & Vázquez R. ( 2014) Mathematical analysis of variational isogeometric methods. Acta Numer. , 23, 157– 287. Google Scholar CrossRef Search ADS   Berdinsky D., Kim T.-W., Bracco C., Cho D., Mourrain B., Oh M.-J. & Kiatpanichgij S. ( 2014) Dimensions and bases of hierarchical tensor-product splines. J. Comput. Appl. Math. , 257, 86– 104. Google Scholar CrossRef Search ADS   Bornemann P. B. & Cirak F. ( 2013) A subdivision-based implementation of the hierarchical b-spline finite element method. Comput. Methods Appl. Mech. Engrg. , 253, 584– 598. Google Scholar CrossRef Search ADS   Bressan A. ( 2011) Isogeometric regular discretization for the Stokes problem. IMA J. Numer. Anal. , 31, 1334– 1356. Google Scholar CrossRef Search ADS   Bressan A. & Sangalli G. ( 2013) Isogeometric discretizations of the Stokes problem: stability analysis by the macroelement technique. IMA J. Numer. Anal. , 33, 629– 651. Google Scholar CrossRef Search ADS   Brezzi F. & Fortin M. ( 1991) Mixed and Hybrid Finite Element Methods . Springer Series in Computational Mathematics, vol. 15. New York: Springer, pp. x+350. Google Scholar CrossRef Search ADS   Buchegger F., Jüttler B. & Mantzaflaris A. ( 2016) Adaptively refined multi-patch B-splines with enhanced smoothness. Appl. Math. Comput. , 272, 159– 172. Buffa A. & Giannelli C. ( 2016) Adaptive isogeometric methods with hierarchical splines: error estimator and convergence. Math. Models Methods Appl. Sci. , 26, 1– 25. Google Scholar CrossRef Search ADS   Buffa A., Giannelli C., Morgenstern P. & Peterseim D. ( 2016) Complexity of hierarchical refinement for a class of admissible mesh configurations. Comput. Aided Geom. Design , 47, 83– 92. Google Scholar CrossRef Search ADS   Buffa A., Rivas J., Sangalli G. & Vázquez R. ( 2011) Isogeometric discrete differential forms in three dimensions. SIAM J. Numer. Anal. , 49, 818– 844. Google Scholar CrossRef Search ADS   Buffa A., Sangalli G. & Vázquez R. ( 2014) Isogeometric methods for computational electromagnetics: B-spline and T-spline discretizations. J. Comput. Phys. , 257, 1291– 1320. Google Scholar CrossRef Search ADS   Costabel M., Crouzeix M., Dauge M. & Lafranche Y. ( 2015) The inf-sup constant for the divergence on corner domains. Numer. Methods Partial Differential Equations , 31, 439– 458. Google Scholar CrossRef Search ADS   Costabel M. & Dauge M. ( 2015) On the inequalities of Babuška–Aziz, Friedrichs and Horgan–Payne. Arch. Ration. Mech. Anal. , 217, 873– 898. Google Scholar CrossRef Search ADS   Cottrell J. A., Hughes T. J. R. & Bazilevs Y. ( 2009) Isogeometric Analysis: Toward Integration of CAD and FEA . Chichester, England: John Wiley & Sons. Google Scholar CrossRef Search ADS   Dauge M. ( 1989) Stationary Stokes and Navier-Stokes systems on two- or three-dimensional domains with corners. I. Linearized equations. SIAM J. Math. Anal. , 20, 74– 97. Google Scholar CrossRef Search ADS   Engelman M. & Bercovier M. ( 1979) Estimation de l’erreur pour la résolution du problème de Stokes par un élément fini quadrilatéral conforme sur les vitesses. C. R. Acad. Sci. Paris Sér. A-B , 288, A555– A557. Evans J. A. & Hughes T. J. R. ( 2013a) Isogeometric divergence-conforming B-splines for the Darcy–Stokes–Brinkman equations. Math. Models Methods Appl. Sci. , 23, 671– 741. Google Scholar CrossRef Search ADS   Evans J. A. & Hughes T. J. R. ( 2013b) Isogeometric divergence-conforming B-splines for the steady Navier–Stokes equations. Math. Models Methods Appl. Sci. , 23, 1421– 1478. Google Scholar CrossRef Search ADS   Evans J. A., Scott M. A., Shepherd K., Thomas D. & Vázquez R. ( 2017) Hierarchical B-spline complexes of discrete differential forms. ArXiv preprint arXiv:1708.04195. Forsey D. R. & Bartels R. H. ( 1988) Hierarchical B-spline refinement. SIGGRAPH Comput. Graph. , 22, 205– 212. Google Scholar CrossRef Search ADS   Giannelli C. & Jüttler B. ( 2013a) Bases and dimensions of bivariate hierarchical tensor-product splines. J. Comput. Appl. Math. , 239, 162– 178. Google Scholar CrossRef Search ADS   Giannelli C. & Jüttler B. ( 2013b) Local and adaptive refinement with hierarchical B-splines. Boll. Unione Mat. Ital. (9) , 6, 735– 740. Giannelli C., Jüttler B. & Speleers H. ( 2012) THB-splines: the truncated basis for hierarchical splines. Comput. Aided Geom. Design , 29, 485– 498. Google Scholar CrossRef Search ADS   Giannelli C., Jüttler B. & Speleers H. ( 2014) Strongly stable bases for adaptively refined multilevel spline spaces. Adv. Comput. Math. , 40, 459– 490. Google Scholar CrossRef Search ADS   Hughes T. J. R., Cottrell J. A. & Bazilevs Y. ( 2005) Isogeometric analysis: CAD, finite elements, NURBS, exact geometry and mesh refinement. Comput. Methods Appl. Mech. Engrg. , 194, 4135– 4195. Google Scholar CrossRef Search ADS   Johannessen K. A., Kumar M. & Kvamsdal T. ( 2015) Divergence-conforming discretization for Stokes problem on locally refined meshes using LR B-splines. Comput. Methods Appl. Mech. Engrg. , 293, 38– 70. Google Scholar CrossRef Search ADS   Kuru G., Verhoosel C. V., van der Zee K. G. & van Brummelen E. H. ( 2014) Goal-adaptive isogeometric analysis with hierarchical splines. Comput. Methods Appl. Mech. Engrg. , 270, 270– 292. Google Scholar CrossRef Search ADS   Mokriš D., Jüttler B. & Giannelli C. ( 2014) On the completeness of hierarchical tensor-product B-splines. J. Comput. Appl. Math. , 271, 53– 70. Google Scholar CrossRef Search ADS   Rüberg T. & Cirak F. ( 2012) Subdivision-stabilised immersed b-spline finite elements for moving boundary flows. Comput. Methods Appl. Mech. Engrg. , 209/212, 266– 283. Google Scholar CrossRef Search ADS   Speleers H. & Manni C. ( 2016) Effortless quasi-interpolation in hierarchical spaces. Numer. Math. , 132, 155– 184. Google Scholar CrossRef Search ADS   Stenberg R. ( 1990) Error analysis of some finite element methods for the Stokes problem. Math. Comp. , 54, 495– 508. Google Scholar CrossRef Search ADS   Verfürth R. ( 1984) Error estimates for a mixed finite element approximation of the Stokes equations. RAIRO Anal. Numér. , 18, 175– 182. Google Scholar CrossRef Search ADS   Vuong A.-V., Giannelli C., Jüttler B. & Simeon B. ( 2011) A hierarchical approach to adaptive local refinement in isogeometric analysis. Comput. Methods Appl. Mech. Engrg. , 200, 3554– 3567. Google Scholar CrossRef Search ADS   © The authors 2017. Published by Oxford University Press on behalf of the Institute of Mathematics and its Applications. All rights reserved. This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/about_us/legal/notices) For permissions, please e-mail: journals. permissions@oup.com
TI - Inf–sup stability of isogeometric Taylor–Hood and Sub-Grid methods for the Stokes problem with hierarchical splines
JF - IMA Journal of Numerical Analysis
DO - 10.1093/imanum/drx031
DA - 2018-04-01
UR - https://www.deepdyve.com/lp/oxford-university-press/inf-sup-stability-of-isogeometric-taylor-hood-and-sub-grid-methods-for-0duM6W5gJe
SP - 955
EP - 975
VL - 38
IS - 2
DP - DeepDyve
ER -