BPX preconditioners for isogeometric analysis using analysis-suitable T-splines

TL;DR

This paper introduces and analyzes two optimal multilevel BPX preconditioners for isogeometric analysis on locally refined T-meshes, demonstrating their efficiency and theoretical optimality through numerical experiments.

Contribution

The paper develops two novel BPX preconditioners tailored for analysis-suitable T-splines on locally refined meshes, with proven optimal complexity.

Findings

Both preconditioners achieve optimal complexity.

Numerical experiments confirm theoretical efficiency.

Preconditioners outperform existing methods in tests.

Abstract

We propose and analyze optimal additive multilevel solvers for isogeometric discretizations of scalar elliptic problems for locally refined T-meshes. Applying the refinement strategy in Morgenstern and Peterseim (2015, Comput. Aided Geom. Design, 34, 50-66) we can guarantee that the obtained T-meshes have a multilevel structure, and that the associated T-splines are analysis-suitable, for which we can define a dual basis and a stable projector. Taking advantage of the multilevel structure, we develop two BPX preconditioners: the first on the basis of local smoothing only for the functions affected by a newly added edge by bisection, and the second smoothing for all the functions affected after adding all the edges of the same level. We prove that both methods have optimal complexity, and present several numerical experiments to confirm our theoretical results, and also to compare the practical performance of the proposed preconditioners.

Tables

Table 1. Table 1: Condition numbers for the square domain. N.P: no preconditioner. Jac: Jacobi smoother. GS: Gauss-Seidel smoother

	Biquadratic				Bicubic				Biquartic
Levels	Dofs	N.P.	Jac.	GS	Dofs	N.P.	Jac.	GS	Dofs	N.P.	Jac.	GS
2	85	14.8	8.0	2.7	137	48.3	33.8	6.6	234	373.9	167.7	23.9
3	135	23.7	10.6	3.5	215	33.6	56.7	7.6	339	272.1	294.5	35.4
4	216	46.5	14.7	5.3	325	64.6	74.2	13.2	495	373.1	393.7	46.5
5	344	74.2	16.7	6.0	496	75.2	88.4	12.9	705	299.6	494.7	58.8
6	569	139.3	20.6	7.5	768	181.0	103.6	19.7	1047	378.4	572.2	67.0
7	961	234.5	20.4	7.6	1233	217.6	109.6	20.4	1575	313.8	642.8	82.6
8	1690	447.7	23.2	8.7	2045	540.6	123.6	24.9	2505	737.3	694.3	86.0
9	3042	798.9	22.1	8.5	3530	686.1	124.3	24.8	4101	918.9	743.1	103.6
10	5643	1560.9	24.6	9.4	6286	1785.7	133.2	27.9	7071	2286.9	770.1	108.7
11	10643	2907.2	23.1	9.1	11539	2371.1	133.9	27.4	12531	3030.4	805.1	118.1
12	20444	5770.1	25.5	9.8	21639	6382.3	140.8	29.9	23037	7868.2	816.5	121.1
13	39652	11040.6	24.0	9.4	41340	8729.9	140.4	29.4	43137	10855.7	843.2	127.1
14	77677	22109.3	26.2	9.9	79952	23983.8	145.9	31.2	82539	28931.9	846.6	128.1
15	152949	42979.9	24.6	9.6	156197	33390.6	145.0	30.6	159567	40886.2	867.2	132.7

Table 2. Table 2: Condition numbers for the curved L-shaped domain. N.P: no preconditioner. Jac: Jacobi smoother. GS: Gauss-Seidel smoother

	Biquadratic				Bicubic				Biquartic
Levels	Dofs	N.P.	Jac.	GS	Dofs	N.P.	Jac.	GS	Dofs	N.P.	Jac.	GS
2	275	60.7	10.2	2.8	436	122.6	47.1	6.7	735	1060.8	233.4	23.9
3	430	80.0	16.2	4.0	676	147.9	82.0	9.9	1057	744.6	394.5	36.5
4	682	202.4	22.6	5.6	1016	368.5	113.8	13.4	1537	1048.1	581.4	48.9
5	1074	270.5	27.5	6.6	1538	478.7	145.9	16.7	2177	870.5	756.6	61.6
6	1764	654.7	32.8	8.0	2370	1117.5	172.9	20.1	3221	1963.0	910.3	74.5
7	2954	914.5	36.3	8.5	3780	1509.3	198.4	22.7	4821	2506.6	1047.1	86.0
8	5168	2203.6	41.1	9.6	6244	3510.3	220.3	25.6	7641	5633.6	1168.9	97.9
9	9250	3238.5	43.6	9.7	10726	5013.1	240.2	27.4	12457	7710.7	1276.9	106.5
10	17104	7864.6	47.8	10.6	19046	11862.7	258.5	29.4	21421	17773.4	1375.3	116.2
11	32154	12027.0	49.8	10.5	34856	17793.4	275.0	30.4	37853	25944.3	1461.5	122.0
12	61656	29424.9	53.5	11.3	65256	42832.1	290.1	31.9	69473	61381.7	1541.2	129.2
13	119378	46178.5	54.9	11.1	124458	66460.2	303.5	32.4	129873	93799.5	1610.5	132.5

Table 3. Table 3: Condition number obtained with the micro technique, after adding all the lines of the same generation (between parentheses)

Biquadratic			Bicubic			Biquartic
Levels	Jac	GS	Levels	Jac	GS	Levels	Jac	GS
6 (2)	13.0	4.6	8 (2)	57.8	13.2	9 (2)	346.3	53.6
11 (3)	22.3	8.6	14 (3)	112.6	23.7	16 (3)	653.6	91.2
20 (4)	28.7	10.0	24 (4)	145.9	31.1	28 (4)	877.8	138.3
28 (5)	34.3	14.8	33 (5)	182.5	42.9	38 (5)	1097.8	166.9
43 (6)	41.4	14.8	49 (6)	198.5	45.9	56 (6)	1220.7	214.4
57 (7)	44.0	19.1	64 (7)	220.6	56.2	72 (7)	1379.3	234.5
84 (8)	48.4	18.6	92 (8)	232.6	56.6	102 (8)	1504.3	266.3
110 (9)	47.9	22.3	119 (9)	247.2	67.1	130 (9)	1604.3	285.0
161 (10)	51.8	20.6	171 (10)	252.2	63.5	184 (10)	1704.0	303.2
211 (11)	50.3	23.6	222 (11)	260.7	72.2	236 (11)	1758.4	316.8
310 (12)	54.0	21.4	322 (12)	264.1	66.3	338 (12)	1826.7	326.1

Table 4. Table 4: Condition numbers with the alternative refinement. N.P: no preconditioner. Jac: Jacobi smoother. GS: Gauss-Seidel smoother

	Biquadratic				Bicubic				Biquartic
Levels	Dofs	N.P.	Jac.	GS	Dofs	N.P.	Jac.	GS	Dofs	N.P.	Jac.	GS
2	184	23.4	6.4	2.4	213	33.5	30.0	5.9	244	243.3	134.6	19.8
3	569	103.5	10.2	3.6	620	126.6	48.8	9.0	673	390.8	285.8	34.8
4	1934	402.8	11.3	4.4	2027	608.9	60.5	12.0	2122	853.1	368.3	47.0
5	7059	1530.0	12.0	4.9	7234	2542.9	65.2	14.0	7411	4318.7	409.0	57.2
6	26904	5999.2	12.3	5.0	27241	10270.6	67.3	15.0	27580	18742.2	423.4	63.5
7	104989	23738.0	12.5	5.1	105648	41233.4	68.9	15.6	106309	77464.3	427.3	67.0

Equations

$$\text{SOLVE}⟶\text{ESTIMATE}⟶\text{MARK}⟶\text{REFINE},$$

$$-\Delta u=f\text{mbox}in\Omega ,u=0\text{mbox}on\partial \Omega ,$$

$$a(u,v)=(f,v)\forall v\in \mathcal{V}$$

$$a(u,v)={\int }_{\Omega }\nabla u\cdot \nabla vdx,(f,v)={\int }_{\Omega }fvdx,$$

$$(Au,v)=a(u,v),\forall u,v\in \mathcal{V}$$

$$Au=b$$

$$BAu=Bb.$$

$$\Vert u-u^k{\Vert }_A\le 2{\left(\frac{\sqrt{\kappa (BA)}-1}{\sqrt{\kappa (BA)}+1}\right)}^k\Vert u-u^0{\Vert }_A,$$

$$\mathcal{V}=\sum _{i=0}^J{\mathcal{V}}_i,$$

$$B:=\sum_{i=0}^{J}I_{i}R_{i}I^{t}_{i}\quad\big{(}=\sum_{i=0}^{J}R_{i}Q_{i}\big{)},$$

$$\sum _{i=0}^J\Vert v_i{\Vert }_A^2\le C_1\Vert v{\Vert }_A^2.$$

$$\left|\sum _{i=0}^J\sum _{j=i+1}^J(u_i,v_j{)}_A\right|\le C_2{\left(\sum _{i=0}^J\Vert u_i{\Vert }_A^2\right)}^{1/2}{\left(\sum _{i=0}^J\Vert v_i{\Vert }_A^2\right)}^{1/2}.$$

$$C_{4}^{-1}\|u_{i}\|_{A}^{2}\leq(R_{i}^{-1}u_{i},u_{i})\leq C_{3}\|u_{i}\|_{A}^{2}{\color[rgb]{0,0,0}\quad\text{ for any }u_{i}\in{\cal V}_{i}.}$$

$$\kappa (BA)\le (1+2C_2)C_1{C}_3{C}_4.$$

$$0={\xi }_0=\cdots ={\xi }_p<{\xi }_{p+1}\le \cdots \le {\xi }_{n-1}<{\xi }_n=\cdots ={\xi }_{n+p}=1,$$

$$B_{i,p}(\zeta )=B[{\Xi }_{i,p}](\zeta ),\zeta \in (0,1).$$

$$\text{mbox}supp(B_{i,p})=[{\xi }_i,{\xi }_{i+p+1}].$$

$${\tilde{I}}_j:=({\xi }_{j-p},{\xi }_{j+p+1}).$$

$${\mathcal{S}}_p(\Xi ):=\text{mbox}span\{B_{i,p},i=0,\dots ,n-1\}.$$

$${\lambda }_{j,p}(B_{i,p})={\delta }_{ij},0\le i,j\le n-1,$$

$$|{\lambda }_{j,p}(f)|\lesssim |{\xi }_{j+p+1}-{\xi }_j{|}^{-1/q}\Vert f{\Vert }_{L^q({\xi }_j,{\xi }_{j+p+1})}.$$

$${\lambda }_{i,p}(f)=\lambda [{\Xi }_{i,p}](f).$$

$${\Pi }_{p,\Xi }:L^2([0,1])\to {\mathcal{S}}_p(\Xi ),{\Pi }_{p,\Xi }(f)=\sum _{j=0}^{n-1}{\lambda }_{j,p}(f)B_{j,p}.$$

$$\Vert {\Pi }_{p,\Xi }(f){\Vert }_{L^2(I_k)}\le C\Vert f{\Vert }_{L^2({\tilde{I}}_k)},$$

$$|{\Pi }_{p,\Xi }(f){|}_{H^1(I_k)}\le C|f{|}_{H^1({\tilde{I}}_k)},$$

$${\mathbf{\Xi }}_{\mathbf{i},\mathbf{p}}=\{{\Xi }_{i_1,p_1},{\Xi }_{i_2,p_2}\}.$$

$$\left\{B_{\mathbf{i},\mathbf{p}}(\mathbf{\zeta })=B[{\Xi }_{i_1,p_1}]({\zeta }_1)B[{\Xi }_{i_2,p_2}]({\zeta }_2),\text{mbox}forall\mathbf{i}\in \mathbf{I}\right\}.$$

$$S_{\mathbf{p}}(\mathbf{\Xi })=\text{mbox}span\{B_{\mathbf{i},\mathbf{p}}(\mathbf{\zeta }),\mathbf{i}\in \mathbf{I}\}.$$

$${\mathcal{M}}^{\mathcal{B}}=\{Q_{\mathbf{j}}=I_{1,j_1}\times I_{2,j_2},\text{mbox}for0\le j_d\le N_d-1\}.$$

$${\tilde{Q}}_{\mathbf{j}}={\tilde{I}}_{1,j_1}\times {\tilde{I}}_{2,j_2},$$

Full text

1 Introduction

The analysis and development of adaptive schemes is one of the most active areas of research in the context of isogeometric analysis (IGA), a recent methodology for the solution of partial differential equations with high continuity splines. The main idea of adaptive methods is to obtain a good accuracy of the solution with less computational effort by applying local mesh refinement, hence adaptive IGA schemes require the use of different spline spaces that break the tensor-product structure of B-splines. The most popular alternatives in the IGA research community are T-splines [37], LR-splines [18, 28] and hierarchical splines [39].

In particular, in this paper we focus on T-splines, introduced by Sederberg et al. in [37] and applied in IGA for the first time in [2, 21]. T-splines are constructed from a T-mesh, a rectangular tiling with hanging nodes, and T-spline blending functions are defined from their local knot vectors, which are computed from the tiling. The mathematical research on T-splines has been very active in recent years, and it has led to the introduction of analysis-suitable (or dual compatible) T-splines [31, 4, 5, 30], a subset of T-splines with good approximation properties that provide local linear independence.

A standard adaptive scheme based on mesh refinement can be written in a loop of the form [14]

[TABLE]

and suitable strategies for all the steps of the adaptive scheme are needed in order to guarantee its efficiency. In this paper we focus on the solution of the linear system arising in the SOLVE module for T-splines, and study the optimality of a suitable BPX preconditioner. Several domain decomposition preconditioners have been recently studied in the IGA context: overlapping Schwarz methods [8], balancing domain decomposition by constraint (BDDC) methods [7, 3] and dual primal isogeometric tearing and interconnecting methods [29, 34]. Multilevel preconditioners for IGA have been extensively analyzed in the tensor-product setting: the BPX preconditioner in [13], and multigrid preconditioners [22, 25, 26, 24, 23, 20]. Moreover, preconditioners based on fast solvers for Sylvester-like equations have been proposed in the recent paper [35].

The T-splines obtained with the refinement strategy in [33] are analysis-suitable by construction, and present a multilevel structure, which makes them very appealing to apply multilevel preconditioners. In the present paper, we first present theoretically optimal multilevel preconditioners for IGA on locally refined T-meshes, extending the results in [13] to T-splines.

For the study of the optimality of the BPX preconditioner we follow [14], writing the preconditioner in the framework of the parallel subspace correction (PSC) method. In this framework, the optimality follows from two basic properties: a stable space decomposition, and the strengthened Cauchy-Schwarz inequality. The proof of these two properties is the core of this work, and as a consequence we obtain that the BPX preconditioner gives a uniformly bounded condition number, which is independent of the mesh size $h$, but depends on the degree $p$ as in [13].

The construction of the BPX preconditioner as in [14] is performed by adding a new edge to the T-mesh by bisection, and the new level is defined by the functions appearing or modified by the insertion of this edge. An alternative construction is also proposed, adding all the edges of the same generation at once, and defining the functions of the new level as those appearing or modified after the insertion of all edges. The theoretical optimality for this alternative construction, that we name macro decomposition, follows from the previous one, but the numerical results show an improved performance.

The paper is organized as follows. In Section 2, we introduce the framework of parallel subspace correction (PSC) method and present its convergence theory based on the two properties mentioned above. We briefly review the basics of univariate/multivariate B-splines in Section 3. In Section 4, we give a new definition of T-meshes by bisections and discuss ${\bf p}$-admissible T-meshes and their fundamental properties as in [33]. In Section 5, we construct a space decomposition on ${\bf p}$-admissible T-meshes and then prove that the two aforementioned properties are satisfied. In Section 6 we obtain the optimality result for the BPX preconditioner, and also deal with the macro decomposition, showing that is also an efficient space decomposition for the purpose of implementation. Some numerical results that validate our theory are presented in Section 7.

2 Preliminaries

2.1 Problem setting

We are interested in the second order Laplacian with homogeneous Dirichlet boundary conditions,

[TABLE]

where $\partial\Omega$ denotes the boundary of $\Omega$ and $f\in L^{2}(\Omega)$. The isogeometric approximation to the solution of (1) is the function $u\in\mathcal{V}$ such that

[TABLE]

where

[TABLE]

and ${\cal V}$ is the isogeometric discrete space. Defining a linear operator $A:\mathcal{V}\rightarrow\mathcal{V}$ by

[TABLE]

and also $b\in\mathcal{V}$ by $(b,v)=(f,v),\ \forall v\in\mathcal{V}$, we have to solve the linear operator equation

[TABLE]

for some $u\in\mathcal{V}$.

In the rest of the paper, we will adopt the following compact notation. Given two real numbers $a,b$ we write $a\lesssim b$, when $a\leq Cb$ for a generic constant $C$ independent of the knot vectors and the mesh size $h$ but depending on the spline degree $p$ and the geometric map ${\mathbf{F}}$ (defined below), and we write $a\approx b$ when $a\lesssim b$ and $b\lesssim a$.

2.2 Preconditioned conjugate gradient method

Let $B$ be a symmetric positive definite (SPD) operator. Applying it to both sides of $Au=b$, we get an equivalent equation

[TABLE]

The preconditioned conjugate gradient (PCG) method can be viewed as a conjugate gradient method applied to solving $BAu=Bb$ where $B$ is called a preconditioner (see, e.g., [14, 41] for an extensive description).

Let $u^{k},\ k=0,1,2,\ldots$, be the solution sequence of the PCG algorithm. Then the following error estimate is well-known:

[TABLE]

which implies that the PCG method converges faster with a smaller condition number $\kappa(BA)$.

2.3 The method of parallel subspace corrections

The method of parallel subspace corrections (PSC) provides a particular construction of the iteration operator $B$. The starting point is a suitable decomposition of $\mathcal{V}$:

[TABLE]

where $\mathcal{V}_{i}$ are subspaces of $\mathcal{V}$. The model problem (1) can be split into sub-problems in each $\mathcal{V}_{i}$ with smaller size. Throughout this paper, we use the following operators, for $i=0,1,\ldots,J$:

• •

$Q_{i}:\mathcal{V}\rightarrow\mathcal{V}_{i}$ the projection in the inner product $(\cdot,\cdot):=(\cdot,\cdot)_{L^{2}(\Omega)}$;

• •

$I_{i}:\mathcal{V}_{i}\rightarrow\mathcal{V}$ the natural inclusion;

• •

$P_{i}:\mathcal{V}\rightarrow\mathcal{V}_{i}$ the projection in the inner product $(\cdot,\cdot)_{A}:=(A\cdot,\cdot)$;

• •

$A_{i}:\mathcal{V}_{i}\rightarrow\mathcal{V}_{i}$ the restriction of $A$ to the subspace $\mathcal{V}_{i}$;

• •

$R_{i}:\mathcal{V}_{i}\rightarrow\mathcal{V}_{i}$ an approximation of $A_{i}^{-1}$.

This method performs the correction on each subspace in parallel, with the operator $B$ defined by

[TABLE]

where $I_{i}^{t}$ denote the adjoint of $I_{i}$ with respect to $(\cdot,\cdot)$. It is readily to check that $Q_{i}A=A_{i}P_{i}$ and $Q_{i}=I_{i}^{t}$ with $(I^{t}_{i}u,v_{i})=(u,I_{i}v_{i})$.

The convergence analysis of (PSC) is based on the following two important properties [14]:

(A1) Stable Decomposition For any $v\in\mathcal{V}$, there exists a decomposition $v=\sum_{i=0}^{J}v_{i},\ v_{i}\in\mathcal{V}_{i},\ i=0,\ldots,J$ such that

[TABLE]

(A2) Strengthened Cauchy-Schwarz (SCS) inequality For any $u_{i},v_{i}\in\mathcal{V}_{i},i=0,\ldots,J$

[TABLE]

For a space decomposition satisfying both properties, one can prove the following result on the preconditioned linear system [14]:

Theorem 2.1

*Let $\mathcal{V}=\sum_{i=0}^{J}\mathcal{V}_{i}$ be a space decomposition satisfying ***(A1) **and (A2), and let $R_{i}$ be SPDs for $i=0,\ldots,J$ such that

[TABLE]

Then $B$ defined by (2) is SPD and

[TABLE]

The goal of the paper is the construction of a preconditioner like (2) for T-splines, and the proof that **(A1) **and (A2) are satisfied.

3 Splines

In this section we recall the definition and main properties of B-splines, mainly to fix the notation. For a more extensive description on the use of splines in isogeometric analysis, the reader is referred to [27, 15], and to [6] for a mathematical perspective.

3.1 Univariate B-splines

Given two positive integers $p$ and $n$, we say that $\Xi:=\{\xi_{0},\ldots,\xi_{n+p}\}$ is a $p$-open knot vector if

[TABLE]

where repeated knots are allowed. From the knot vector $\Xi$, univariate B-spline basis functions of degree $p$ are defined recursively using the Cox-De Boor formula (see [16]). The definition of each B-spline ${B}_{i,p},\ i=0,\ldots,n-1$, is determined by a $p+2$ local knot vector $\Xi_{i,p}=\{\xi_{i},\ldots,\xi_{i+p+1}\}$. Whenever necessary, we will stress it by adopting the following notation:

[TABLE]

Thus, the basis function ${B}_{i,p}$ has support

[TABLE]

Throughout the paper, we assume that the maximum multiplicity of the internal knots is less than or equal to the degree $p$, that is, the B-spline functions are at least continuous.

Let us select from $\Xi$ a subset $\{\xi_{i_{k}},k=0,\ldots,N\}$ of non-repeated knots, or breakpoints, with $\xi_{i_{0}}=0$, $\xi_{i_{N}}=1$. We point out that the local mesh size of the element $I_{k}=(\xi_{i_{k}},\xi_{i_{k+1}})$ is called $h_{k}=\xi_{i_{k+1}}-\xi_{i_{k}}$ for $k=0,\ldots,N-1$. Moreover, to the element $I_{k}=(\xi_{i_{k}},\xi_{i_{k+1}})$, that can be written as $(\xi_{j},\xi_{j+1})$ for a certain $j$, we associate the support extension $\widetilde{I}_{j}$ defined by

[TABLE]

We denote by

[TABLE]

the univariate B-splines space spanned by those B-splines of degree $p$. The functions ${B}_{i,p}$ are a partition of unity, as shown in [36].

Following [36, Theorem 4.41], we define suitable functionals $\lambda_{j,p}=\lambda[\xi_{j},\ldots,\xi_{j+p+1}]$, for $0\leq j\leq n-1$, which are dual to the B-splines basis functions, that is

[TABLE]

where $\delta_{ij}$ is the Kronecker delta. The following estimate of the functionals $\lambda_{j,p}$ will be useful in the sequel.

Lemma 3.1

If $f\in L^{q}(\xi_{j},\xi_{j+p+1})$, with $1\leq q\leq+\infty$ then

[TABLE]

Proof. The proof can be found in [36, Theorem 4.41]. $\square$

We note that these dual functionals are locally defined and depend only on the corresponding local knot vector, namely,

[TABLE]

Let $\Pi_{p,\Xi}$ be the projection that is built from the dual basis as detailed in [36, Theorem 12.6], that is,

[TABLE]

Assumption 3.1

$\{\xi_{i_{0}},\xi_{i_{1}},\ldots,\xi_{i_{N}}\}$* is locally quasi-uniform, that is, there exists a constant $\theta\geq 1$ such that the mesh sizes $h_{k}=\xi_{i_{k+1}}-\xi_{i_{k}}$ satisfy the relation $\theta^{-1}\leq h_{k}/h_{k+1}\leq\theta$, for $k=0,\ldots,N-2$.*

Proposition 3.2

For any non-empty knot span $I_{k}=(\xi_{i_{k}},\xi_{i_{k+1}})$, we have

[TABLE]

where the constant $C$ depends only on the degree $p$. Moreover, if Assumption 3.1 holds, we have

[TABLE]

where the constant $C$ depends only on $p$ and $\theta$ and $H^{1}$ denotes the classical Sobolev norm.

Proof. The proof can be found in [6, Proposition 2.2]. $\square$

3.2 Multivariate splines

Multivariate B-splines can be constructed by means of tensor products. We discuss here the bivariate case, the higher-dimensional case being analogous.

For $d=1,2$, assume that $n_{d}\in\mathbb{N}$, the degree $p_{d}$ and the $p_{d}$-open knot vector $\Xi_{d}=\{\xi_{d,0},\ldots,\xi_{d,n_{d}+p_{d}}\}$ are given. We set the polynomial degree vector ${\mathbf{p}}=(p_{1},p_{2})$ and ${\bf\Xi}=\{\Xi_{1},\Xi_{2}\}$. We introduce a set of multi-indices ${\bf I}=\{{\bf i}=(i_{1},i_{2}):0\leq i_{d}\leq n_{d}-1\}$ and for each multi-index ${\bf i}=(i_{1},i_{2})$, we define the local knot vector

[TABLE]

Then we can introduce the set of multivariate B-splines

[TABLE]

The spline space in the parametric domain ${\Omega}=(0,1)^{2}$ is then

[TABLE]

We also introduce the set of non-repeated interface knots for each direction $\{\xi_{d,i_{0}},\ldots,\xi_{d,i_{N_{d}}}\}$, $d=1,2$, which determine the intervals $I_{d,j_{d}}=(\xi_{d,i_{j_{d}}},\xi_{d,i_{j_{d}+1}})$, for $0\leq j_{d}\leq N_{d}-1$. These intervals lead to the Cartesian grid ${\mathcal{M^{B}}}$ (or simply ${\cal M}$) in the unit domain ${\Omega}=(0,1)^{2}$, also called the Bézier mesh:

[TABLE]

For a generic element $Q_{\bf j}$, we also define its support extension

[TABLE]

with $\widetilde{I}_{d,j_{d}}$ the univariate support extension by (5).

3.2.1 Multivariate quasi-interpolants

We note that when the univariate quasi-interpolants are defined from a dual basis, as in (6), then the multivariate quasi-interpolant is also defined from a dual basis. Indeed, we have

[TABLE]

where each dual functional is defined from the univariate dual bases by the expression

[TABLE]

4 T-splines

The main drawback of B-splines is their tensor-product structure, which prevents local refinement as required by adaptive methods. One of the alternatives is the use of T-splines [37], a superset of B-splines that allows for local refinement. The use of T-splines in IGA was first explored in [2, 21], and it has led to a growing interest for the analysis of their mathematical properties. In this section we are collecting mathematical results from [31, 4, 5] (linear-independence, dual basis and projectors), [30, 12] (nestedness and space characterization) and [6, Section 7] (local linear independence), following mainly the notation in [6].

We restrict ourselves to T-splines where refinement is always performed by bisection, and the multiplicity is never reduced. A more general setting can be considered, but it would only add technical difficulties without adding more insight.

4.1 T-mesh defined by bisection

As in the previous section, let us assume that we are given the degrees $p_{d}$, the integers $n_{d}$ and the open knot vectors $\Xi_{d}$, and let us denote $m_{d}=n_{d}+p_{d}$ for $d=1,2$. For simplicity we assume that the internal knots in $\Xi_{d}$ are not repeated and equally spaced, so the element size in each parametric direction can be denoted by $h_{1}$ and $h_{2}$. Our starting point is the index mesh ${\cal T}_{0}$ in the index domain $[0,m_{1}]\times[0,m_{2}]$, defined as the Cartesian grid of unit squares

[TABLE]

which is associated to the tensor-product B-spline space $S_{\bf p}(\boldsymbol{\Xi})$.

To define T-splines by bisection, we start introducing, for any integer $\ell\geq 0$ and for $d=1,2$, the set of rational indices

[TABLE]

and we notice that ${\cal I}^{\ell}_{d}\subset{\cal I}^{\ell^{\prime}}_{d}$ if $\ell\leq\ell^{\prime}$. We also define the ordered knot vectors

[TABLE]

in a recursive way: starting from $\Xi_{d}^{0}=\Xi_{d}$, for $\ell>0$ and for any new index $k\in{\cal I}^{\ell}_{d}\setminus{\cal I}^{\ell-1}_{d}$, we define the knot

[TABLE]

which is well defined because $k-\frac{1}{2^{\ell}},k+\frac{1}{2^{\ell}}\in{\cal I}_{d}^{\ell-1}$. Clearly, $\Xi_{d}^{\ell}\subset\Xi_{d}^{\ell+1}$ for $\ell\geq 0$, and the interval size is $h_{d,\ell}=h_{d}/2^{\ell}$. Notice that in this procedure we do not introduce new knots between the repeated knots of the open knot vector.

We also define, for an arbitrary rectangular element in the index domain $\tau=[k_{1},k_{1}+t_{1}]\times[k_{2},k_{2}+t_{2}]$, with indices $k_{d},k_{d}+t_{d}\in{\cal I}^{\ell}_{d}$, the following bisection operators (see [33, Definition 2.5]):

[TABLE]

Notice that the bisection operators split the element in two adding a vertical and a horizontal edge, respectively, only when the corresponding length in the parametric domain is greater than zero, that is, when $\xi_{d,k_{d}+t_{d}}-\xi_{d,k_{d}}>0$.

Starting from the Cartesian grid ${\cal T}_{0}$, we define a T-mesh ${\cal T}={\cal T}_{N}$ by successive applying bisection of elements, in the form

[TABLE]

where we use the formal addition (see [33, Definition 2.6] and also [14, Section 3])

[TABLE]

and the bisection operator can be either ${\rm bisect}_{x}$ or ${\rm bisect}_{y}$. Moreover, we define the finest level $\Lambda$ as the minimum integer such that the $d$th coordinate of all vertices in the T-mesh ${\cal T}$ belongs to ${\cal I}_{d}^{\Lambda}$, for $d=1,2$.

Remark 4.1

We notice that a T-mesh defined with this procedure automatically satisfies the first condition in the definition of admissible mesh in [6, Definition 7.10], that is, the first $p_{d}+1$ lines closer to each boundary are completely contained in the mesh. Instead, the second condition of not having T-junctions in the so-called frame region is not satisfied, because T-junctions may appear on the interface between the frame and the active region. In any case, these T-junctions do not affect the definition of the T-spline functions, and the mesh can still be considered “admissible”.

4.2 Analysis suitable T-splines

To construct the blending functions associated to a T-mesh ${\cal T}$, we have to define first the set of anchors, that we denote by ${\cal A}_{\bf p}({\cal T})$, and that depends on the degree, see [6, Definition 7.13]. These are either the set of vertices ($p_{1}$ and $p_{2}$ odd), elements ($p_{1}$ and $p_{2}$ even), horizontal edges ($p_{1}$ even, $p_{2}$ odd) or vertical edges ($p_{1}$ odd, $p_{2}$ even) in the active region, which is the rectangle

[TABLE]

see the examples in Figure 1. We also define the horizontal (resp. vertical) skeleton of the mesh, and denote it by ${\rm hSk}({\cal T})$ (resp. ${\rm vSk}({\cal T})$), as the union of all horizontal (resp. vertical) edges and vertices. The union of ${\rm hSk}({\cal T})$ and ${\rm vSk}({\cal T})$ will be called skeleton. Then, for each anchor we construct an ordered horizontal index vector of $p_{1}+2$ indices, tracing a horizontal line from the anchor and collecting the closest $\lfloor(p_{1}+2)/2\rfloor$ indices leftwards and rightwards where the line intersects the vertical skeleton of the mesh, plus the index of the anchor if the degree is odd, see [6, Definition 7.14] and Figure 1. A vertical index vector of $p_{2}+2$ indices is constructed in an analogous way tracing a vertical line passing through the anchor.

Then, given an anchor ${\bf A}\in{\cal A}_{\bf p}({\cal T})$ with horizontal and vertical index vectors

[TABLE]

we define the local knot vectors

[TABLE]

and from these local knot vectors we define the associated bivariate function

[TABLE]

We will denote by $S({\cal A}_{\bf p}({\cal T}))={\rm span}\{B_{{\bf A},{\bf p}}:{\bf A}\in{\cal A}_{\bf p}({\cal T})\}$ the space generated by the T-splines.

A sufficient condition to guarantee linear independence of the T-spline functions is the dual compatibility condition, introduced in [4, 5] and equivalent to the analysis-suitability condition in [30]. The former was generalized to arbitrary dimension in [6]. We say that two local knot vectors overlap if both of them can be written as subvectors, with consecutive indices, of the same global knot vector. Then we say that the T-mesh is dual compatible if for each pair of anchors ${\bf A}^{\prime},{\bf A}^{\prime\prime}\in{\cal A}_{\bf p}({\cal T})$, with ${\bf A}^{\prime}\neq{\bf A}^{\prime\prime}$, there exists a direction $d\in\{1,2\}$ such that the local knot vectors $\Xi_{{\bf A}^{\prime},d,p_{d}}$ and $\Xi_{{\bf A}^{\prime\prime},d,p_{d}}$ are different and overlap. See [6, Section 7.1] for details.

Assuming that the T-mesh is dual compatible, then the functionals

[TABLE]

form a dual basis for $\{{B}_{\bf A,p}:{\bf A}\in\mathcal{A}_{\bf p}(\mathcal{T})\}$ (see [6, Proposition 7.3]). Moreover, we can build the projection operator ${\bf\Pi}_{\bf p}^{\cal T}\equiv{\bf\Pi}_{\mathcal{A}_{\bf p}(\mathcal{T})}:L^{2}({\Omega})\rightarrow S(\mathcal{A}_{\bf p}(\mathcal{T}))$ by

[TABLE]

For the analysis of the projector properties, we make use of the Bézier mesh, that we define as in [6, Section 7.3]. The Bézier mesh is different from the T-mesh, and plays a similar role to the mesh in finite elements. We start recalling that T-junctions are internal vertices of the T-mesh with valence equal to three, that can be grouped in horizontal ($\vdash,\dashv$) and vertical ($\top,\bot$) T-junctions. For a T-junction of type $\dashv$ with index coordinates $(\bar{\imath},\bar{\jmath})$, we define the extension as the minimal horizontal line that, passing through the T-junction, intersects $\lfloor p_{1}/2\rfloor$ (closed) vertical edges to its left, and $\lceil p_{1}/2\rceil$ to its right. The extensions in the other cases are defined in a similar way, using $p_{2}$ for vertical T-junctions. Then, we define ${\rm ext}_{\bf p}({\cal T})$ the extended T-mesh in the index domain adding to ${\cal T}$ all the T-junction extensions. The Bézier mesh ${\cal M^{T}}$ (or simply ${\cal M}$) associated to ${\cal T}$ is then defined as the collection of non-empty elements in the domain $\Omega$ of the form

[TABLE]

as depicted in Figure 2.

For a Bézier element $Q\in{\cal M^{T}}$, and in general for any subdomain $Q\subset\Omega$, we define the support extension as the union of the supports of the functions whose support intersects $Q$, that is

[TABLE]

where int($C$) denotes the interior of a set $C$. Moreover, we define $\bar{Q}$ as the smallest rectangle in $\Omega$ containing $\widetilde{Q}$. The following result holds (see [6, Proposition 7.7]):

Proposition 4.1

Let $\mathcal{T}$ be a dual compatible T-mesh. Then there exists a constant $C$, depending only on ${\bf p}$, such that for any Bézier element $Q\in{\cal M^{T}}$ the projector (11) satisfies

[TABLE]

Finally, we notice that for each anchor the index vectors (10) define a local Cartesian grid of $(p_{1}+1)(p_{2}+1)$ cells, called tiled floor in [4, 5]. Moreover, we also define the parametric tiled floor, as the set of non-empty cells

[TABLE]

We remark that in general, the cells of the tiled floor do not coincide with the elements of the T-mesh, and the cells of the parametric tiled floor do not coincide with the elements of the Bézier mesh.

4.3 Analysis suitable T-splines by bisection

In order to apply BPX preconditioners to analysis-suitable T-splines, it is necessary to define a suitable refinement procedure that provides a multilevel structure. In this section we adopt the refinement strategy introduced and analysed in [33], and present a new local quasi-uniformity result necessary for the analysis of multilevel preconditioners. The idea in that work is to refine by bisection alternating the refinement direction, and whenever a new edge is added, a recursive algorithm is called to refine the elements in the neighborhood, ensuring that the condition of analysis-suitability is preserved. One of the advantages of this refinement algorithm is that it allows to associate a level (or generation) to each element and function, as required by multilevel methods.

We start setting the generation $g(\tau)=0$ for all the elements of the Cartesian grid $\tau\in{\cal T}_{0}$. Then, the T-mesh ${\cal T}$ is defined as in (8)–(9), choosing the bisection

[TABLE]

and setting the generation $g(\tau^{\prime})=g(\tau)+1$ for $\tau^{\prime}\in{\rm bisect}(\tau)$ (see [33, Definition 2.6]). Moreover, we say that the bisection $b_{\tau}$ has generation $g(b_{\tau})=g(\tau)+1$. Without loss of generality (see also the proof of Theorem 3.6 in [33]), in the following we will assume that the bisections in (8) are ordered by their generation, that is, if $k>k^{\prime}$ then $g(b_{\tau_{k}})\geq g(b_{\tau_{k^{\prime}}})$, and the same relation holds for the elements generated by the two bisections. In the following, we will denote by $L$ the finest generation (or level), that is, the generation of the last bisection.

Remark 4.2

We recall that elements with zero length in one parametric direction are not bisected in that direction, see Section 4.1. However, when applying the bisection operator their generation is increased by one, in such a way that the next time the bisection operator is applied, they will be refined in the other direction.

To define the generation of the functions, we denote by $\Phi_{0}=\{B_{{\bf A},{\bf p}}:{\bf A}\in{\cal A_{\bf p}}({\cal T}_{0})\}$ the functions of the tensor product spline space in the coarsest mesh, and then define for $k=1,\ldots,N$ the collection of T-spline functions newly appeared or modified after the bisection $b_{\tau_{k-1}}$ as

[TABLE]

An example of the definition of $\Phi_{k}$ is shown in Figure 3. Functions in $\Phi_{0}$ have generation 0, whereas functions in $\Phi_{k}$ have the same generation of the bisection $b_{\tau_{k-1}}$. To alleviate notation, in the following we will denote by $\ell_{k}$ the generation (or level) of functions in $\Phi_{k}$, that is

[TABLE]

with the convention that $g(\tau_{-1})=-1$. Notice that the subscript $k$ varies from [math] to $N$, while $\ell_{k}$ takes values from 0 to $L$.

To each $\Phi_{k}$, we associate a subspace

[TABLE]

We also define the support of functions in $\Phi_{k},1\leq k\leq N$ and its support extension as

[TABLE]

and

[TABLE]

respectively.

Notice that, since we alternate the directions of refinement, and recalling the notation of Section 4.1, the local knot vectors of a function of generation $\ell$ are contained in

[TABLE]

For convenience, we also introduce the corresponding set of rational indices, ${\cal\bf I}^{\ell}=\{{\cal I}_{1}^{\lceil\ell/2\rceil},{\cal I}_{2}^{\lfloor\ell/2\rfloor}\}$, and the Bézier mesh in the parametric domain ${\cal M}_{\ell}^{\cal B}$. For each generation $\ell$ we have that the mesh size is $h_{\ell}\approx 2^{-\ell/2}$. This important relation between generation and mesh size can be also represented using additional notation as

[TABLE]

Apart from the generations, that are necessary to provide the multilevel structure, we also need some definitions to adapt the refinement algorithm from [33] to the case of having an open knot vector. Given two points ${\bf x},{\bf x^{\prime}}\in\mathbb{R}^{2}$ we define the distance between them componentwise as the vector

[TABLE]

Moreover, let us define, for a point in the index domain ${\bf x}=(x_{1},x_{2})\in[0,m_{1}]\times[0,m_{2}]$, its translated version ${\bf\widetilde{x}}=(\widetilde{x}_{1},\widetilde{x}_{2})$ given by

[TABLE]

Given a bisection T-mesh ${\cal T}$ and an element $\tau\in{\cal T}$, we denote its middle point as ${\bf x}_{\tau}=(x_{1,\tau},x_{2,\tau})$. We can define the distance of a point ${\bf x}$ to the element $\tau$, and the distance between two elements $\tau$ and $\tau^{\prime}$ as

[TABLE]

respectively. Then, we define the ${\bf p}$-neighborhood of $\tau$ as (see [33, Definition 2.4] and Figure 4)

[TABLE]

where

[TABLE]

We also need to define the set (see [33, Corollary 2.15])

[TABLE]

Given a bisection T-mesh ${\cal T}$ and an element $\tau\in{\cal T}$, we say that the bisection of $\tau$ is ${\bf p}$-admissible or simply admissible, if all $\tau^{\prime}\in{\cal G}_{\bf p}(\tau)$ satisfy $g(\tau^{\prime})\geq g(\tau)$. Moreover, we say that a bisection T-mesh ${\cal T}$ is ${\bf p}$-admissible, if it can be obtained as in (8) with a sequence of admissible bisections, see the example of Figure 4. It has been proved in [33, Theorem 3.6] that admissible T-meshes are also analysis suitable, and therefore the dual basis and the projector of the previous section can be built.

Remark 4.3

The result in [33] does not take into account the repeated knots of the open knot vector. However, the same ideas apply using the definitions above, because the bisection of zero measure elements only adds new lines in the “safe” direction, without causing intersection of T-junction extensions. The use of the translated points in practice forces that a line arriving at a repeated knot, which is on the boundary of the parametric domain, should continue until the boundary of the index domain.

Besides the analysis suitability condition, we need a local quasi-uniformity result, for which it is necessary to use the following auxiliary lemmas.

Lemma 4.2

Let ${\bf A}\in{\cal A}_{\bf p}({\cal T})$ be an anchor of a ${\bf p}$-admissible T-mesh ${\cal T}$. Then, any cell in the parametric tiled floor of ${\bf A}$ contains at most two Bézier elements of ${\cal M^{T}}$.

Proof. Given a cell of the tiled floor, from [5, Lemma 3.2] it does not contain any vertex of ${\cal T}$ in its interior, and any line of the extended mesh in its interior belongs to a T-junction extension. Since the mesh is analysis suitable, only vertical or horizontal lines can be found in its interior, but not both. Since a cell in the parametric tiled floor corresponds to a cell in the tiled floor, the result holds because we are refining by bisection alternating the refinement directions. $\square$

Lemma 4.3

Given a ${\bf p}$-admissible T-mesh ${\cal T}$ and an element $\tau\in{\cal T}$, for any $\tau^{\prime}\in{\cal G}_{\bf p}(\tau)$ it holds $g(\tau^{\prime})\geq g(\tau)-1$.

Proof. See [33, Lemma 2.14]

Lemma 4.4

Let ${\bf A}\in{\cal A}_{\bf p}({\cal T})$ be an anchor of a ${\bf p}$-admissible T-mesh ${\cal T}$, associated to a function of generation $\ell$. Then the length of the cells of the parametric tiled floor of ${\bf A}$ in each parametric direction is equal to

[TABLE]

Proof. The lemma is proved by induction, similarly to [33, Theorem 3.6]. The result clearly holds for the functions in the tensor-product space associated to the T-mesh ${\cal T}_{0}$. Assuming it is true for ${\cal T}_{k}$, we have to prove it for ${\cal T}_{k+1}={\cal T}_{k}+b_{\tau_{k}}$. Since the bisections can be ordered by their generation, it holds that $\max_{\tau\in{\cal T}_{k}}g(\tau)=g(\tau_{k})+1$. We assume that $g(\tau_{k})$ is even, the odd case is proved analogously exchanging the role of the vertical and horizontal directions.

Introducing ${\cal T}^{u}_{\ell}$ the Cartesian grid obtained after applying all possible bisections of generation $\ell$, using the same arguments as in [33, Theorem 3.6], we have the following result about the edges of the T-mesh in the $U_{\bf p}(\tau_{k})$ region

[TABLE]

We have to check that the result is true for the functions introduced or modified by the bisection $b_{\tau_{k}}$, that is, for functions in $\Phi_{k+1}$ of generation $\ell_{k+1}=g(\tau_{k})+1$. Given one of these functions, from the length of the index vectors and (16)–(17) its associated anchor must intersect $U_{\bf p}(\tau_{k})$. Moreover, using [33, Corollary 2.15] this anchor is associated to a geometrical entity (vertex, edge or element) that belongs to an element $\tau^{\prime}\in{\cal G}_{\bf p}(\tau_{k})$. Using first Lemma 4.3 and then the definition of ${\bf p}$-admissible T-mesh, for any element $\tau^{\prime\prime}\in{\cal G}_{\bf p}(\tau^{\prime})$ it holds $g(\tau^{\prime\prime})\geq g(\tau^{\prime})-1\geq g(\tau_{k})-1$.

Since the horizontal index vector is collected from the vertical lines of elements in ${\cal G}_{\bf p}(\tau^{\prime})$, in the horizontal direction the intervals of the index vector have lengths at most $2^{-g(\tau_{k})/2}$ and $2^{-(g(\tau_{k})/2+1)}$, and recalling that $\ell_{k+1}=g(\tau_{k})+1$ is odd, we have the result for the horizontal direction.

For the vertical index vector, following an analogous reasoning we obtain that the intervals in the vertical direction also have lengths $2^{-g(\tau_{k})/2}$ and $2^{-(g(\tau_{k})/2+1)}$. Moreover, since the anchor of any function in ${\Phi}_{k+1}$ is horizontally aligned with the new edge, using (16) the length of the vertical intervals is always $2^{-g(\tau_{k})/2}=2^{-(\ell_{k+1}-1)/2}$.

Finally, the result is proved passing from the index vectors to the local knot vectors, using the length of the intervals in the knot vectors $\Xi_{d}^{\ell}$ from Section 4.1. $\square$

Corollary 4.5

In a ${\bf p}$-admissible T-mesh, for any basis function ${B}_{{\bf A},{\bf p}}\in\Phi_{k}$, and for any Bézier element $Q\subset{\rm supp}(B_{{\bf A},{\bf p}})$, it holds that $h_{Q}\approx h_{\ell_{k}}\approx 2^{-\ell_{k}/2}$.

Proposition 4.6

*For a ${\bf p}$-admissible T-mesh ${\cal T}$, and for any Bézier element $Q\in\cal M^{T}$ of the associated Bézier mesh, it holds that $h_{Q}\approx h_{\widetilde{Q}}\approx h_{\bar{Q}}$.

Proof. We only need to prove that $h_{\bar{Q}}\lesssim h_{\widetilde{Q}}\lesssim h_{Q}$, because $Q\subset\widetilde{Q}\subset\bar{Q}$. We start proving $h_{\widetilde{Q}}\lesssim h_{Q}$. Since the T-mesh is analysis-suitable, in the element $Q$ there are at most $(p_{1}+1)(p_{2}+1)$ basis functions that do not vanish [6, Proposition 7.6], and for each of these functions the element $Q$ is contained in a cell of its parametric tiled floor. Moreover, from Lemma 4.2 this cell contains at most two elements of the Bézier mesh. Finally, each function contains at most $(p_{1}+1)(p_{2}+1)$ cells in its tiled floor, and the result of Lemma 4.4 states that the size of all these cells is comparable. Combining all these results we obtain $h_{\widetilde{Q}}\lesssim h_{Q}$ with a constant that depends on the degree ${\bf p}$.

Let us denote by $h_{1,Q}$ and $h_{2,Q}$ the length of $Q$ and in each parametric direction, and by $h_{1,\widetilde{Q}}$ and $h_{2,\widetilde{Q}}$ the length of $\widetilde{Q}$, given by the difference in the first (second) parametric direction between the rightmost (uppermost) and leftmost (downmost) points in $\widetilde{Q}$. From the previous results and the relative position of the Bézier element inside the tiled floor of each function (see the proof of [6, Proposition 7.6]), we have that

[TABLE]

and since $\bar{Q}$ is the minimum rectangle that contains $\widetilde{Q}$, it also holds that $h_{\bar{Q}}\lesssim h_{Q}$.

$\square$

Remark 4.4

For simplicity we have assumed that the internal knots in the knot vectors $\Xi_{d}^{\ell}$ are equally spaced. The results of this section, and in particular Lemma 4.4, can be extended to the case of local quasi-uniform knot vectors under Assumption 3.1.

Proposition 4.7

For a ${\bf p}$-admissible T-mesh defined as in (8), and for all $u\in L^{2}(\Omega)$, and $\ell\geq 1$, we have that

[TABLE]

where $\omega_{k}$ and $\widetilde{\omega}_{k}$ are as in (14) and (15), respectively.

Proof. The result follows from a careful counting of the number of times a Bézier element $Q$ (or its children, if it is bisected) can appear in the sets $\omega_{k}$ and $\widetilde{\omega}_{k}$ of the same generation $\ell_{k}$. Studying the worst case scenario, when all the elements in the neighboring of $Q$ are bisected, we see from the graphical explanation in Figure 5 that the element is contained in at most $(2p_{1}+1)(2\lceil p_{2}/2\rceil+1)$ (or $(2p_{2}+1)(2\lceil p_{1}/2\rceil+1)$) sets $\omega_{k}$.

A similar counting, that can be understood from Figure 6, gives that the element $Q$ is contained in at most $(4p_{1}+1)(4\lceil p_{2}/2\rceil+2\lfloor p_{2}/2\rfloor+1)$ (or $(4p_{2}+1)(4\lceil p_{1}/2\rceil+2\lfloor p_{1}/2\rfloor+1)$) sets $\widetilde{\omega}_{k}$. The result holds with a constant that depends on the degree ${\bf p}$.

We remark that in the figures we are mixing Bézier elements and elements in the T-mesh. In this worst case scenario, and staying away from the boundary, there is a one-to-one correspondence between both. $\square$

Remark 4.5

As already mentioned at the beginning of this section, the results can be extended to T-splines with lower continuity, that is, with repeated knots in the T-mesh. It would be necessary to adapt the definition of T-mesh by bisection in Section 4.1 to take into account the multiplicity. The definition of analysis-suitable T-splines also applies to the case of reduced regularity, while the definition of p-admissible T-meshes applies as well, taking into account that the neighborhood is defined in the T-mesh. We note that reducing the regularity reduces the support extension $\widetilde{Q}$. Even if the proofs could be extended applying the same ideas, allowing for repeated knots would require a much more intricate notation, since one could not rely on the definition of the index sets ${\cal I}_{d}^{\ell}$ of each level.

Remark 4.6

The constants appearing in the proofs of the previous results, and in particular in Propositions 4.1, 4.6 and 4.7, are independent of the mesh size but dependent on the degree. It is also important to note that the constants in (18) depend not only on the degree, but more precisely in the number of times that each single element appears in $\omega_{k}$ (and $\widetilde{\omega}_{k}$), that is, in the overlaps between the differents $\omega_{k}$ of the same generation. This can have an important impact in the performance of the preconditioner, as we will see in Section 7.

5 Space decomposition on AS T-meshes

After having introduced all the preliminary results, we are now in a position to prove the main results of the paper, that is, a space decomposition that satisfies **(A1) **and (A2).

Let ${\cal T}$ be a ${\bf p}$-admissible T-mesh obtained by successive bisections from the index T-mesh ${\cal T}_{0}$. We give a decomposition of T-splines space $\mathcal{V}:=S(\mathcal{A}_{\bf p}(\mathcal{T}))$ using these successive bisections. By (13), we have the following space decomposition:

[TABLE]

In each subspace $\mathcal{V}_{k}$, thanks to Lemma 4.4 and [1], we have the following inverse inequality

[TABLE]

where we recall that $h_{\ell_{k}}$ is the mesh size at level $\ell_{k}$, which is the level of functions in $\Phi_{k}$.

5.1 Stable decomposition

Theorem 5.1** (Space decomposition over AS T-meshes)**

For any $v\in\mathcal{V}$, there exist $v_{k}\in{\cal V}_{k},\ k=0,\ldots,N$ such that $v=\sum_{k=0}^{N}v_{k}$ and

[TABLE]

Proof. First we consider an auxiliary decomposition over uniformly refined spaces of tensor-product B-splines. Let us recall that ${\cal T}_{\ell}^{u}$ indicates, for $0\leq\ell\leq L$, the Cartesian mesh obtained after applying all possible bisections of generation $\ell$, $S_{\bf p}({\bf\Xi}^{\ell})$ is the associated space of tensor-product B-splines of degree ${\bf p}$, and ${\bf\Pi}_{\bf p,\Xi^{\ell}}:S_{\bf p}({\bf\Xi}^{L})\rightarrow S_{\bf p}({\bf\Xi}^{\ell})$ is the multivariate quasi-interpolant in (7). We state (without proof) the following well-known stable decomposition for the space $S_{\bf p}({\bf\Xi}^{L})=\sum_{\ell=0}^{L}S_{\bf p}({\bf\Xi}^{\ell})$ (see [13] for details)

Lemma 5.2** (Stable decomposition for quasi-uniform meshes)**

For any $\bar{v}\in S_{\bf p}({\bf\Xi}^{L})$, let $\bar{v}_{\ell}=({\bf\Pi}_{{\bf p},{\bf\Xi}^{\ell}}-{\bf\Pi}_{{\bf p},{\bf\Xi}^{\ell-1}})\bar{v}$ for $\ell=0,\ldots,L$, setting ${\bf\Pi}_{{\bf p},{\bf\Xi}^{-1}}:=0$. Then $\bar{v}=\sum_{\ell=0}^{L}\bar{v}_{\ell}$ is a stable decomposition in the sense that

[TABLE]

We now give a multilevel decomposition of $v$ using a sequence of projection operators ${\bf\Pi}_{\bf p}^{{\cal T}_{k}}:\mathcal{V}\rightarrow S(\mathcal{A}_{\bf p}(\mathcal{T}_{k})),\ k=0,1,\ldots,N$, which are defined in (11).

First let $v=\sum_{\ell=0}^{L}\bar{v}_{\ell}$ with $\bar{v}_{\ell}:=({\bf\Pi}_{\bf p,\Xi^{\ell}}-{\bf\Pi}_{{\bf p},{\bf\Xi}^{\ell-1}})v\in S_{\bf p}({\bf\Xi}^{\ell})$. Then the slicing operator ${\bf\Pi}_{\bf p}^{{\cal T}_{k}}-{\bf\Pi}_{\bf p}^{{\cal T}_{k-1}}$ verifies an important property:

[TABLE]

with the convention that ${\bf\Pi}_{\bf p}^{{\cal T}_{-1}}:=0$ and thus $v=\sum_{k=0}^{N}v_{k}$. If $\ell_{k}=g(\tau_{k-1})+1$ is the generation of the functions in $\Phi_{k}$, obtained after the bisection of the element $\tau_{k-1}$ in (8), then for $l\leq\ell_{k}-1=g(\tau_{k-1})$ it holds that $({\bf\Pi}_{\bf p}^{{\cal T}_{k}}-{\bf\Pi}_{\bf p}^{{\cal T}_{k-1}})\bar{v}_{l}=0$, which is shown in details in Appendix A.1. Thus for $1\leq k\leq N$, we have

[TABLE]

Also, Proposition 4.1 shows that

[TABLE]

From the definition of the regions $\omega_{k}$, applying first the inequality above and then inequality (18), we infer that for $\ell\geq 1$

[TABLE]

Also we notice that when $\ell=0$ the space is tensor product, then $v_{0}=\bar{v}_{0}$, and we trivially have

[TABLE]

Applying the discrete Hardy inequality [14, Lemma 4.3] with $s=\gamma^{2}=1/2$ to $a_{\ell}=\|\bar{v}_{\ell}\|^{2}$ and $b_{\ell}=\sum_{k:\ell_{k}=\ell}\|v_{k}\|^{2}$, we obtain

[TABLE]

Finally, Lemma 5.2 leads to

[TABLE]

Combining this inequality with (20) gives the desired estimate (21). $\square$

5.2 Strengthened Cauchy-Schwarz (SCS) inequality

The proof of the SCS inequality for T-splines relies on the same result for tensor-product B-splines. Although multilevel methods for B-splines are now rather classical, we have not explicitly found the result for tensor-product splines in the literature. For this reason, we start this section proving the SCS inequality in the tensor-product case.

Lemma 5.3** (SCS inequality for B-splines on globally quasi-uniform meshes)**

For $u_{i}\in S_{\bf p}({\bf\Xi}^{i})$ and $u_{j}\in S_{\bf p}({\bf\Xi}^{j})$ with $j\geq i$, we have

[TABLE]

where $\gamma<1$ is a constant such that $h_{i}\approx\gamma^{i}$.

Proof. Recalling that $S_{\bf p}({\bf\Xi}^{i})=\mbox{span}\{{B}_{\bf i,p}({\boldsymbol{\zeta}}),~{}{\bf i}\in{\bf I}^{i}\}$ and $S_{\bf p}({\bf\Xi}^{j})=\mbox{span}\{{B}_{\bf j,p}({\boldsymbol{\zeta}}),~{}{\bf j}\in{\bf I}^{j}\}$, any $u_{j}\in S_{\bf p}({\bf\Xi}^{j})$ can be written as

[TABLE]

which can also be classified as follows: for each Bézier element $Q^{i}\in{\cal M}_{i}^{\cal B}$ we define the sets of indices

[TABLE]

of basis functions with support completely outside of $Q^{i}$, completely contained in $Q^{i}$, or with just part of it in $Q^{i}$, where we recall that the element $Q^{i}$ is open. Then, we have

[TABLE]

Denoting $a_{Q^{i}}(u,v):=\int_{Q^{i}}\nabla u\cdot\nabla v~{}dx$, we clearly have $a_{Q^{i}}(u_{i},u_{\rm out})=0$, and thus

[TABLE]

Define $\Gamma:=\bigcup_{{\bf j}\in{\bf I}^{j}_{\rm B}}{\rm supp}(B_{\bf j,p})\cap Q^{i}$, which is the union of smaller elements $Q^{j}\in{\cal M}_{j}^{\cal B}$ contained in $Q^{i}$ and in the support of a function that is not completely contained in $Q^{i}$. It is easy to see that $|\Gamma|\approx h_{i}h_{j}$, with the implicit constant depending on the degree.

First we estimate the term corresponding to $u_{\rm B}$. Since $\mbox{supp}(u_{\rm B})\cap Q^{i}=\Gamma$, we have

[TABLE]

Using that $\Gamma\subset Q^{i}$ with $|\Gamma|\approx h_{i}h_{j}$, the inverse inequality [11, Lemma 4.5.3], and finally that $h_{j}\approx\gamma^{j}$, we obtain

[TABLE]

The inverse inequality on globally quasi-uniform meshes, followed by [7, Proposition 5.1 and Corollary 5.1] shows that

[TABLE]

where a detailed proof of the last inequality can be found in Appendix A.2. Combining the above two inequalities leads to

[TABLE]

Now, for the estimate of the term corresponding to $u_{\rm in}$ we have

[TABLE]

It follows from the inverse inequality that

[TABLE]

Noting that $\|u_{\rm in}\|_{L^{2}(Q^{i})}\lesssim\|u_{j}\|_{L^{2}(Q^{i})}$ that can be shown in a similar way to Appendix A.2, we get

[TABLE]

Combining the two estimates (23) and (24), together with $0<\gamma<1$, yields

[TABLE]

Summing over $Q^{i}$ and using the Cauchy-Schwarz inequality, we obtain

[TABLE]

which ends the proof. $\square$

Theorem 5.4** (SCS inequality for AS T-meshes)**

For any $u_{i},v_{i}\in\mathcal{V}_{i},\ 0\leq i\leq N$, we have

[TABLE]

Proof. The proof is similar to the proof of Theorem 4.6 in [14], and is divided into four steps.

1. For a fixed $i\in\{1,2,\ldots,N\}$, we denote

[TABLE]

Observe that $w_{\ell}^{i}\in S_{\bf p}({\bf\Xi}^{\ell})$ and $\ell=\ell_{j}\geq\ell_{i}$ by the assumption that if $k>k^{\prime}$, then $\ell_{k}\geq\ell_{k^{\prime}}$, while $u_{i}\in S_{\bf p}({\bf\Xi}^{\ell_{i}})$.

For any Bézier element $Q\subset{\omega}_{i}$, we apply Lemma 5.3 over $Q$ to $u_{i}$ and $w_{\ell}^{i}$ to obtain

[TABLE]

Then

[TABLE]

Moreover, from the definition of $w_{\ell}^{i}$

[TABLE]

Thus we obtain

[TABLE]

2. Fix $u_{i}$ and consider

[TABLE]

Thus we get

[TABLE]

3. We sum over $i$ with fixed generation $\ell_{i}=k$:

[TABLE]

where in the third inequality we used the usual Cauchy-Schwarz inequality and in the last inequality we used Proposition 4.7.

4. Finally we sum over all the generations $0\leq k\leq L$ to get

[TABLE]

where in the last inequality we used the inequality (see [43, Lemma 5])

[TABLE]

Since $\sum^{L}_{k=0}\sum_{i:\ell_{i}=k}=\sum^{N}_{i=0}$, using a scaled Poincaré inequality, namely,

[TABLE]

we have

[TABLE]

$\square$

6 BPX preconditioners on locally quasi-uniform AS T-meshes

After the proof of the stable decomposition **(A1) **and the SCS inequality (A2), we can apply the parallel subspace correction method. We introduce two different space decompositions: the first one, that we called micro decomposition, is based on the same subspaces in the previous two sections; the second one, that we denote macro decomposition, is based on a smaller number of subspaces, that collect all the bisections of the same generation at once.

6.1 Micro decomposition

We apply parallel subspace correction methods to the space decomposition

[TABLE]

where the subspaces ${\cal V}_{k}$ are defined in (13), and thus obtain BPX preconditioners on locally quasi-uniform AS T-meshes.

Theorem 6.1

If the preconditioner $B$ of (2) is based on the space decomposition (25) and the SPD smoothers satisfying (3), then we have

[TABLE]

Proof. It follows immediately from Theorems 5.1, 5.4 and 2.1. $\square$

Remark 6.1

We observe that we can use standard smoothers such as Jacobi and symmetric Gauss-Seidel iterations satisfying (3). Recently, Hofreither et al. in [24] proposed a new smoother based on the mass matrix and a boundary correction, that was later improved in [23]. The results in those papers show that on a uniform mesh, the multigrid with the new smoother is more robust, in the sense that convergence is independent of both the mesh size and the spline degree. The extension of this smoother to T-splines is one of our future research topics. We notice however that this would eliminate the dependency on the degree for $K_{3}$ and $K_{4}$ in (4), but not for $K_{1}$ and $K_{2}$.

6.2 Macro decomposition

As an alternative to the previous decomposition, we introduce a macro space decomposition where each level contains all the elements of the same corresponding generation. In order to define the macro space decomposition, to each generation $\ell$ we associate a subspace ${\mathcal{W}}_{\ell}$ that contains the functions added or modified after inserting all the lines of generation $\ell$. For instance, in the T-meshes of Figure 7, ${\mathcal{W}}_{0}$ would be the tensor-product space in the black mesh, and the subspaces ${\mathcal{W}}_{1}$, ${\mathcal{W}}_{2}$ and ${\mathcal{W}}_{3}$ would consist of the T-splines added or modified after the bisection sequence for adding all the red vertical edges, all the blue horizontal edges, and all the pink vertical edges, respectively.

To be more precise, starting from the Cartesian mesh ${\cal T}_{0}$, we define for $1\leq\ell\leq L$ the T-mesh

[TABLE]

which is the finest T-mesh of generation $\ell$. Similarly to the sets $\Phi_{k}$ in Section 4.3, we can now define the sets

[TABLE]

and the subspace

[TABLE]

noting that this subspace contains all the subspaces of the micro decomposition of the corresponding generation, in fact, ${\cal W}_{\ell}=\sum_{k:\ell_{k}=\ell}{\cal V}_{k}$.

Defining those ${\mathcal{W}}_{\ell}$ gives the decomposition

[TABLE]

We can use the general framework presented in Sections 4-5 to obtain the uniformly bounded condition number of the BPX preconditioner on the macro space decomposition (26).

Remark 6.2

It is worth to mention that since our macro structure is similar to the locally refined grids proposed in [9, 10, 42], BPX optimality on the space decomposition (26) can also be obtained in an analogous way to the proofs in those references.

6.3 Application to IGA

All the theoretical results presented in this paper are proved in the unit domain ${\Omega}=(0,1)^{2}$. The results can be easily extended to the IGA setting, where the physical domain $\Omega$ is defined as the image of the unit parametric domain through a parametrization ${\bf F}$, that is, $\Omega={\bf F}((0,1)^{2})$. The parametrization ${\bf F}$ can be defined as a T-spline on a coarse mesh, and both ${\bf F}$ and its inverse should be regular, see [13] for details.

7 Numerical results

We have performed some numerical tests to support our theoretical results with numerical evidence. We solve the model problem (1) with a suitable $f$ that makes it necessary to refine towards a corner, noticing however that the condition number does not depend on $f$. Since we lack a true implementation of T-splines, we have run the tests using the Octave/Matlab software GeoPDEs [17, 38] for some particular T-meshes, defining a tensor-product space for each level, and then collecting the active functions of different levels. Our implementation is clearly inefficient for T-splines, hence we do not present computational times.

7.1 Macro tests

In the first two tests we follow the macro space decomposition in Section 6.2, that is, the subspace ${\cal W}_{k}$ contains the functions that have been added or modified after inserting all the lines of generation $k$.

7.1.1 Square domain

In the first test the domain is the unit square $\Omega_{0}=(0,1)^{2}$, which is refined near the origin with the following procedure: the level zero mesh is a Cartesian grid, and then for each level we refine, in the T-mesh, the maximal square subregion such that the T-mesh remains ${\bf p}$-admissible in the sense of [33]. The number of elements of the initial Cartesian grid has been chosen to obtain a similar number of degrees of freedom for each degree. In particular, it is equal to $7^{2}$, $8^{2}$ and $10^{2}$ for biquadratic, bicubic and biquartic splines, respectively. We show in Figure 7 the refined meshes for the biquadratic and the bicubic case, after all the bisections of the fourth generation.

For the computations with the BPX preconditioner, the inclusion operator $I_{k}$ can be computed through knot insertion, while the restriction operator is simply its transpose. We have compared the results with two different smoothers $R_{k}$: Jacobi and symmetric Gauss-Seidel, with one single iteration each. The condition number of the preconditioned system is estimated computing the minimum and maximum eigenvalues with Lanczos’ method. The stopping criterion is set to an initial tolerance of $10^{-6}$ for coarse meshes. Table 1 summarizes the obtained results, where the number of degrees of freedom corresponds to the matrix size, that is, after applying Dirichlet boundary conditions.

As predicted by the theory, the numerical results show that the condition number remains bounded by a constant during $h$-refinement. As already observed in [13], the constant deteriorates with the degree $p$, and as explained in [24, 23] (see also [19]) this is due to the bad behavior of Jacobi and symmetric Gauss-Seidel smoothers with splines. The smoother recently proposed in [23] for tensor-product B-splines can be probably adapted to AS T-splines by bisection, although we expect a considerable loss of efficiency of the preconditioner, since the construction of the smoother in [23] strongly relies on the Kronecker tensor-product structure of the matrices. We also remark that, although the condition number is lower using symmetric Gauss-Seidel as a smoother, the iteration cost in our implementation is higher, thus in our tests the Jacobi smoother performs better than the Gauss-Seidel one.

7.1.2 Curved L-shaped domain

As a second test we have considered the curved L-shaped domain in Figure 8. The domain is defined with three patches, each one the image of the unit square through a different parametrization ${\bf F}$, in such a way that the reentrant corner is always the image of the origin. We have constructed, in the parametric domain, the same meshes as in the previous test, which are then mapped to the physical domain through ${\bf F}$. The relative orientation of the three patches guarantees that the meshes match on the interfaces, and basis functions on the interfaces are glued together with $C^{0}$ continuity.

We have run the same kind of tests as before, the results are summarized in Table 2. The numerical results show that the presence of the parametrization does not greatly affect the condition number, as long as the parametrization remains regular. These results are in agreement with those in [13] for tensor-product B-splines.

7.2 Comparison with other variants

7.2.1 Micro test

Our implementation is based on tensor-product spaces, therefore it is not general enough to test the preconditioner using the micro technique, that is, defining the decomposition spaces ${\cal V}_{k}$ as in (13) adding one edge at each step. However, in order to better understand how the micro decomposition could work we have tried a modified version of the method: instead of adding all the bisections of the same generation at the same time, as in the macro technique, we group aligned bisection edges of the same generation into one single level. For instance, in the example of Figure 7(c) for the bicubic case, and recalling that we start from generation zero, the first, second and third generations contain seven, six and ten new levels, respectively.

In Table 3 we show the condition numbers of the system obtained after all the lines of the same generation have been obtained, that can be compared with the ones obtained with the macro technique in Table 1. Notice that we do not show the number of degrees of freedom and the condition number of the unpreconditioned system, since they are the same as in Table 1. Although the condition number remains bounded, as predicted by the theory, the results show a loss of efficiency of the preconditioner. This seems to be a consequence of the behavior of the constants in the inequalities (18), that depend on the overlaps of the sets $\omega_{k}$ (and $\widetilde{\omega}_{k}$) of the same generation, see also Remark 4.6. Gathering all the bisections of the same generation together, as in the macro technique, reduces the number of overlaps and the value of these constants, while considering the bisections separately, as in the micro technique, increases the number of overlaps. We expect even a worse behavior for the “pure” micro decomposition of Section 6.1. We also remark that this behavior does not appear in the case of finite elements [14], because of the reduced support of $C^{0}$ finite element functions.

7.2.2 An alternative refinement

In all the previous tests the refinement is performed as in [33] to obtain ${\bf p}$-admissible meshes, that is, the elements are refined alternatively in the vertical and the horizontal direction, depending on their generation. For this new test we have tried a different kind of refinement, suitable to refine towards a corner, that preserves the multilevel structure and the analysis-suitability of the mesh. This refinement gives a higher proportion of square elements with respect to ${\bf p}$-admissible meshes, although more elongated elements can also appear. Starting from a Cartesian grid, the elements of the T-mesh in a square subregion next to the corner are split in four, bisecting them vertically and horizontally simultaneously. Then, some other elements have to be bisected, either horizontally or vertically, to maintain the analysis-suitability property. To simplify the computations, in this second step we refine together a set of aligned elements, although this is not really necessary to maintain the analysis-suitable property. All the functions added or modified during the refinement of the square region are considered to have the same generation, independently of whether the refinement is horizontal, vertical or both. A detail of the meshes obtained after three refinement steps, for the biquadratic case, is shown in Figure 9. Lines added at the same step are drawn with the same color.

For the numerical tests we use the macro technique as in Section 7.1, that is, the space ${\cal W}_{k}$ is defined after adding all the lines of the same generation. In Table 4 we show the results obtained applying the BPX preconditioner with this kind of refinement. The results show a better behavior of the BPX preconditioner compared to the meshes in Section 7.1, and in fact the numbers are very similar to the ones obtained for tensor-product B-splines in [13]. Indeed, this seems a good motivation to look for alternative refinement techniques for T-splines, that are less restrictive than the one in [33]. This may be possible after the definition of truncated T-splines [40], for which has been proved linear independence with a less restrictive condition than the one of analysis-suitable T-splines.

8 Conclusions

We have presented the analysis of an additive multilevel preconditioner for T-splines. The T-meshes are assumed to be ${\bf p}$-admissible, as in the construction by [33], which guarantees analysis-suitability of the T-splines. The optimality of the preconditioner under $h$-refinement is proved following the ideas by [14], and relies on a stable space decomposition (Theorem 5.1) and the strengthened Cauchy-Schwarz inequality (Theorem 5.4). In the method by [14] a new level is defined each time we add one edge to the mesh, and only the functions modified by this edge are added to the level. We have studied an alternative construction, where the new level is defined after adding all the edges of the same generation, and considering the newly added and modified functions as having the same level. The theoretical analysis for this alternative decomposition follows from the previous one, but the numerical results show a better performance. We also remark that, since the refinement strategy in [33] has been generalized in [32] to the three dimensional case with trivariate odd degree AS T-tsplines, our proposed preconditioners can be also generalized to the three dimensional case accordingly.

Acknowledgments

The authors would like to thank Prof. Annalisa Buffa for many fruitful discussions on the subject of the paper. The work of Durkbin Cho was supported by Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Science, ICT & Future Planning (2015R1A1A1A05001109). The work of Rafael Vázquez was partially supported by the European Research Council through the FP7 ERC Consolidator Grant no. 616563 HIGEOM (PI: Giancarlo Sangalli) and the FP7 ERC Advanced Grant no. 694515 CHANGE (PI: Annalisa Buffa), and by European Union’s Horizon 2020 research and innovation programme through the grant no. 680448 CAxMan. This support is gratefully acknowledged.

This is a pre-copyedited, author-produced PDF of an article accepted for publication in IMA Journal of Numerical Analysis following peer review. The version of record is available online at: https://doi.org/10.1093/imanum/dry032.

Appendix A Mathematical proofs of auxiliary results

A.1 Difference of T-spline projectors applied to coarse B-splines

The first auxiliary result was used in the proof of Theorem 5.1. Letting $\bar{v}_{\ell}:=({\bf\Pi}_{\bf p,\Xi^{\ell}}-{\bf\Pi}_{{\bf p},{\bf\Xi}^{\ell-1}})v\in S_{\bf p}({\bf\Xi}^{\ell}),\ 0\leq\ell\leq L$, we need to prove that for $l\leq\ell_{k}-1$

[TABLE]

which is a consequence of the following lemma.

Lemma A.1

Let ${\cal T}_{k}={\cal T}_{k-1}+b_{\tau_{k-1}}$ be a ${\bf p}$-admissible T-mesh generated by bisection of $\tau_{k-1}$, and let $\ell_{k}=g(\tau_{k-1})+1$. Then it holds that

[TABLE]

for any function $v\in S_{\bf p}({\bf\Xi}^{\ell})$, with $\ell<\ell_{k}$.

Proof. Since $S_{\bf p}({\bf\Xi}^{i})\subset S_{\bf p}({\bf\Xi}^{j})$ for $i\leq j$, we only need to prove the result for a function $v\in S_{\bf p}({\bf\Xi}^{\ell_{k}-1})$. Let us first define, similarly to $\Phi_{k}$, the set

[TABLE]

which is the collection of functions removed after the bisection $b_{\tau_{k-1}}$. Since the bisection only affects basis functions in $\Psi_{k-1}$ (removed) and $\Phi_{k}$ (added), from the definition of the dual functionals it is clear that

[TABLE]

From the nestedness of the T-spline spaces $S({\cal A}_{\bf p}({\cal T}_{k-1}))\subset S({\cal A}_{\bf p}({\cal T}_{k}))$ (see [33, Corollary 5.8]), it holds that $\mbox{span}(\Psi_{k-1})\subset\mbox{span}(\Phi_{k})$, and by the local linear independence of analysis-suitable T-splines, we have

[TABLE]

Moreover, from Lemma 4.4, and in particular from (16)-(17), we know that functions in $\Phi_{k}$ see all the lines of the uniform Cartesian mesh ${\cal T}^{u}_{\ell_{k}-1}$ in the vicinity of $\tau_{k-1}$, and using the previous property for the supports, we know that the same is true for functions in $\Psi_{k-1}$.

Let us define, starting from the uniform Cartesian mesh ${\cal T}^{\prime}_{0}={\cal T}^{u}_{\ell_{k}-1}$ the auxiliary family of ${\bf p}$-admissible T-meshes

[TABLE]

where $k_{0}$ is the minimum integer such that $g(\tau_{k_{0}})=g(\tau_{k-1})$ in the construction (8) of the T-mesh ${\cal T}_{k+1}$. In other words, starting from the Cartesian grid of level $\ell_{k}-1$, we perform the bisections of the same generation until we reach $b_{\tau_{k-1}}$. From what we have just seen, it is clear that $\Phi^{\prime}_{k-k_{0}}=\Phi_{k}$ and $\Psi^{\prime}_{k-k_{0}-1}=\Psi_{k-1}$, and by the nestedness of ${\bf p}$-admissible T-splines it holds $S_{\bf p}({\bf\Xi}^{\ell_{k}-1})=S({\cal A}_{\bf p}({\cal T}^{\prime}_{0}))\subset S({\cal A}_{\bf p}({\cal T}^{\prime}_{j}))$ for any $j=0,\ldots,k-k_{0}$. Thus we have, for $v\in S_{\bf p}({\bf\Xi}^{\ell_{k}-1})$

[TABLE]

which ends the proof. $\square$

A.2 Appropriate discrete norms

The second auxiliary result, that we present in this subsection, is $\|u_{\rm B}\|_{L^{2}(\Gamma)}\lesssim\|u_{j}\|_{L^{2}(Q^{i})}$ in the proof of Lemma 5.3.

Let $Q^{i}\in{\cal M}^{\cal B}_{i}$ be a Bézier element of the Cartesian grid of generation $i$. For $j\geq i$ we define the set of indices

[TABLE]

which is also the union of indices ${\bf I}^{j}_{\rm in}$ and ${\bf I}^{j}_{\rm B}$ defined as in the proof of Lemma 5.3, namely, ${\bf I}^{j}(Q^{i})={\bf I}^{j}_{\rm in}\cup{\bf I}^{j}_{\rm B}$. By local linear independence of B-splines, a function $z\in S_{\bf p}({\bf\Xi}^{j})$ restricted to $Q^{i}$ can be written as

[TABLE]

Let also $Q^{j}\in{\cal M}^{\cal B}_{j}$ such that $Q^{j}\subset Q^{i}$, and define ${\bf I}^{j}(Q^{j}):=\{{\bf j}\in{\bf I}^{j}:Q^{j}\subset{\rm supp}(B_{\bf j,p})\}$. Obviously, ${\bf I}^{j}(Q^{j})\subset{\bf I}^{j}(Q^{i})$. We define a local discrete norm $|\cdot|_{Q^{j}}$ as follows: for any $z\in S_{\bf p}({\bf\Xi}^{j})|_{Q^{i}}$,

[TABLE]

and also define the global discrete norm $|\cdot|_{Q^{i}}$ on $Q^{i}$ as follows: for any $z\in S_{\bf p}({\bf\Xi}^{j})|_{Q^{i}}$,

[TABLE]

Proposition A.1

Let $Q^{i}\in{\cal M}^{\cal B}_{i}$ and let, for $j\geq i$, $Q^{j}\in{\cal M}^{\cal B}_{j}$ such that $Q^{j}\subset Q^{i}$. Then for any $z\in S_{\bf p}({\bf\Xi}^{j})|_{Q^{i}}$ it holds

[TABLE]

and

[TABLE]

Proof. The proof can be found in [7, Proposition 5.1 and Corollary 5.1]. $\square$

By Proposition A.1, it is clear that

[TABLE]