# Isogeometric Tearing and Interconnecting Solvers for Linearized Elasticity
in multi-patch Isogeometric Analysis

Jarle Sogn  Department of Mathematics, University of Oslo, Postboks 1053,
Blindern, 0316 Oslo, Norway  Stefan Takacs  Institute of Numerical
Mathematics, Johannes Kepler University Linz, Altenberger Str. 69, 4040 Linz,
Austria

###### Abstract

We consider the linearized elasticity equation, discretized with multi-patch
Isogeometric Analysis. A standard discretization error analysis is based on
Korn’s inequality, which degrades for certain geometries, such as long and
thin cantilevers. This phenomenon is known as geometry locking. We observe
that high-order methods, like Isogeometric Analysis is beneficial in such a
setting. The main focus of this paper is the construction and analysis of a
domain decomposition solver, namely an Isogeometric Tearing and
Interconnecting (IETI) solver, where we prove that the convergence behavior
does not depend on the constant of Korn’s inequality for the overall domain,
but only on the corresponding constants for the individual patches. Moreover,
our analysis is explicit in the choice of the spline degree. Numerical
experiments are provided which demonstrates the efficiency of the proposed
solver.

###### keywords:

Linear elasticity , Isogeometric Analysis , Convergence analysis , Domain
decomposition solver

††journal: arXiv

##  1 Introduction

In this paper, we consider a discretization and a solution strategy for the
linearized elasticity problem. For small deformations $\mathbf{u}$ of a body
in static equilibrium, the equation

$-2\mu\;\text{div}\,\epsilon(\mathbf{u})-\lambda\;\nabla\text{div}\,\mathbf{u}=\mathbf{f}$
(1)

holds in a variational sense, where
$\epsilon(\mathbf{u})=\frac{1}{2}(\nabla\mathbf{u}+\nabla^{\top}\mathbf{u})$
is the symmetric gradient, $\lambda,\mu>0$ are Lamé coefficients modelling the
material properties and $\mathbf{f}$ is the external volumetric force acting
on the body, cf. [10]. In this model, the computational domain
$\Omega\subset\mathbb{R}^{d}$ models the shape of the undeformed object with
boundary $\partial\Omega$. As boundary conditions, we use Dirichlet boundary
conditions to prescribe the displacement on some part $\Gamma_{D}$ of the
boundary with positive measure and Neumann boundary conditions on the
remainder $\Gamma_{N}:=\partial\Omega\setminus\Gamma_{D}$ of the boundary to
prescribe acting forces:

$\displaystyle\mathbf{u}$ $\displaystyle=\mathbf{g}_{D}$
$\displaystyle\mbox{on}\quad\Gamma_{D},$ (2) $\displaystyle
2\mu\;\epsilon(\mathbf{u})\;\mathbf{n}+\lambda\;\text{div}\,\mathbf{u}\;\mathbf{n}$
$\displaystyle=\mathbf{g}_{N}$ $\displaystyle\mbox{on}\quad\Gamma_{N},$ (3)

where $\mathbf{g}_{D}$ denotes the prescribed displacement and
$\mathbf{g}_{N}$ the prescribed boundary force density.

Several challenges may arise when discretizing and solving Equation (1). The
ratio $\mu/\lambda$ takes small values for nearly incomressible materials. The
limit case $\mu/\lambda\to 0$ can be understood as the case of completely
incompressible materials. In these cases, the well-posedness (in the sense of
the Lax-Milgram lemma) deteriorates, which might lead to a poor approximation
of the solution with standard Galerkin discretizations. This phenomenon is
known as material locking and geometry locking. A common remedy is to
introduce an auxiliary variable and derive a mixed variational formulation,
cf. [5, 24, 2, 3]. The focus of this paper is not on material locking; thus,
we assume $\mu/\lambda$ is well bounded away from zero.

Another challenge is related to the symmetric gradient. For showing well-
posedness in the standard Sobolev space $H^{1}$, one typically uses Korn’s
inequality, which introduces a constant that usually depends on the geometry
of the object in question. This happens for example for long and thin
cantilevers. Also such domains might lead to a poor approximation of the
solution with standard Galerkin discretizations, which is known as geometry
locking. We will discuss this issue in detail and observe that high-order
methods help to mitigate related issues, compared to standard discretizations
with piecewise linear functions.

We use Isogeometric Analysis (IgA) as discretization method. IgA was
introduced in the seminal paper [12] as a new technique for discretizing
partial differential equations (PDEs) that improves the integration of
simulation and computer aided design (CAD). In IgA, both the computational
domain and the solution of the PDE are represented as linear combination of
tensor-product B-splines or non-uniform rational B-splines (NURBS). Only
simple domains can be parameterized using a single geometry mapping. More
complicated domains are usually composed of several patches, where each patch
is parameterized independently (multi-patch domains). For more information on
IgA, see the review paper [4] and the references therein.

The main focus of this paper is the intoduction and the analysis of a fast
domain decomposition solver for the linearized elasticity problem, discretized
using IgA on multi-patch domains. We consider FETI-DP methods, which were
originally introduced in [8] and has been used to solve linear elasticity,
see, e.g., [13, 9, 14]. For multi-patch IgA domains, the patches themselves
serve as a natural choice for the substructures. FETI-DP was first adapted to
IgA in [15] and named the Dual-Primal Isogeometric Tearing and Interconnecting
(IETI-DP) method. In recent years, there has been a lot of research in IETI-DP
methods, see, e.g., [11, 20, 23, 27, 21, 17].

Isogeometric domain decomposition methods for mixed linear elasticity systems
have been studied in [27, 18], where convergence estimates are provided. We
study the primal formulation and use the results from [16, 20] to prove a
convergence estimate which is explicit with respect to both the spline degree
and grid-size (Theorem 7). The analysis shows that the performance of the
solver does not depend on the shape of the overall computational domain, but
only on the shapes of the individual patches. This means that one can obtain a
robust IETI solver by setting up an adequate subdivision of the computational
domain into patches.

This paper is organized as follows. In Section 2, we introduce the variational
formulation of the model problem and its discretization. The well-posedness
and the discretization errors are discussed in the subsequent Section 3. The
proposed solver is introduced in Section 4. The following Section 5 is
dedicated to its analysis. We close with numerical results and concluding
remarks in Section 6.

##  2 Problem formulation and discretization

In this section we formulate the variational formulation of the elasticity
problem (1) – (3) and its discretization. Let $\Omega\subset\mathbb{R}^{2}$ be
an open, connected and bounded domain with Lipschitz boundary
$\partial\Omega$. $L^{2}(\Omega)$ and $H^{s}(\Omega)$ denote the standard
Lebesgue and Sobolev spaces with standard norms $\|\cdot\|_{L^{2}(\Omega)}$,
$\|\cdot\|_{H^{s}(\Omega)}$, seminorms $|\cdot|_{H^{s}(\Omega)}$, and scalar
products $(\cdot,\cdot)_{L^{2}(\Omega)}$. Analogously, we define these spaces,
norms, seminorms and scalar products on other domains and on manifolds in
$\mathbb{R}^{2}$, like boundaries. For ease of notation, we assume that
Dirichlet boundary conditions are homogeneous. So, the displacement has to be
in $H^{1}_{0,D}(\Omega):=\left\\{v\in
H^{1}(\Omega):v|_{\Gamma_{D}}=0\right\\}$. Here and in what follows a
restriction to (a part of) the boundary is to be understood in the sense of a
trace operator. By standard arguments, the following problem is the weak form
of the elasticity equations.

######  Problem 1.

Given $\mathbf{f}\in[L^{2}(\Omega)]^{2}$ and
$\mathbf{g}_{N}\in[L^{2}(\Gamma_{N})]^{2}$, find
$\mathbf{u}\in[H^{1}_{0,D}(\Omega)]^{2}$ such that

$a(\mathbf{u},\mathbf{v})=\langle F,\mathbf{v}\rangle\quad\mbox{for
all}\quad\mathbf{v}\in[H^{1}_{0,D}(\Omega)]^{2},$ (4)

where

$\displaystyle a(\mathbf{u},\mathbf{v})$
$\displaystyle:=2\mu(\epsilon(\mathbf{u}),\epsilon(\mathbf{v}))_{L^{2}(\Omega)}+\lambda\;(\text{div}\,\mathbf{u},\text{div}\,\mathbf{v})_{L^{2}(\Omega)},\mbox{
and }$ $\displaystyle\langle F,\mathbf{v}\rangle$
$\displaystyle:=(\mathbf{f},\mathbf{v})_{L^{2}(\Omega)}+(\mathbf{g}_{N},\mathbf{v})_{L^{2}(\Gamma_{N})}.$

We discretize the elasticity equations by means of multi-patch Isogeometric
Analysis. We assume that the computational domain $\Omega$ is composed of $K$
non-overlapping patches $\Omega^{(k)}$, i.e., the domains $\Omega^{(k)}$ are
open, connected and bounded with Lipschitz boundary such that

$\displaystyle\overline{\Omega}=\bigcup_{k=1}^{K}\overline{\Omega^{(k)}}\quad\text{and}\quad\Omega^{(k)}\cap\Omega^{(\ell)}=\emptyset\quad\text{for
all}\quad k\neq\ell,$

where $\overline{T}$ denotes the closure of the set $T$. We require that the
patches form an admissible decomposition, i.e., that there are no T-junctions.
This means that for any two patch indices $k\not=\ell$, the set
$\partial{\Omega^{(k)}}\cap\partial{\Omega^{(\ell)}}$ is either a common edge
$\Gamma^{(k,\ell)}:=\partial{\Omega^{(k)}}\cap\partial{\Omega^{(\ell)}}$, a
common vertex, or empty (cf. [20, Ass. 2]). Moreover, we assume that both the
Neumann boundary and the Dirichlet boundary consist of whole edges.

For each patch index $k$, the set $\mathcal{N}_{\Gamma}(k)$ contains the
indices $\ell$ of patches $\Omega^{(\ell)}$ that share an edge with
$\Omega^{(k)}$. The common vertices of two or more patches – that are not
located on the Dirichlet boundary – are denoted by
$\mathbf{x}_{1},\ldots,\mathbf{x}_{J}$. For each $j=1,\ldots,J$, the set
$\mathcal{N}_{\mathbf{x}}(j)$ contains the indices of all patches
$\Omega^{(k)}$ such that $x_{j}\in\partial\Omega^{(k)}$. We assume that there
is a constant $C_{1}$ such that

$|\mathcal{N}_{\mathbf{x}}(j)|\leq C_{1}\quad\mbox{for all}\quad j=1,\ldots,J$
(5)

(cf. [20, Ass. 3]). The set of all vertices associated to one patch is denoted
by

$\mathcal{N}_{\mathbf{x}}^{-1}(k):=\\{j:k\in\mathcal{N}_{\mathbf{x}}(j)\\}.$

Each patch $\Omega^{(k)}$ is parameterized by a geometry mapping

$\displaystyle\mathbf{G}_{k}:\widehat{\Omega}:=(0,1)^{2}\rightarrow\Omega^{(k)}:=\mathbf{G}_{k}(\widehat{\Omega})\subset\mathbb{R}^{2},$
(6)

which can be continuously extended to the closure of the parameter domain
$\widehat{\Omega}$. In IgA, the geometry mapping is typically represented
using B-splines or NURBS. As usual, the computational methods do not depend on
such a representation. We only assume that the geometry mappings are not too
much distorted, i.e., there is a constant $C_{2}>0$ such that

$\|\nabla\mathbf{G}_{k}\|_{L^{\infty}(\widehat{\Omega})}\leq
C_{2}\,H_{k}\quad\text{and}\quad\|(\nabla\mathbf{G}_{k})^{-1}\|_{L^{\infty}(\widehat{\Omega})}\leq
C_{2}\,H_{k}^{-1}$ (7)

holds for all $k=1,\ldots,K$, where $H_{k}$ is the diameter of the patch
$\Omega^{(k)}$ (cf. [20, Ass. 1]).

On the parameter domain $\widehat{\Omega}=(0,1)^{2}$, we choose a B-spline
space of some freely chosen spline degree
$\mathpzc{p}\in\mathbb{N}:=\\{1,2,3,\ldots\\}$. The spline space depends on a
freely chosen $\mathpzc{p}$-open knot vector
$\Xi^{(k,\delta)}:=(\xi^{(k,\delta)}_{0},\cdots,\xi^{(k,\delta)}_{N^{(k,\delta)}})$
for each patch $k$ and each spacial direction $\delta\in\\{1,2\\}$. We allow
repeated (inner) knots in order to reduce smoothness; since we require a
$H^{1}$-conforming discretization, no more than $\mathpzc{p}$ interior knots
may be repeated. Based on these knot vectors, the corresponding B-splines are
obtained by the Cox-de Boor formula (cf. [7, (2.1) and (2.2)]); these
B-splines span the spline space $S^{(k,\delta,\mathpzc{p})}\subset
H^{1}(0,1)$. Based on these univariate splines, we introduce the corresponding
tensor product spline space as
$S^{(k,\mathpzc{p})}:=S^{(k,1,\mathpzc{p})}\otimes S^{(k,2,\mathpzc{p})}$ as
discretization space on the parameter domain $\widehat{\Omega}$, equipped with
the standard tensor product basis.

The function spaces on the physical patches $\Omega^{(k)}$ are defined via the
_pull-back principle_ :

$\mathbf{V}^{(k)}:=\left\\{v:v\circ\mathbf{G}_{k}\in
S^{(k,\mathpzc{p})},\;v\big{|}_{\Gamma_{D}\cap\partial\Omega^{(k)}}=0\right\\}^{2}.$

The grid size $\widehat{h}_{k}$ on the parameter domain and the grid size
$h_{k}$ on the physical patch are defined by

$\widehat{h}_{k}:=\max\left\\{\xi_{i}^{(k,\delta)}-\xi_{i-1}^{(k,\delta)}\,:\,i=1,\ldots,N^{(k,\delta)},\,\delta=1,2\right\\}\quad\mbox{and}\quad
h_{k}:=H_{k}\widehat{h}_{k},$ (8)

where the definition of the latter is motivated by (7). Analogously, we define
the minimum grid size

$\widehat{h}_{\min,k}:=\min\left\\{\xi_{i}^{(k,\delta)}-\xi_{i-1}^{(k,\delta)}\neq
0\,:\,i=1,\ldots,N^{(k,\delta)},\,\delta=1,2\right\\}$

and $h_{\min,k}:=H_{k}\widehat{h}_{\min,k}$. We assume the grid to be quasi-
uniformity, i.e., there is a constant $C_{3}$ such that

$\frac{\widehat{h}_{k}}{\widehat{h}_{\min,k}}=\frac{h_{k}}{h_{\min,k}}\leq
C_{3}\text{ for }k=1,\ldots,K,$ (9)

cf. [20, Ass. 4].

Finally, we introduce a global discretization space. We assume that the
discretization is fully matching (see [20, Ass. 5]), this means that for any
two patches $\Omega^{(k)}$ and $\Omega^{(\ell)}$, sharing an edge
$\Gamma^{(k,\ell)}$, the traces agree,
$\mathbf{V}^{(k)}\big{|}_{\Gamma^{(k,\ell)}}=\mathbf{V}^{(\ell)}\big{|}_{\Gamma^{(k,\ell)}}$
and are represented using the same basis. So, as global discretization space,
we use

$\mathbf{V}=\left\\{\mathbf{v}\in[H^{1}_{0,D}(\Omega)]^{2}\,:\,\mathbf{v}|_{\Omega^{(k)}}\in\mathbf{V}^{(k)}\text{
for }k=1,\ldots,K\right\\}.$

For non-matching interfaces or for patches with T-junctions, discontinuous
Galerkin methods could be considered, see, e.g., [21]. This, however goes
beyond this paper.

The discretized problem reads as follows.

######  Problem 2.

Given $\mathbf{f}\in[L^{2}(\Omega)]^{2}$ and
$\mathbf{g}_{N}\in[L^{2}(\Gamma_{N})]^{2}$, find $\mathbf{u}_{h}\in\mathbf{V}$
such that

$a(\mathbf{u}_{h},\mathbf{v}_{h})=\langle
F,\mathbf{v}_{h}\rangle\quad\mbox{for all}\quad\mathbf{v}_{h}\in\mathbf{V},$
(10)

where $a(\cdot,\cdot)$ and $\langle F,\cdot\rangle$ are as in Problem 1.

######  Notation 1.

Here and in what follows, $c$ denotes a generic positive constant that only
depends on the constants $C_{1}$, $C_{2}$ and $C_{3}$ from (5), (7) and (9).
We write $a\lesssim b$ if there is such a generic constant such that $a\leq
c\,b$. We write $a\eqsim b$ if $a\lesssim b\lesssim a$.

##  3 Well posedness and discretization error

To show well-posedness both of Problem 1 and the discretized Problem 2, we use
the Lax-Milgram lemma [6, Theorem 2.7.7], which is applicable since both
$[H^{1}_{0,D}(\Omega)]^{2}$ and $\mathbf{V}$ are Hilbert spaces. Thus,
existence and uniqueness of solutions to the Problems 1 and 2 and their
continuous dependence on the data is guaranteed if $a(\cdot,\cdot)$ is
coercive and bounded and if $\langle F,\cdot\rangle$ is bounded; the latter
follows directly from the Cauchy-Schwarz inequality and the trace theorem. For
the further analysis, the estimates on the bilinear form are more relevant;
also here, the boundedness is straight forward:

$a(\textbf{u},\textbf{v})\leq(2\mu+2\lambda)\|\epsilon(\mathbf{u})\|_{L^{2}(\Omega)}\|\epsilon(\mathbf{v})\|_{L^{2}(\Omega)}\leq(2\mu+2\lambda)|\mathbf{u}|_{H^{1}(\Omega)}|\mathbf{v}|_{H^{1}(\Omega)}.$
(11)

The crucial point is the coercivity. Assuming

$\alpha_{\Omega,\Gamma_{D}}:=\inf_{\mathbf{u}\in
H^{1}_{0,D}(\Omega)}\frac{\|\epsilon(\mathbf{u})\|_{L^{2}(\Omega)}}{|\mathbf{u}|_{H^{1}(\Omega)}}>0,$
(12)

we obtain

$a(\textbf{u},\textbf{u})\geq
2\mu\|\epsilon(\mathbf{u})\|_{L^{2}(\Omega)}^{2}\geq
2\mu\alpha_{\Omega,\Gamma_{D}}^{2}|\mathbf{u}|_{H^{1}(\Omega)}^{2}.$

Using a standard Poincaré–Friedrichs inequality, cf. [6, Proposition 5.3.4],
we know that

$C_{F}\|\mathbf{u}\|_{H^{1}(\Omega)}\leq|\mathbf{u}|_{H^{1}(\Omega)}\leq\|\mathbf{u}\|_{H^{1}(\Omega)}\quad\mbox{for
all}\quad\mathbf{u}\in\mathbf{V},$

where the constant $C_{F}>0$ depends on $\Omega$ and $\Gamma_{D}$. Since the
Poincaré–Friedrichs inequality guarantees that $|\cdot|_{H^{1}(\Omega)}$ is
actually a norm on $H^{1}_{0,D}(\Omega)$, we state the following estimates
using the _norm_ $|\cdot|_{H^{1}(\Omega)}$.

An estimate of the kind of (12) is obtained using Korn’s second inequality.
The following Lemma 1 is a handy statement that allows to the construction of
various kinds of Korn’s inequalities, which we will use throughout this paper.

######  Remark 1.

Although the constants of both the Korn and the Poincaré–Friedrichs inequality
depend on the size of $\Omega$, the dependence of the Poincaré–Friedrichs
inequality is less severe since the $L^{2}$-norm and the $H^{1}$-seminorm
scale differently with the diameter of $\Omega$. This is different for Korn’s
inequality: $\|\epsilon(\cdot)\|_{L^{2}(\Omega)}$ and
$|\cdot|_{H^{1}(\Omega)}$ both are invariant under scaling (or scale like
$(\mathrm{diam}\,\Omega)^{d/2-1}$ for $d\neq 2$).

The following Lemma is used to show (12).

######  Lemma 1.

Let $\ell(\mathbf{u}):[H^{1}(\Omega)]^{2}\to\mathbb{R}$ be a bounded linear
operator such that $\ell(\mathbf{r})\neq 0$ for

$\mathbf{r}(x,y):=\begin{pmatrix}-y\\\ x\end{pmatrix}.$ (13)

Then, there is a constant $\alpha_{\Omega,\ell}>0$ that only depends on
$\Omega$ and $\ell$ such that

$\|\epsilon(\mathbf{v})\|_{L^{2}(\Omega)}+|\ell(\mathbf{v})|\geq\alpha_{\Omega,\ell}|\mathbf{v}|_{H^{1}(\Omega)}\quad\mbox{for
all}\quad\mathbf{v}\in[H^{1}(\Omega)]^{2}.$ (14)

A poof of this Lemma is given at the end of this section. Since the proof is
not constructive, we cannot give a lower bound for $\alpha_{\Omega,\ell}$ in
general. Note that the conditions of the Lemma are satisfied for
$\ell_{\Gamma_{D}}(\mathbf{v}):=(\mathbf{v},\mathbf{r})_{L^{2}(\Gamma_{D})}$
with $\mathbf{r}$ as in (13); thus we obtain (12) with a constant
$\alpha_{\Omega,\Gamma_{D}}>0$ that depends on $\Omega$ and $\Gamma_{D}$.

######  Remark 2.

A simple scaling argument shows that the constant $\alpha_{\Omega,\Gamma_{D}}$
does not depend on the diameter of $\Omega$, but only on its shape.

Since the assumptions of the Lax-Milgram lamma hold, we know existence and
uniqueness of a solution $\mathbf{u}$ to Problem 1 and $\mathbf{u}_{h}$ to
Problem 2. Using Ceá’s lemma (in its version [6, Remark 2.8.5] for symmetric
bilinear forms), we can bound the discretization error
$|\mathbf{u}-\mathbf{u}_{h}|_{H^{1}(\Omega)}$ by the approximation error as
follows:

$|\mathbf{u}-\mathbf{u}_{h}|_{H^{1}(\Omega)}\leq\alpha_{\Omega,\Gamma_{D}}^{-1}\,\sqrt{\frac{\mu+\lambda}{\mu}}\inf_{\mathbf{v}_{h}\in\mathbf{V}}|\mathbf{u}-\mathbf{v}_{h}|_{H^{1}(\Omega)}.$

This can be further estimated using an approximation error estimate.

######  Lemma 2.

Let $q\in\\{1,\cdots,\mathpzc{p}\\}$ and assume that the geometry functions
$\mathbf{G}_{k}$ are smooth enough such that there is a constant
$\widetilde{C}_{2}$ such that

$|\mathbf{v}\circ\mathbf{G}_{k}|_{H^{q+1}(\widehat{\Omega})}\leq\widetilde{C}_{2}H_{k}^{q}\|\mathbf{v}\|_{H^{q+1}(\Omega^{(k)})}\quad\mbox{
for all }\quad\mathbf{v}\in[H^{q+1}(\Omega)]^{2}.$ (15)

Then, the estimate

$\inf_{\mathbf{v}_{h}\in\mathbf{V}}|\mathbf{v}-\mathbf{v}_{h}|_{H^{1}(\Omega)}^{2}\lesssim\widetilde{C}_{2}^{2}\pi^{-2q}\sum_{k=1}^{K}h_{k}^{2q}\|\mathbf{v}\|_{H^{q+1}(\Omega^{(k)})}^{2}\;\mbox{holds
for all}\;\mathbf{v}\in[H^{q+1}(\Omega)]^{2}.$ (16)

This lemma is a slight extension of [25, Theorem 3.6]; for completeness, a
proof is given at the end of this section. Using this lemma, we immediately
obtain the discretization error estimate

$|\mathbf{u}-\mathbf{u}_{h}|_{H^{1}(\Omega)}^{2}\lesssim\widetilde{C}_{2}^{2}\alpha_{\Omega,\Gamma_{D}}^{-2}\,\frac{\mu+\lambda}{\mu}\pi^{-2q}\sum_{k=1}^{K}h_{k}^{2q}\|\mathbf{u}\|_{H^{q+1}(\Omega^{(k)})}^{2}.$
(17)

We observe that the convergece deteriorates if $\alpha_{\Omega,\Gamma_{D}}\to
0$ or $\mu/\lambda\to 0$. The latter happens in case of almost incompressible
materials; incompressible materials correspond to the limit case. This effect
is known as material locking. To avoid material locking, one usually considers
a saddle point formulation of the elasticity equations, however, this is not
in the scope of this paper.

The deterioration of the convergence in case $\alpha_{\Omega,\Gamma_{D}}$ is
too small, is usually known as geometry locking. In the following, we show
cases where this can happen and show that the actual convergence might still
be good, even if $\alpha_{\Omega,\Gamma_{D}}\to 0$. These discussions are done
for a spacial case. Consider a rectangular domain
$\Omega_{K}:=(0,K)\times(0,1)$ with $K\in\\{1,2,3,\cdots\\}$, consisting of
square patches $\Omega^{(k)}:=(k-1,k)\times(0,1)$ for $k=1,\ldots,K$, i.e.,

$\overline{\Omega_{K}}=\overline{(0,K)\times(0,1)}=\bigcup_{k=1}^{K}\overline{\Omega^{(k)}}=\bigcup_{k=1}^{K}\overline{(k-1,k)\times(0,1)},$

see Figure 1. The patches are parameterized in the canonical way by
$\mathbf{G}_{k}(\widehat{x},\widehat{y}):=(\widehat{x}+k-1,\widehat{y})^{\top}$.
We consider the linearized elasticity equation with homogeneous Dirichlet
boundary conditions on $\Gamma_{D}:=\Gamma^{(0)}$ and Neumann boundary
conditions on $\Gamma_{N}:=\partial\Omega\setminus\Gamma_{D}$, the remainder
of the boundary.

$\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}x$$\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}y$$\Omega^{(1)}$$\Omega^{(2)}$$\Omega^{(3)}$$\Omega^{(4)}$$\Omega^{(5)}$$\Gamma_{D}=\Gamma^{(0)}$$\Gamma^{(1)}$$\Gamma^{(2)}$$\Gamma^{(3)}$$\Gamma^{(4)}$$\Gamma^{(5)}$

Figure 1: Rectangular domain, consisting of $K=5$ square patches

The following Lemma states that the constant

$\alpha_{\Omega_{K},\Gamma_{D}}:=\inf_{\mathbf{u}\in
H^{1}_{0,D}(\Omega_{K})}\frac{\|\epsilon(\mathbf{u})\|_{L^{2}(\Omega_{K})}}{|\mathbf{u}|_{H^{1}(\Omega_{K})}}$

for the corresponding Korn inequality grows linearly with the number of
patches.

######  Lemma 3.

For all $K\in\\{1,2,3,\ldots\\}$, we have
$\alpha_{\Omega_{K},\Gamma_{D}}\eqsim K^{-1}$.

###### Proof.

An upper bound for $\alpha_{\Omega_{K},\Gamma_{D}}$ is obtained for the choice
$\mathbf{v}(x,y):=((2y-1)x,-x^{2})^{\top}$, which obviously satisfies the
Dirichlet conditions. For this choice, we can explicitly compute
$\alpha_{\Omega_{K},\Gamma_{D}}^{-1}=\sqrt{1+2K^{2}}\geq K$. (Note that for
$\mathpzc{p}\geq 2$, the function $\mathbf{v}(x,y)$ is in the spline space
$\mathbf{V}$, which means that the constant $\alpha_{\Omega_{K},\Gamma_{D}}$
cannot be improved by restricting to the discrete space $\mathbf{V}$.)

For a lower bound, we use Korn’s second inequality for the parameter domain
$\widehat{\Omega}:=(0,1)^{2}$. Let
$\mathbf{v}(x,y)=(v_{1}(x,y),v_{2}(x,y))^{\top}$ be arbitrary and define
$\mathbf{w}(x,y)=(w_{1}(x,y),w_{2}(x,y))^{\top}:=(v_{1}(Kx,y),\allowbreak
K^{-1}v_{2}(Kx,y))^{\top}$. Then, we have using simple scaling arguments

$\displaystyle\frac{|\mathbf{v}|_{H^{1}(\Omega)}^{2}}{\|\epsilon(\mathbf{v})\|_{L^{2}(\Omega)}^{2}}\leq
1+\frac{\|\partial_{y}v_{1}\|_{L^{2}(\Omega)}^{2}+\|\partial_{x}v_{2}\|_{L^{2}(\Omega)}^{2}}{\|\partial_{x}v_{1}\|_{L^{2}(\Omega)}^{2}+\tfrac{1}{2}\|\partial_{y}v_{1}+\partial_{x}v_{2}\|_{L^{2}(\Omega)}^{2}+\|\partial_{y}v_{2}\|_{L^{2}(\Omega)}^{2}}$
$\displaystyle\quad\leq
1+\frac{\|\partial_{y}w_{1}\|_{L^{2}(\widehat{\Omega})}^{2}+\|\partial_{x}w_{2}\|_{L^{2}(\widehat{\Omega})}^{2}}{K^{-2}\|\partial_{x}w_{1}\|_{L^{2}(\widehat{\Omega})}^{2}+\tfrac{1}{2}\|\partial_{y}w_{1}+\partial_{x}w_{2}\|_{L^{2}(\widehat{\Omega})}^{2}+K^{2}\|\partial_{y}w_{2}\|_{L^{2}(\widehat{\Omega})}^{2}}$
$\displaystyle\quad\lesssim
K^{2}\frac{|\mathbf{w}|_{H^{1}(\widehat{\Omega})}^{2}}{\|\epsilon(\mathbf{w})\|_{L^{2}(\widehat{\Omega})}^{2}}\leq
K^{2}\alpha_{\widehat{\Omega}}^{2},$

where $\alpha_{\widehat{\Omega}}>0$ is the constant resulting from the use of
Korn’s inequality on the parameter domain $\widehat{\Omega}$. Since this
constant is independent of the actual geometry (particularly the number of
patches $K$), we have $\alpha_{\widehat{\Omega}}\eqsim 1$. This finishes the
proof. ∎

Based on these statements, the discretization error estimate (17) degrades for
$K\to\infty$. For $\lambda=0$, the following lemma shows that the proposed
isogeometric discretization does not suffer from geometry locking. This
estimate is restricted to $\mathpzc{p}\geq 2$, which hints that it make sense
to consider discretizations other than the standard linear finite elements.

######  Lemma 4.

Let $\lambda=0$ and $\mathpzc{p}\geq 2$. For all $K\in\\{1,2,3,\ldots\\}$, we
have

$|\mathbf{u}-\mathbf{u}_{h}|_{H^{1}(\Omega)}^{2}\lesssim\inf_{\mathbf{v}_{h}\in\mathbf{V}}|\mathbf{u}-\mathbf{v}_{h}|_{H^{1}(\Omega)}^{2},$
(18)

where $\mathbf{u}$ and $\mathbf{u}_{h}$ are the solutions to the Problem 1 and
the discretized Problem 2, respectively. Under the assumptions of Lemma 2, we
have further

$|\mathbf{u}-\mathbf{u}_{h}|_{H^{1}(\Omega)}^{2}\lesssim\widetilde{C}_{2}^{2}\pi^{-2q}\sum_{k=1}^{K}h_{k}^{2q}\|\mathbf{u}\|_{H^{q+1}(\Omega^{(k)})}^{2}.$

###### Proof.

Let $\ell(w_{1},w_{2}):=\int_{0}^{1}(2y-1)\,w_{1}(1,y)\,\mathrm{d}y$. Due to
the trace inequality, it is well-defined and bounded in
$[H^{1}(\widehat{\Omega})]^{2}$. Moreover, we immediately observe that
$\ell(\mathbf{r})\neq 0$ for $\mathbf{r}$ as in (13). Thus, Lemma 1 yields

$|\mathbf{w}|_{H^{1}(\widehat{\Omega})}\leq\alpha_{\widehat{\Omega},\ell}^{-1}\left(\|\epsilon(\mathbf{w})\|_{L^{2}(\widehat{\Omega})}+|\ell(\mathbf{w})|\right),$

where $\alpha_{\widehat{\Omega},\ell}>0$ is the constant from Lemma 1. Since
this constant is independent of the actual geometry (particularly the number
of patches $K$), we have $\alpha_{\widehat{\Omega},\ell}\eqsim 1$. Since this
statement is invariant under translations, we obtain also

$|\mathbf{w}|_{H^{1}(\Omega)}^{2}\lesssim\left(\|\epsilon(\mathbf{w})\|_{L^{2}(\Omega)}^{2}+\sum_{k=1}^{K}\left(\int_{0}^{1}(2y-1)\,w_{1}(k,y)\,\mathrm{d}y\right)^{2}\right)$
(19)

for all $\mathbf{w}(\mathbf{x})=(w_{1}(\mathbf{x}),w_{2}(\mathbf{x}))^{\top}$.
Define

$\psi_{k}(x,y):=\begin{cases}((2y-1)x,-x^{2})^{\top}&\mbox{for }x<k,\\\
((2y-1)k,k^{2}-2kx)^{\top}&\mbox{for }x\geq k\end{cases}$

and observe that

$\epsilon(\psi_{k})=\begin{pmatrix}2y-1&0\\\ 0&0\end{pmatrix}\chi_{x<k},$

where $\chi_{x<k}$ is the characteristic function, i.e., it takes the value
$1$ for $x<k$ and vanishes otherwise. Since $\psi_{k}\in\mathbf{V}$, the
Galerkin orthogonality implies

$\int_{\Omega}\epsilon(\mathbf{u}-\mathbf{u}_{h}):\epsilon(\psi_{k})\,\mathrm{d}\mathbf{x}=0\quad\mbox{for
all}\quad k\in\\{1,\ldots,K\\}.$

Let $\mathbf{w}:=\mathbf{u}-\mathbf{u}_{h}=(w_{1},w_{2})^{\top}$. Using the
explicit representation of $\psi_{k}$ and the Dirichlet conditions on
$\Gamma_{D}=\Gamma^{(0)}$, we obtain

$0=\int_{0}^{k}\int_{0}^{1}(2y-1)\partial_{x}w_{1}(x,y)\,\mathrm{d}y\,\mathrm{d}x=\int_{0}^{1}(2y-1)\,w_{1}(k,y)\,\mathrm{d}y.$

In combination with (19), this yields
$|\mathbf{u}-\mathbf{u}_{h}|_{H^{1}(\Omega)}^{2}\lesssim\|\epsilon(\mathbf{u}-\mathbf{u}_{h})\|_{L^{2}(\Omega)}^{2}.$
Using Galerkin orthogonality, we immediately obtain (18). The second results
directly follows using Lemma 2. ∎

We are interested in a fast solver for this problem, which does not degrade
for domains with small Korn constant. For the solution, we aim for a domain
decomposition approach, where the convergence rates only depend on the Korn
constants for the individual patches, which serve as substructures for the
domain decomposition solver; specifically, we require that

$\|\epsilon(\mathbf{v})\|_{L^{2}({\Omega^{(k)}}\raisebox{5.60002pt}{})}\geq\alpha_{k}|\mathbf{v}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}\;\mbox{for
all}\;\mathbf{v}\in[H^{1}({{\Omega^{(k)}}\raisebox{8.00003pt}{}})]^{2}\;\mbox{with}\;\int_{{\Omega^{(k)}}\raisebox{5.60002pt}{}}\text{curl}\,\mathbf{v}\;\mathrm{d}\mathbf{x}=0,$
(20)

which follows from [1, eq. (3.0.1) and (3.1.1)]. Certainly, the constant
$\alpha_{k}$ only depends on the shape of $\Omega^{(k)}$ and is independent on
the number of patches and the shape of the overall domain $\Omega$.

Finally, we give proofs for Lemmas 1 and 2.

Proof of Lemma 1. Let $\mathbf{v}\in[H^{1}(\Omega)]^{2}$ arbitrary but fixed.
[1, eq. (3.0.1) and (3.1.1)] states that there is a constant
$\alpha_{\Omega}>0$ such that

$\|\epsilon(\mathbf{w})\|_{L^{2}(\Omega)}\geq\alpha_{\Omega}|\mathbf{w}|_{H^{1}(\Omega)}\;\mbox{for
all}\;\mathbf{w}\in[H^{1}({\Omega})]^{2}\;\mbox{with}\;\int_{\Omega}\text{curl}\,\mathbf{w}\;\mathrm{d}\mathbf{x}=0.$
(21)

For $\mathbf{r}$ as in (13), we have $\text{curl}\,\mathbf{r}=2$. For
$\rho:=\tfrac{1}{2}\int_{\Omega}\text{curl}\,\mathbf{v}\;\mathrm{d}\mathbf{x}$,
we have $\mathbf{v}=\rho\,\mathbf{r}+\mathbf{w}$, where
$\int_{\Omega}\text{curl}\,\mathbf{w}\;\mathrm{d}\mathbf{x}=0$. Using the
triangle inequality, the linearity of $\ell$ and its boundedness, we obtain

$\displaystyle|\mathbf{v}|_{H^{1}(\Omega)}$
$\displaystyle\leq|\mathbf{w}|_{H^{1}(\Omega)}+|\rho\,\mathbf{r}|_{H^{1}(\Omega)}=|\mathbf{w}|_{H^{1}(\Omega)}+\frac{|\mathbf{r}|_{H^{1}(\Omega)}}{|\ell(\mathbf{r})|}|\ell(\rho\mathbf{r})|$
$\displaystyle\leq|\mathbf{w}|_{H^{1}(\Omega)}+\frac{|\mathbf{r}|_{H^{1}(\Omega)}}{|\ell(\mathbf{r})|}|\ell(\textbf{w})|+\frac{|\mathbf{r}|_{H^{1}(\Omega)}}{|\ell(\mathbf{r})|}|\ell(\textbf{v})|$
$\displaystyle\leq\left(1+\frac{|\mathbf{r}|_{H^{1}(\Omega)}}{|\ell(\mathbf{r})|}\|\ell\|_{[H^{1}(\Omega)]^{*}}\right)|\mathbf{w}|_{H^{1}(\Omega)}+\frac{|\mathbf{r}|_{H^{1}(\Omega)}}{|\ell(\mathbf{r})|}|\ell(\mathbf{v})|,$

where $\|\ell\|_{[H^{1}(\Omega)]^{*}}$ denotes the dual norm of $\ell$. Using
(21), $\epsilon(\mathbf{v}-\mathbf{w})=\epsilon(\rho\mathbf{r})=0$, and
$\alpha_{\Omega}\leq 1$, we further have

$\displaystyle|\mathbf{v}|_{H^{1}(\Omega)}$
$\displaystyle\leq\alpha_{\Omega}^{-1}\left(1+\frac{|\mathbf{r}|_{H^{1}(\Omega)}}{|\ell(\mathbf{r})|}\max\\{1,\|\ell\|_{[H^{1}(\Omega)]^{*}}\\}\right)\left(\|\epsilon(\mathbf{v})\|_{L^{2}(\Omega)}+|\ell(\mathbf{v})|\right),$

which finishes the proof for
$\alpha_{\Omega,\ell}^{-1}=\alpha_{\Omega}\left(1+\frac{|\mathbf{r}|_{H^{1}(\Omega)}}{|\ell(\mathbf{r})|}\max\\{1,\|\ell\|_{[H^{1}(\Omega)]^{*}}\\}\right)^{-1}$.
$\square$

Proof of Lemma 2. This lemma is a slight generalization of [25, Theorem 3.6].
As in [25], we construct the function $\mathbf{v}_{h}$ patchwise. We choose
$\mathbf{v}_{h}\circ\mathbf{G}_{k}:=\widehat{\Pi}_{k}(\mathbf{u}\circ\mathbf{G}_{k})$
as in [25, Eq. (3.11)], where $\widehat{\Pi}_{k}$ is the projector from [25,
Eq. (3.9)]. This method indeed gives a function $\mathbf{v}_{h}\in\mathbf{V}$,
see [25, Lemma 3.4]. For deriving the error estimates themselves, we do not
use the $\mathpzc{p}$-robust estimates from [26], but the more recent
$\mathpzc{p}$-robust estimates from [19, Theorem 3.1]. Together with (7), and
(8), we obtain

$\displaystyle\inf_{\mathbf{v}_{h}\in\mathbf{V}}|\mathbf{u}-\mathbf{v}_{h}|_{H^{1}(\Omega)}^{2}$
$\displaystyle\lesssim\sum_{k=1}^{K}|(I-\widehat{\Pi}_{k})(\mathbf{u}\circ\mathbf{G}_{k})|_{H^{1}}^{2}\leq\pi^{-2q}\sum_{k=1}^{K}\widehat{h}_{k}^{2q}|\mathbf{u}\circ\mathbf{G}_{k}|_{H^{q+1}(\widehat{\Omega})}^{2},$

form which the desired result follows using (15). $\square$

##  4 A IETI-DP solver for the elasticity problem

In this section, we outline the setup of the proposed IETI-DP method for
solving the discretized elasticity problem (10). For the IETI-DP method, we
have to assemble the variational problem locally, i.e., we have

$\displaystyle{a^{(k)}}\raisebox{8.00003pt}{}({\mathbf{u}^{(k)}}\raisebox{8.00003pt}{},{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{})$
$\displaystyle:=2\mu(\epsilon({\mathbf{u}^{(k)}}\raisebox{8.00003pt}{}),\epsilon({\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}))_{L^{2}(\Omega^{(k)})}+\lambda(\text{div}\,{\mathbf{u}^{(k)}}\raisebox{8.00003pt}{},\text{div}\,{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{})_{L^{2}(\Omega^{(k)})}$
$\displaystyle\langle{F^{(k)}}\raisebox{8.00003pt}{},{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}\rangle$
$\displaystyle:=(\mathbf{f},{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{})_{L^{2}(\Omega^{(k)})}+(\mathbf{g}_{N},{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{})_{L^{2}({\Gamma_{N}^{(k)}}\raisebox{5.60002pt}{})},$
(22)

where
${\Gamma_{N}^{(k)}}\raisebox{8.00003pt}{}:=\Gamma_{N}\cap\partial{\Omega^{(k)}}\raisebox{8.00003pt}{}$.
If the problems would not be coupled, we would have
${\mathbf{u}^{(k)}}\raisebox{8.00003pt}{}\in{\mathbf{V}^{(k)}}\raisebox{8.00003pt}{}$
with

${a^{(k)}}\raisebox{8.00003pt}{}({\mathbf{u}^{(k)}}\raisebox{8.00003pt}{},{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{})\stackrel{{\scriptstyle!}}{{=}}\langle{F^{(k)}}\raisebox{8.00003pt}{},{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}\rangle\quad\mbox{for
all}\quad{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}\in{\mathbf{V}^{(k)}}\raisebox{8.00003pt}{}.$

By fixing bases, we obtain linear systems

${A^{(k)}}\raisebox{8.00003pt}{}{\underline{\mathbf{u}}^{(k)}}\raisebox{8.00003pt}{}\stackrel{{\scriptstyle!}}{{=}}\underline{\mathbf{f}}^{(k)}.$
(23)

Here and in what follows, underlined quantities refer to the coefficient
representations of the corresponding functions. We represent the spaces
$\mathbf{V}^{(k)}$ as direct sum
$\mathbf{V}^{(k)}=\mathbf{V}^{(k)}_{\Gamma}\oplus\mathbf{V}^{(k)}_{\mathrm{I}}$,
where $\mathbf{V}^{(k)}_{\mathrm{I}}$ is spanned by the basis functions which
vanish on the interfaces and $\mathbf{V}^{(k)}_{\Gamma}$ is spanned by the
remaining functions, i.e., the functions that are active on the interfaces
(which includes the primal degrees of freedom). Assuming a corresponding
ordering of the basis functions, we have

$\begin{pmatrix}{A_{\Gamma\Gamma}^{(k)}}\raisebox{8.00003pt}{}&{A_{\Gamma\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}\\\
{A_{\mathrm{I}\Gamma}^{(k)}}\raisebox{8.00003pt}{}&{A_{\mathrm{I}\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}\\\
\end{pmatrix}\begin{pmatrix}\underline{\mathbf{u}}_{\Gamma}^{(k)}\\\
\underline{\mathbf{u}}_{\mathrm{I}}^{(k)}\end{pmatrix}\stackrel{{\scriptstyle!}}{{=}}\begin{pmatrix}\underline{\mathbf{f}}_{\Gamma}^{(k)}\\\
\underline{\mathbf{f}}_{\mathrm{I}}^{(k)}\end{pmatrix}.$ (24)

These local systems correspond to pure Neumann problems, unless the
corresponding patch contributes to the Dirichlet boundary. Therefore, these
systems are not uniquely solvable since the constant is in the null space of
${A^{(k)}}\raisebox{8.00003pt}{}$. To ensure that the system matrices of the
patch-local problems are non-singular, we introduce the function values on
each of the corners of the patch are primal degrees of freedom, i.e.,

$\mathbf{u}^{(k)}(\textbf{x}_{j})=\mathbf{u}^{(\ell)}(\textbf{x}_{j})\qquad\mbox{for
all}\qquad j=1,\ldots,J\mbox{ and }k,\ell\in\mathcal{N}_{\textbf{x}}(j).$

The convergence analysis that we develop carries over automatically to the
case where additional primal degrees of freedom are introduced, like the edge
averages of $\mathbf{u}$ or its normal component; we consider such choices in
the numerical experiments presented in Section 6. For the theoretical
considerations, we require only the corner values to be chosen as primal
degrees of freedom.

The corresponding Lagrangian multipliers are denoted by
${\mu^{(k)}}\raisebox{8.00003pt}{}$. We obtain

$\bar{A}^{(k)}\bar{\underline{\mathbf{u}}}^{(k)}:=\begin{pmatrix}{A_{\Gamma\Gamma}^{(k)}}\raisebox{8.00003pt}{}&{A_{\Gamma\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}&{C^{(k)}}\raisebox{8.00003pt}{}^{\top}\\\
{A_{\mathrm{I}\Gamma}^{(k)}}\raisebox{8.00003pt}{}&{A_{\mathrm{I}\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}&0\\\
{C^{(k)}}\raisebox{8.00003pt}{}&0&0\\\
\end{pmatrix}\begin{pmatrix}{\underline{\mathbf{u}}_{\Gamma}^{(k)}}\raisebox{8.00003pt}{}\\\
{\underline{\mathbf{u}}_{\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}\\\
{\underline{\mu}^{(k)}}\raisebox{8.00003pt}{}\end{pmatrix}\stackrel{{\scriptstyle!}}{{=}}\begin{pmatrix}{\underline{\mathbf{f}}_{\Gamma}^{(k)}}\raisebox{8.00003pt}{}\\\
{\underline{\mathbf{f}}_{\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}\\\
0\end{pmatrix}=:\bar{\underline{\mathbf{f}}}^{(k)},$ (25)

where ${C^{(k)}}\raisebox{8.00003pt}{}$ evaluates the primal degrees of
freedom associated to the patch ${\Omega^{(k)}}\raisebox{8.00003pt}{}$, such
that
${C^{(k)}}\raisebox{8.00003pt}{}{\underline{\mathbf{u}}_{\Gamma}^{(k)}}\raisebox{8.00003pt}{}=0$.
The continuity between the patches is enforced by the matrices
${B^{(k)}}\raisebox{8.00003pt}{}$. The term
$\sum_{k=1}^{K}{B^{(k)}}\raisebox{8.00003pt}{}\underline{\mathbf{u}}_{\Gamma}^{(k)}$
evaluates to a vector containing the differences of the coefficients of any
two matching basis functions. Here, we do not include the vertex values of
(see Figure 2) since these are primal degrees of freedom. The relation

$\sum_{k=1}^{K}{B^{(k)}}\raisebox{8.00003pt}{}\underline{\mathbf{u}}_{\Gamma}^{(k)}=0$

guarantees the continuity across the patches.

$\Omega^{(1)}$$\Omega^{(3)}$$\Omega^{(2)}$$\Omega^{(4)}$ Figure 2: Enforcing
continuity across the patches. The corners are excluded.

Moreover, we introduce the primal problem, i.e., the global problem for the
primal degrees of freedom. We use a $A^{(k)}$-orthogonal basis for the primal
degrees of freedom. This basis is represented in terms of the basis functions
of the basis for
$\mathbf{V}_{\Gamma}^{(k)}\times\mathbf{V}_{\mathrm{I}}^{(k)}$ using the
matrix $\Psi^{(k)}$, which are the solution of the system

$\begin{pmatrix}{A_{\Gamma\Gamma}^{(k)}}\raisebox{8.00003pt}{}&{A_{\Gamma\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}&{C^{(k)}}\raisebox{8.00003pt}{}^{\top}\\\
{A_{\mathrm{I}\Gamma}^{(k)}}\raisebox{8.00003pt}{}&{A_{\mathrm{I}\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}&0\\\
{C^{(k)}}\raisebox{8.00003pt}{}&0&0\\\
\end{pmatrix}\begin{pmatrix}{\Psi_{\Gamma}^{(k)}}\raisebox{8.00003pt}{}\\\
{\Psi_{\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}\\\
{\Phi^{(k)}}\raisebox{8.00003pt}{}\\\ \end{pmatrix}=\begin{pmatrix}0\\\ 0\\\
{R^{(k)}}\raisebox{8.00003pt}{}\end{pmatrix},$ (26)

where ${R^{(k)}}\raisebox{8.00003pt}{}$ is a boolean matrix that selects the
primal degrees of freedom which are active on the patch $\Omega^{(k)}$. We
have ${R^{(k)}}\raisebox{8.00003pt}{}\in\mathbb{R}^{N^{(k)}_{\Pi}\times
N_{\Pi}}$, where $N^{(k)}_{\Pi}$ is the number of primal degrees of freedom
corresponding to the $k$’th patch (thus $N^{(k)}_{\Pi}=8$ if the patch does
not contribute to the Dirichlet boundary) and $N_{\Pi}$ is the overall number
of primal degrees of freedom.

We define the system matrix, right-hand side and jump matrix for the primal
problem as

$A_{\Pi}:=\sum^{K}_{k=1}{\Psi^{(k)}}\raisebox{8.00003pt}{}^{\top}{A^{(k)}}\raisebox{8.00003pt}{}{\Psi^{(k)}}\raisebox{8.00003pt}{},\quad\underline{\mathbf{f}}_{\Pi}:=\sum^{K}_{k=1}{\Psi^{(k)}}\raisebox{8.00003pt}{}^{\top}{\underline{\bar{\mathbf{f}}}^{(k)}}\raisebox{8.00003pt}{}\quad\text{and}\quad
B_{\Pi}:=\sum^{K}_{k=1}{B^{(k)}}\raisebox{8.00003pt}{}{\Psi^{(k)}}\raisebox{8.00003pt}{}.$

Finally we are able to write down the overall IETI-DP system:

$\begin{pmatrix}{\bar{A}^{(1)}}\raisebox{8.00003pt}{}&&&&{\bar{B}^{(1)}}\raisebox{8.00003pt}{}^{\top}\\\
&\ddots&&&\vdots\\\
&&{\bar{A}^{(K)}}\raisebox{8.00003pt}{}&&{\bar{B}^{(K)}}\raisebox{8.00003pt}{}^{\top}\\\
&&&{A}_{\Pi}&{B}_{\Pi}^{\top}\\\
{\bar{B}^{(1)}}\raisebox{8.00003pt}{}&\ldots&{\bar{B}^{(K)}}\raisebox{8.00003pt}{}&{B}_{\Pi}\\\
\end{pmatrix}\begin{pmatrix}{\underline{\bar{\mathbf{x}}}^{(1)}}\raisebox{8.00003pt}{}\\\
\vdots\\\ {\underline{\bar{\mathbf{x}}}^{(K)}}\raisebox{8.00003pt}{}\\\
\underline{{\mathbf{x}}}_{\Pi}\\\
\underline{\lambda}\end{pmatrix}=\begin{pmatrix}{\underline{\bar{\mathbf{f}}}^{(1)}}\raisebox{8.00003pt}{}\\\
\vdots\\\ {\underline{\bar{\mathbf{f}}}^{(K)}}\raisebox{8.00003pt}{}\\\
\underline{\mathbf{f}}_{\Pi}\\\ 0\end{pmatrix}.$ (27)

For solving this linear system, we take the Schur complement with respect to
the Lagrange multipliers $\underline{\lambda}$. This means that we solve

$F\underline{\lambda}=\underline{g},$ (28)

where

$F:={B}_{\Pi}{A}^{-1}_{\Pi}{B}_{\Pi}^{\top}+\sum^{K}_{k=1}{\bar{B}^{(k)}}\raisebox{8.00003pt}{}{\bar{A}^{(k)}}\raisebox{8.00003pt}{}^{-1}{\bar{B}^{(k)}}\raisebox{8.00003pt}{}^{\top}\;\mbox{and}\;\underline{g}:={B}_{\Pi}{A}^{-1}_{\Pi}{\underline{\mathbf{f}}}_{\Pi}+\sum^{K}_{k=1}\bar{B}^{(k)}{\bar{A}^{(k)}}\raisebox{8.00003pt}{}^{-1}\bar{\underline{\mathbf{f}}}^{(k)}$

The patch-local solutions are then recovered by

$\underline{{\mathbf{x}}}_{\Pi}={A}_{\Pi}^{-1}\left(\underline{{\mathbf{f}}}_{\Pi}-B_{\Pi}^{\top}\underline{\lambda}\right)\quad\mbox{and}\quad{\underline{\bar{\mathbf{x}}}^{(k)}}\raisebox{8.00003pt}{}={\bar{A}^{(k)}}\raisebox{8.00003pt}{}^{-1}\left(\underline{\bar{\mathbf{f}}}^{(k)}-{\bar{B}^{(k)}}\raisebox{8.00003pt}{}^{\top}\underline{\lambda}\right).$

The final solution is then obtained by distributing
$\underline{{\mathbf{x}}}_{\Pi}$ to the patches.

We solve the system (28) with a conjugate gradient solver, preconditioned with
the scaled Dirichlet preconditioner

$M_{\mathrm{sD}}:=\sum_{k=1}^{K}{B^{(k)}}\raisebox{8.00003pt}{}\mathcal{D}_{k}^{-1}{S_{A}^{(k)}}\raisebox{8.00003pt}{}\mathcal{D}_{k}^{-1}{B^{(k)}}\raisebox{8.00003pt}{}^{\top},$

where $\mathcal{D}_{k}:=2I$ is set up based on the principle of multiplicity
scaling and
${S_{A}^{(k)}}\raisebox{8.00003pt}{}:={A_{\Gamma\Gamma}^{(k)}}\raisebox{8.00003pt}{}-{A_{\Gamma\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}{A_{\mathrm{I}\mathrm{I}}^{(k)}}\raisebox{8.00003pt}{}^{-1}{A_{\mathrm{I}\Gamma}^{(k)}}\raisebox{8.00003pt}{}$
is the local Schur complement based on the splitting (24).

######  Remark 3.

In this paper we focus on solving (28) for the symmetric positive definite
matrix $F$. This requires the exact representation of the inverses of
${\bar{A}^{(k)}}\raisebox{8.00003pt}{}$. If these are unavailable, then
efficient approximations of the inverses can be used by instead solving (27).
For more information see, e.g., [22, 17].

##  5 Condition number analysis for the IETI solver

The convergence rates of the conjugate gradient solver can be estimated based
on the condition number of the preconditioned system $M_{\mathrm{sD}}\bar{F}$,
based on the framework introduced in [16] and based on the results from [20],
where a IETI-DP solver for the Poisson problem has been analyzed in detail.

In this paper, we are interested in a IETI-DP solver for the linearized
elasticity problem. The extension of the results from [20] to the vector-
valued is not an issue since all proofs could be reused verbatim. The main
difference is that we have to show the coercivity of the bilinear form
associated to the linearized elasticity problem, which again requires the use
of the Korn inequality. In general, not every patch touches the Dirichlet
boundary $\Gamma_{D}$. Thus, we need a framework that allows the application
of the Korn inequality.

Let $\mathbf{V}_{d}:=\mathbf{V}^{(1)}\times\cdots\times\mathbf{V}^{(K)}$. The
space $\widetilde{\mathbf{V}}\subseteq\mathbf{V}_{d}$ models all functions
whose corner values are continuous:

$\widetilde{\mathbf{V}}:=\left\\{\mathbf{v}\in\mathbf{V}_{d}:\mathbf{v}^{(k)}(\textbf{x}_{j})=\mathbf{v}^{(\ell)}(\textbf{x}_{j})\mbox{
for all }k,\ell\in\mathcal{N}_{\textbf{x}}(j)\mbox{ with
}j\in\\{1,\ldots,J\\}\right\\},$

where $\mathbf{v}=(\mathbf{v}^{(1)},\cdots,\mathbf{v}^{(K)})$.

In the following, we show that the $H^{1}$-seminorm can be bounded from above
by the semi-norm
$\|\epsilon(\cdot)\|_{L^{2}({\Omega^{(k)}}\raisebox{5.60002pt}{})}$ if we
guarantee that the corner values vanish. This, of course, is only possible for
spaces, where the notation of the corner values is well-defined, like the
discrete space ${\mathbf{V}^{(k)}}\raisebox{8.00003pt}{}$. Let
$\mathcal{Q}_{k}:=\\{\mathbf{v}:\mathbf{v}\circ\mathbf{G}_{k}\mbox{ is
bilinear}\\}$ be the space of all isogeometric functions with bilinear pre-
image. Since a bilinear function is uniquely given by its values at the
corners of the unit square $\widehat{\Omega}$, the following projector is
well-defined:

${\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}:{\mathbf{V}^{(k)}}\raisebox{8.00003pt}{}\to\mathcal{Q}_{k},\quad\mbox{such
that}\quad({\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{})(\textbf{x}_{\ell})={\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}(\textbf{x}_{\ell})\quad\mbox{for
all}\quad\ell\in\mathcal{N}_{\textbf{x}}^{-1}(k).$

The following lemma gives upper bounds both for the $H^{1}$-seminorm and the
$L^{\infty}_{0}$-seminorm, given by

$|{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}|_{L^{\infty}_{0}({\Omega^{(k)}}\raisebox{5.60002pt}{})}:=\inf_{\mu_{k}\in\mathbb{R}^{2}}\|{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-\mu_{k}\|_{L^{\infty}({\Omega^{(k)}}\raisebox{5.60002pt}{})},$

of
${\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-{\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}$,
which vanishes on the corners by construction.

######  Lemma 5.

For all $k=1,\ldots,K$, the inequalities

$\displaystyle|{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-{\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}$
$\displaystyle\lesssim\alpha_{k}^{-1}\Lambda_{k}^{1/2}\|\epsilon({\mathbf{v}^{(k)}}\raisebox{8.00003pt}{})\|_{L^{2}({\Omega^{(k)}}\raisebox{5.60002pt}{})},$
$\displaystyle|{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-{\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}|_{L_{0}^{\infty}({\Omega^{(k)}}\raisebox{5.60002pt}{})}$
$\displaystyle\lesssim\alpha_{k}^{-1}\Lambda_{k}^{1/2}\|\epsilon({\mathbf{v}^{(k)}}\raisebox{8.00003pt}{})\|_{L^{2}({\Omega^{(k)}}\raisebox{5.60002pt}{})}$

hold for all $\mathbf{v}\in\widetilde{\mathbf{V}}$, where $\alpha_{k}$ is the
constant from (20), $\Lambda_{k}:=1+\log\mathpzc{p}+\log\frac{H_{k}}{h_{k}}$.

###### Proof.

For any $\mathbf{v}\in\widetilde{\mathbf{V}}$, there is a decomposition
${\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}={\mathbf{w}^{(k)}}\raisebox{8.00003pt}{}+\rho_{k}\mathbf{r}$,
with
$\int_{{\Omega^{(k)}}\raisebox{5.60002pt}{}}\text{curl}\,{\mathbf{w}^{(k)}}\raisebox{8.00003pt}{}\,\mathrm{d}x=0$,
$\mathbf{r}$ as in (13) and $\rho_{k}\in\mathbb{R}$. Using the triangle
inequality, we have

$|{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-{\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}\leq|{\mathbf{w}^{(k)}}\raisebox{8.00003pt}{}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}+|{\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-\rho_{k}\mathbf{r}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}.$
(29)

Since both
${\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}$
and $\mathbf{r}$ are bilinear functions, the term
${\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-\rho_{k}\mathbf{r}$
is uniquely determined by the function values at the corners. Since only the
$H^{1}$-seminorm is considered, only deviations from the average matter; thus

$|{\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-\rho_{k}\mathbf{r}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}\lesssim\psi_{k}:=\inf_{\mu_{k}\in\mathbb{R}^{2}}\sum_{\ell\in\mathcal{N}_{\mathbf{x}}^{-1}(k)}|({\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{})(\mathbf{x}_{\ell})-\rho_{k}\mathbf{r}(\mathbf{x}_{\ell})-\mu_{k}|.$
(30)

Since
${\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}$
and ${\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}$ coincide at the corners, we
obtain

$\psi_{k}=\inf_{\mu_{k}\in\mathbb{R}^{2}}\sum_{\ell\in\mathcal{N}_{\mathbf{x}}^{-1}(k)}|{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}(\mathbf{x}_{\ell})-\rho_{k}\mathbf{r}(\mathbf{x}_{\ell})-\mu_{k}|=\inf_{\mu_{k}\in\mathbb{R}^{2}}\sum_{\ell\in\mathcal{N}_{\mathbf{x}}^{-1}(k)}|{\mathbf{w}^{(k)}}\raisebox{8.00003pt}{}(\mathbf{x}_{\ell})-\mu_{k}|.$

Using [20, Lemma 4.14] and the usual relation between parameter domain and
physical domain and the Poincaré–Friedrichs inequality, we obtain further

$\psi_{k}\lesssim\Lambda_{k}^{1/2}\inf_{\mu_{k}\in\mathbb{R}^{2}}\|{\mathbf{w}^{(k)}}\raisebox{8.00003pt}{}-\mu_{k}\|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}\\\
\lesssim\Lambda_{k}^{1/2}|{\mathbf{w}^{(k)}}\raisebox{8.00003pt}{}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}.$
(31)

Together with (30) and (29), this yields

$\displaystyle|{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-{\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}$
$\displaystyle\lesssim\Lambda_{k}^{1/2}|{\mathbf{w}^{(k)}}\raisebox{8.00003pt}{}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}.$

Finally, (20) gives the desired first estimate. For the second result, use
again the triangle inequality, the fact that the $L_{0}^{\infty}$-seminorm of
the bilinear function
${\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-\rho_{k}\mathbf{r}$
can be bounded by its values at the corner, [20, Lemma 4.14], and (30) and
(31) yield

$\displaystyle|{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}-{\Pi_{F}^{(k)}}\raisebox{8.00003pt}{}{\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}|_{L_{0}^{\infty}({\Omega^{(k)}}\raisebox{5.60002pt}{})}\lesssim|{\mathbf{w}^{(k)}}\raisebox{8.00003pt}{}|_{L_{0}^{\infty}({\Omega^{(k)}}\raisebox{5.60002pt}{})}+\psi_{k}\lesssim\Lambda_{k}^{1/2}|{\mathbf{w}^{(k)}}\raisebox{8.00003pt}{}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})},$

from which the second estimate follows using (20). ∎

We build our main condition number estimate on the framework from [16]. To
obtain a concise notation, let
$\mathbf{W}_{d}:=\mathbf{W}^{(1)}\times\cdots\times\mathbf{W}^{(K)}$ with
$\mathbf{W}^{(k)}:=\\{\mathbf{v}|_{\partial\Omega^{(k)}}\,:\,\mathbf{v}\in\mathbf{V}^{(k)}\\}$
be the skeleton space and let

$\widetilde{\mathbf{W}}:=\left\\{\mathbf{w}\in\mathbf{W}_{d}:\mathbf{w}^{(k)}(\textbf{x}_{j})=\mathbf{w}^{(\ell)}(\textbf{x}_{j})\;\forall\;k,\ell\in\mathcal{N}_{\textbf{x}}(j)\mbox{
with }j\in\\{1,\ldots,J\\}\right\\},$

where $\mathbf{w}=(\mathbf{w}^{(1)},\cdots,\mathbf{w}^{(K)})$, be the
corresponding space which is continuous at the corners. Here and in what
follows, for any function, say ${\mathbf{w}^{(k)}}\raisebox{8.00003pt}{}$ in
${\mathbf{W}^{(k)}}\raisebox{8.00003pt}{}$, the corresponding underlined
symbol, here ${\underline{\mathbf{w}}^{(k)}}\raisebox{8.00003pt}{}$, denotes
the representation of the corresponding functions with respect to the basis
for the space ${\mathbf{W}^{(k)}}\raisebox{8.00003pt}{}$. The overall vector
$\underline{\mathbf{w}}$ is a block-vector containing the local contributions
${\underline{\mathbf{w}}^{(k)}}\raisebox{8.00003pt}{}$. Analogously, $B$ is
block-row matrix containing the local matrices
${B^{(k)}}\raisebox{8.00003pt}{}$ and $\mathcal{D}$ and $S_{A}$ are block-
diagonal matrices containing the local matrices $\mathcal{D}_{k}$ and
${S_{A}^{(k)}}\raisebox{8.00003pt}{}$, respectively. [16, Theorem 22] states
that the condition number bound depends on the maximum value of the ratio

$\frac{\|\mathcal{D}^{-1}B^{\top}B\underline{\mathbf{v}}\|_{S_{A}}^{2}}{\|\underline{\mathbf{v}}\|_{S_{A}}^{2}},$
(32)

where the vector $\underline{\mathbf{v}}$ is a coefficient representation of a
function $\mathbf{v}\in\widetilde{\mathbf{W}}$. The following Lemma gives an
upper bound for that ratio.

######  Lemma 6.

For all $\mathbf{v}\in\widetilde{\mathbf{W}}$ with coefficient representation
$\underline{\mathbf{v}}$, the estimate

$\|\mathcal{D}^{-1}B^{\top}B\underline{\mathbf{v}}\|_{S_{A}}^{2}\lesssim\frac{\mu+\lambda}{\mu}\,\mathpzc{p}\left(1+\log\mathpzc{p}+\max_{k=1,\ldots,K}\log\frac{H_{k}}{h_{k}}\right)^{2}\left(\max_{k=1,\ldots,K}\alpha_{k}^{-2}\right)\|\underline{\mathbf{v}}\|_{S_{A}}^{2}$

holds.

###### Proof.

Let $\mathbf{w},\mathbf{v}\in\widetilde{\mathbf{W}}$ with coefficient
representations $\underline{\mathbf{w}}$ and $\underline{\mathbf{v}}$ be such
that
$\underline{\mathbf{w}}:=\mathcal{D}^{-1}B^{\top}B\underline{\mathbf{v}}$.
Since the Schur-complement is the energy-minimizing extension and using the
boundedness of $A$, we have
$\|\mathcal{D}^{-1}B^{\top}B\underline{\mathbf{v}}\|_{S_{A}}^{2}=\|\underline{\mathbf{w}}\|_{S_{A}}^{2}$
with

$\displaystyle\|\underline{\mathbf{w}}\|_{S_{A}}^{2}$
$\displaystyle=\sum_{k=1}^{K}\inf_{\begin{array}[]{c}\mbox{\scriptsize${\mathbf{u}^{(k)}}\raisebox{5.60002pt}{}\in{\mathbf{V}^{(k)}}\raisebox{5.60002pt}{}$}\\\
\mbox{\scriptsize$\left.{\mathbf{u}^{(k)}}\raisebox{5.60002pt}{}\right|_{\partial{\Omega^{(k)}}\raisebox{3.92001pt}{}}={\mathbf{w}^{(k)}}\raisebox{5.60002pt}{}$}\end{array}}2\mu\|\epsilon({\mathbf{u}^{(k)}}\raisebox{8.00003pt}{})\|_{L^{2}({\Omega^{(k)}}\raisebox{5.60002pt}{})}^{2}+\lambda\|\text{div}\,{\mathbf{u}^{(k)}}\raisebox{8.00003pt}{}\|_{L^{2}({\Omega^{(k)}}\raisebox{5.60002pt}{})}^{2}$
$\displaystyle\leq
2(\mu+\lambda)\sum_{k=1}^{K}\inf_{\begin{array}[]{c}\mbox{\scriptsize${\mathbf{u}^{(k)}}\raisebox{5.60002pt}{}\in{\mathbf{V}^{(k)}}\raisebox{5.60002pt}{}$}\\\
\mbox{\scriptsize$\left.{\mathbf{u}^{(k)}}\raisebox{5.60002pt}{}\right|_{\partial{\Omega^{(k)}}\raisebox{3.92001pt}{}}={\mathbf{w}^{(k)}}\raisebox{5.60002pt}{}$}\end{array}}|{\mathbf{u}^{(k)}}\raisebox{8.00003pt}{}|_{H^{1}({\Omega^{(k)}}\raisebox{5.60002pt}{})}^{2}.$

This term is the $H^{1}$-energy of the discrete harmonic extension
$\mathcal{H}^{(k)}:\mathbf{W}^{(k)}\rightarrow\mathbf{V}^{(k)}$, which by
definition minimizes
$|\mathcal{H}^{(k)}\mathbf{w}^{(k)}|_{H^{1}(\Omega^{(k)})}$ under the
constraint
$\left(\mathcal{H}^{(k)}\mathbf{w}^{(k)}\right)\Big{|}_{\partial\Omega^{(k)}}=\mathbf{w}^{(k)}$.
This means that we have

$\|\underline{\mathbf{w}}\|_{S_{A}}^{2}\leq
2(\mu+\lambda)\sum_{k=1}^{K}|\mathcal{H}^{(k)}\mathbf{w}^{(k)}|_{H^{1}(\Omega^{(k)})}^{2}.$

Using [20, Theorem 4.2 and Lemma 4.15], we further have

$\displaystyle\|\underline{\mathbf{w}}\|_{S_{A}}^{2}$ $\displaystyle\lesssim
2(\mu+\lambda)\mathpzc{p}\sum_{k=1}^{K}|\mathbf{w}^{(k)}|_{H^{1/2}(\partial\Omega^{(k)})}^{2}$
$\displaystyle\lesssim
2(\mu+\lambda)\mathpzc{p}\sum_{k=1}^{K}\sum_{\ell\in\mathcal{N}_{\Gamma}(k)}\left(|\mathbf{w}^{(k)}|_{H^{1/2}(\Gamma^{(k,\ell)})}^{2}+\Lambda|\mathbf{w}^{(k)}|_{L_{0}^{\infty}(\Gamma^{(k,\ell)})}^{2}\right),$

where
$\Lambda:=\left(1+\log\mathpzc{p}+\max_{k=1,\ldots,K}\log\frac{H_{k}}{h_{k}}\right)$.
Since $\mathbf{v}^{(k)}$ is continuous at the corners, standard arguments
yield

$\mathbf{w}^{(k)}|_{\Gamma^{(k,\ell)}}=\mathbf{v}^{(k)}|_{\Gamma^{(k,\ell)}}-\mathbf{v}^{(\ell)}|_{\Gamma^{(k,\ell)}},$

see, e.g., [20, Proof of Lemma 4.16] for a proof. Using this statement, we
further have

$\|\underline{\mathbf{w}}\|_{S_{A}}^{2}\lesssim
2(\mu+\lambda)\mathpzc{p}\Lambda\sum_{k=1}^{K}\sum_{\ell\in\mathcal{N}_{\Gamma}(k)}\left(|\mathbf{v}^{(k)}|_{H^{1/2}(\Gamma^{(k,\ell)})}^{2}+|\mathbf{v}^{(k)}|_{L_{0}^{\infty}(\Gamma^{(k,\ell)})}^{2}\right).$

We substitute
${\mathbf{v}_{0}^{(k)}}\raisebox{8.00003pt}{}:=(I-\Pi_{F}^{(k)}){\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}$
for ${\mathbf{v}^{(k)}}\raisebox{8.00003pt}{}$. Since
$\Pi_{F}^{(k)}\mathbf{v}^{(k)}|_{\Gamma^{(k,\ell)}}=\Pi_{F}^{(\ell)}\mathbf{v}^{(\ell)}|_{\Gamma^{(k,\ell)}}$,
we know that $\mathbf{v}-\mathbf{v}_{0}$ vanishes on all interfaces. Thus,
$B\underline{\mathbf{v}}_{0}=B\underline{\mathbf{v}}=0$ and consequently

$\|\underline{\mathbf{w}}\|_{S_{A}}^{2}\lesssim\tilde{\Lambda}\sum_{k=1}^{K}\sum_{\ell\in\mathcal{N}_{\Gamma}(k)}\left(|(I-\Pi_{F}^{(k)})\mathbf{v}^{(k)}|_{H^{1/2}(\Gamma^{(k,\ell)})}^{2}+|(I-\Pi_{F}^{(k)})\mathbf{v}^{(k)}|_{L_{0}^{\infty}(\Gamma^{(k,\ell)})}^{2}\right),$

where $\tilde{\Lambda}:=2(\mu+\lambda)\mathpzc{p}\Lambda$. Using [20, Theorem
4.2] and the fact that the $\mathcal{H}^{(k)}$ is the extension operator that
minimizes the $H^{1}$-seminorm, we obtain

$\displaystyle\|\underline{\mathbf{w}}\|_{S_{A}}^{2}$
$\displaystyle\lesssim\tilde{\Lambda}\sum_{k=1}^{K}\left(|\mathcal{H}^{(k)}(I-\Pi_{F}^{(k)})\mathbf{v}^{(k)}|_{H^{1}(\Omega^{(k)})}^{2}+\sum_{\ell\in\mathcal{N}_{\Gamma}(k)}|(I-\Pi_{F}^{(k)})\mathbf{v}^{(k)}|_{L_{0}^{\infty}(\Gamma^{(k,\ell)})}^{2}\right)$
$\displaystyle\leq\tilde{\Lambda}\sum_{k=1}^{K}\left(|\mathcal{H}_{\epsilon}^{(k)}(I-\Pi_{F}^{(k)})\mathbf{v}^{(k)}|_{H^{1}(\Omega^{(k)})}^{2}+\sum_{\ell\in\mathcal{N}_{\Gamma}(k)}|(I-\Pi_{F}^{(k)})\mathbf{v}^{(k)}|_{L_{0}^{\infty}(\Gamma^{(k,\ell)})}^{2}\right),$

where
$\mathcal{H}_{\epsilon}^{(k)}:\mathbf{W}^{(k)}\rightarrow\mathbf{V}^{(k)}$,
which minimizes the energy
$\|\epsilon(\mathcal{H}_{\epsilon}^{(k)}\mathbf{w}^{(k)})\|_{L^{2}(\Omega^{(k)})}$
under the constraint
$\left(\mathcal{H}^{(k)}_{\epsilon}\mathbf{w}^{(k)}\right)\Big{|}_{\partial\Omega^{(k)}}=\mathbf{w}^{(k)}$.
We use the fact that $\mathcal{H}^{(k)}_{\epsilon}$ does not change the values
at the interfaces and
$\Pi_{F}^{(k)}\mathcal{H}^{(k)}_{\epsilon}=\Pi_{F}^{(k)}$ to obtain

$\displaystyle\|\underline{\mathbf{w}}\|_{S_{A}}^{2}$
$\displaystyle\lesssim\tilde{\Lambda}\sum_{k=1}^{K}\left(|(\mathcal{H}_{\epsilon}^{(k)}-\Pi_{F}^{(k)})\mathbf{v}^{(k)}|_{H^{1}(\Omega^{(k)})}^{2}+|(\mathcal{H}_{\epsilon}^{(k)}-\Pi_{F}^{(k)})\mathbf{v}^{(k)}|_{L_{0}^{\infty}(\Omega^{(k)})}^{2}\right)$
$\displaystyle\lesssim\tilde{\Lambda}\sum_{k=1}^{K}\left(|(I-\Pi_{F}^{(k)})\mathcal{H}_{\epsilon}^{(k)}\mathbf{v}^{(k)}|_{H^{1}(\Omega^{(k)})}^{2}+|(I-\Pi_{F}^{(k)})\mathcal{H}_{\epsilon}^{(k)}\mathbf{v}^{(k)}|_{L_{0}^{\infty}(\Omega^{(k)})}^{2}\right).$

Since $(I-\Pi_{F}^{(k)})\mathbf{v}^{(k)}$ vanishes at the corners, Lemma 5
gives

$\displaystyle\|\underline{\mathbf{w}}\|_{S_{A}}^{2}=\|\mathcal{D}^{-1}B^{\top}B\underline{\mathbf{v}}\|_{S_{A}}^{2}$
$\displaystyle\lesssim
2(\mu+\lambda)\mathpzc{p}\,\Lambda^{2}\,\sum_{k=1}^{K}\alpha_{k}^{-2}\left(\|\epsilon(\mathcal{H}_{\epsilon}^{(k)}\mathbf{v}^{(k)})\|_{L^{2}(\Omega^{(k)})}^{2}+\|\epsilon(\mathcal{H}_{\epsilon}^{(k)}\mathbf{v}^{(k)})\|_{L^{2}(\Omega^{(k)})}^{2}\right)$
$\displaystyle\lesssim
2(\mu+\lambda)\mathpzc{p}\,\Lambda^{2}\,\left(\max_{k=1,\ldots,K}\alpha_{k}^{-2}\right)\sum_{k=1}^{K}\|\epsilon(\mathcal{H}_{\epsilon}^{(k)}\mathbf{v}^{(k)})\|_{L^{2}(\Omega^{(k)})}^{2}$
$\displaystyle\leq\frac{\mu+\lambda}{\mu}\mathpzc{p}\,\Lambda^{2}\left(\max_{k=1,\ldots,K}\alpha_{k}^{-2}\right)\,\|\underline{\mathbf{v}}\|_{S_{A}}^{2},$

which finishes the proof. ∎

The primal elasticity formulation considered in this paper is an elliptic
problem, which directly fits into the framework of [16]. Altough, a scalar
valued problem was considered there, the main results carry over directly to
vector-valued problems. So, we can directly apply [16, Theorem 22] to obtain
the following condition number bound.

######  Theorem 7.

Provided that the IETI-DP solver is set up as outlined in the previous
section, the condition number of the preconditioned system satisfies

$\kappa_{\mathrm{ess}}(M_{\mathrm{sD}}F)\leq
c\,\frac{\mu+\lambda}{\mu}\,\mathpzc{p}\left(1+\log\mathpzc{p}+\max_{k=1,\ldots,K}\log\frac{H_{k}}{h_{k}}\right)^{2}\,\left(\max_{k=1,\ldots,K}\alpha_{k}^{-2}\right)\,,$

where $c$ is a constant that only depends on the constants from the
assumptions (5), (7) and (9) (cf. Notation 1), $\alpha_{k}$ are the patch-
local Korn constants from (20) and $\kappa_{\mathrm{ess}}$ denotes the
essential condition number, which is the ratio between the largest eigenvalue
and the smallest non-zero eigenvalue.

###### Proof.

[16, Theorem 22] states that the bound follows directly from Lemma 6. ∎

If standard Krylov space methods are applied to the singular matrix $F$,
preconditioned with a non-singular preconditioner $M_{\mathrm{sD}}$, all
iterations live in the corresponding factor space. The convergence behavior is
dictated by the essential condition number of $M_{\mathrm{sD}}F$, cf. [16,
Remark 23].

##  6 Numerical results

In this section, we present numerical results that illustrate the efficiency
of the proposed IETI-DP solver and the address the issues related to geometry
locking. We consider two domains, a tortoise (Subsection 6.1) and a wrench
(Subsection 6.2). The first domain illustrates that the IETI-DP solver works
also well for domains with non-trivial parmeterizations $\mathbf{G}_{k}$ and
that contain extraordinary vertices (interior vertices touching less than or
more than four patches). We consider the wrench domain as a non-trivial
example of a long and thin domain, susceptible for geometry locking.

###  6.1 Results for the tortoise domain

We consider a tortoise domain, composed of 60 patches, see Figure 3. Each of
the patches is parameterized using B-splines. We consider the elasticity
problem with Lamé parameters $\mu=\lambda=1$ and without external volumetric
forces ($\mathbf{f}=0$). The boundary conditions are chosen as follows. On the
bottom of the feet (depicted in blue color in Figure 3), we have homogeneous
Dirichlet boundary conditions ($\mathbf{g}_{D}=0$). On the top of the tortoise
(depicted in red color in Figure 3), we assume inhomogeneous Neumann boundary
conditions with $\mathbf{g}_{N}=(0,-g_{0})^{\top}$, with $g_{0}=0.00001$. On
the remaining boundary, we set homogeneous Neumann boundary conditions.

Figure 3: A tortoise

The discretization is done as outlined in Section 2. On each patch, we obtain
the discretization space by performing $\ell=1,2,\ldots$ uniform refinement
steps, staring from patchwise global polynomials, i.e., a grid without inner
knots. For solving the overall problem, we use the IETI-DP method as proposed
in Section 4. The reduced system (28) is solved using a conjugate gradient
solver. We use a random initial guess. The overall solver is stopped when the
Euclidean norm of the residual is reduced by a factor of
$\varepsilon=10^{-6}$, compared to the Euclidean norm of the initial residual.
The condition numbers are estimated using the underlying Lanczos iteration of
the conjugate gradient method. The local linear systems are solved using a
sparse LU-solver. All experiments have been implemented using the G+Smo
library111https://github.com/gismo/gismo.

| $\mathpzc{p}=2$ | $\mathpzc{p}=3$ | $\mathpzc{p}=4$ | $\mathpzc{p}=5$ | $\mathpzc{p}=6$
---|---|---|---|---|---
$\ell$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$
$2$ |  $20$ |  $12.2$ |  $23$ |  $16.6$ |  $25$ |  $20.0$ |  $26$ |  $22.6$ |  $26$ |  $25.2$
$3$ |  $24$ |  $17.2$ |  $27$ |  $22.2$ |  $29$ |  $26.1$ |  $30$ |  $29.3$ |  $31$ |  $32.4$
$4$ |  $28$ |  $22.7$ |  $31$ |  $28.6$ |  $33$ |  $33.3$ |  $35$ |  $37.6$ |  $35$ |  $42.0$
$5$ |  $31$ |  $29.2$ |  $34$ |  $36.6$ |  $37$ |  $42.8$ |  $39$ |  $48.8$ |  $39$ |  $53.6$
Table 1: Iteration counts (it) and condition numbers ($\kappa$); tortoise; $\mu=1,\lambda=1$.  | $\mathpzc{p}=2$ | $\mathpzc{p}=3$ | $\mathpzc{p}=4$ | $\mathpzc{p}=5$ | $\mathpzc{p}=6$
---|---|---|---|---|---
$\ell$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$
$2$ |  $13$ |  $4.3$ |  $16$ |  $5.9$ |  $18$ |  $7.6$ |  $18$ |  $9.4$ |  $19$ |  $11.2$
$3$ |  $17$ |  $6.2$ |  $20$ |  $8.7$ |  $21$ |  $11.4$ |  $22$ |  $14.0$ |  $22$ |  $16.2$
$4$ |  $20$ |  $8.7$ |  $23$ |  $12.7$ |  $24$ |  $16.4$ |  $25$ |  $19.5$ |  $26$ |  $22.2$
$5$ |  $23$ |  $12.4$ |  $26$ |  $17.9$ |  $27$ |  $22.4$ |  $29$ |  $26.1$ |  $29$ |  $29.2$
Table 2: Iteration counts (it) and condition numbers ($\kappa$); tortoise; $\mu=1,\lambda=1$; including edge averages of normal component of displacement (33) as primal degrees of freedom.  | $\mathpzc{p}=2$ | $\mathpzc{p}=3$ | $\mathpzc{p}=4$ | $\mathpzc{p}=5$ | $\mathpzc{p}=6$
---|---|---|---|---|---
$\ell$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$
$2$ |  $11$ |  $2.8$ |  $13$ |  $4.5$ |  $15$ |  $5.9$ |  $16$ |  $7.8$ |  $16$ |  $9.5$
$3$ |  $14$ |  $4.7$ |  $16$ |  $6.8$ |  $18$ |  $9.7$ |  $18$ |  $12.2$ |  $18$ |  $14.3$
$4$ |  $17$ |  $6.6$ |  $19$ |  $10.8$ |  $20$ |  $14.4$ |  $21$ |  $17.4$ |  $21$ |  $19.9$
$5$ |  $19$ |  $10.2$ |  $21$ |  $15.7$ |  $24$ |  $20.1$ |  $24$ |  $23.6$ |  $24$ |  $26.5$
Table 3: Iteration counts (it) and condition numbers ($\kappa$); tortoise;
$\mu=1,\lambda=1$; including edge averages of both displacement components
(34) as primal degrees of freedom.

For this setting, we obtain a moderate global Korn constant
$\alpha_{\Omega,\Gamma_{D}}^{-2}\approx 211$; here and in what follows, we
have estimated the constant for $\mathpzc{p}=3$, $\ell=2$ using a singular
value decomposition.

The number of iterations needed to reach the threshold and the estimated
condition number are shown in Table 1. As proposed, we only considered the
corner values as primal degrees of freedom. From the tables we see that the
iteration count and condition number increase slightly as $\ell$ and
$\mathpzc{p}$ increases. This is consistent with the convergence estimate in
Theorem 7.

In Table 2, we present the results, where we enrich the primal degrees of
freedom. We follow the approach presented in [23] for the Stokes equations.
Here, besides the corner values, we also use the averages of the normal
components as primal degrees of freedom, i.e., we require

$\int_{\Gamma^{(k,\ell)}}\mathbf{u}^{(k)}\cdot\mathbf{n}^{(k)}\,\mathrm{d}s+\int_{\Gamma^{(k,\ell)}}\mathbf{u}^{(\ell)}\cdot\mathbf{n}^{(\ell)}\,\mathrm{d}s=0$
(33)

for any two patches $\Omega^{(k)}$ and $\Omega^{(\ell)}$ that share a common
edge $\Gamma^{(k,\ell)}:=\partial\Omega^{(k)}\cap\partial\Omega^{(\ell)}$. We
observe that this choice significantly reduces the iteration count and
condition number.

In Table 3, we present the results, where we enrich the primal degrees of
freedom using standard edge averages

$\int_{\Gamma^{(k,\ell)}}\mathbf{u}^{(k)}\,\mathrm{d}s=\int_{\Gamma^{(k,\ell)}}\mathbf{u}^{(\ell)}\,\mathrm{d}s,$
(34)

which we use besides the corner values. This option yields slightly better
results than the results with normal components only, however this approach
leads to an increase of the number of primal degrees of freedom by a factor of
around $\tfrac{3}{2}$, compared to the version presented in Table 2.

###  6.2 Results for a wrench

In this Subsection, we demonstrate that the convergence of the proposed method
does not suffer in cases where the constant from Korn’s inequality is small.
As computational domain, we choose a wrench with a long and relatively thin
handle, see Figure 4. The patch is compost of $64$ patches.

Figure 4: A wrench

Again, we consider the elasticity problem with Lamé parameters $\mu=\lambda=1$
and with $\mu=1,\lambda=100$ and without external volumetric forces
($\mathbf{f}=0$). We choose a homogeneous Dirichlet boundary condition on the
sides where the wrench would touch the bolt, this is marked with blue color in
Figure 4. The remaining boundary has Neumann boundary conditions. For the
parts that are marked with red color, we impose a downward force, represented
as $\mathbf{g}_{N}=(0,-g_{0})$, with $g_{0}=0.00001$. For the remaining
boundary, we have homogeneous Neumann boundary conditions.

The discretization and solution strategy is the same as for the tortoise
(Subsection 6.1). Compared to the tortoise domain, the estimated value of the
constant for Korn’s inequality is much smaller, since we now get
$\alpha_{\Omega,\Gamma_{D}}^{-2}\approx 1500$.

The iteration count and the estimated condition number are shown in Table 4
and 5. For Table 5, ’normal’ dofs where also used as primals. The global Korn
constant is much smaller for the wrench, however, the iteration counts and
condition numbers are better than the tortoise. This indicates, as predicted
from the theory, that the global Korn constant does not effect the condition
number of the IETI-DP solver. The reason the tortoise has larger iteration
counts and condition numbers, is because the individual geometry mappings are
more distorted than for the wrench.

| $\mathpzc{p}=1$ | $\mathpzc{p}=2$ | $\mathpzc{p}=3$ | $\mathpzc{p}=4$ | $\mathpzc{p}=5$ | $\mathpzc{p}=6$
---|---|---|---|---|---|---
$\ell$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  |
$2$ |  $13$ |  $3.8$ |  $17$ |  $6.0$ |  $19$ |  $7.4$ |  $20$ |  $8.4$ |  $20$ |  $9.4$ |  $20$ |  $10.0$
$3$ |  $16$ |  $5.4$ |  $19$ |  $7.9$ |  $21$ |  $9.4$ |  $22$ |  $10.6$ |  $23$ |  $11.5$ |  $23$ |  $12.4$
$4$ |  $18$ |  $7.2$ |  $22$ |  $10.0$ |  $24$ |  $11.7$ |  $25$ |  $13.1$ |  $26$ |  $14.3$ |  $26$ |  $15.3$
$5$ |  $21$ |  $9.2$ |  $24$ |  $12.5$ |  $26$ |  $14.5$ |  $28$ |  $16.1$ |  $29$ |  $17.4$ |  $29$ |  $18.6$
Table 4: Iteration count (it) and condition numbers ($\kappa$); wrench; $\mu=1,\lambda=1$.  | $\mathpzc{p}=1$ | $\mathpzc{p}=2$ | $\mathpzc{p}=3$ | $\mathpzc{p}=4$ | $\mathpzc{p}=5$ | $\mathpzc{p}=6$
---|---|---|---|---|---|---
$\ell$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  |
$2$ |  $8$ |  $2.0$ |  $11$ |  $3.0$ |  $13$ |  $3.8$ |  $14$ |  $4.6$ |  $14$ |  $5.1$ |  $15$ |  $5.7$
$3$ |  $11$ |  $2.8$ |  $14$ |  $4.2$ |  $15$ |  $5.2$ |  $16$ |  $6.1$ |  $17$ |  $6.7$ |  $17$ |  $7.4$
$4$ |  $13$ |  $3.8$ |  $16$ |  $5.6$ |  $18$ |  $6.9$ |  $19$ |  $7.8$ |  $20$ |  $8.6$ |  $20$ |  $9.3$
$5$ |  $15$ |  $5.2$ |  $18$ |  $7.3$ |  $20$ |  $8.8$ |  $22$ |  $9.9$ |  $22$ |  $10.8$ |  $22$ |  $11.5$
Table 5: Iteration count (it) and condition numbers ($\kappa$); wrench; $\mu=1,\lambda=1$; including edge averages of normal component of displacement (33) as primal degrees of freedom.  | $\mathpzc{p}=1$ | $\mathpzc{p}=2$ | $\mathpzc{p}=3$ | $\mathpzc{p}=4$ | $\mathpzc{p}=5$ | $\mathpzc{p}=6$
---|---|---|---|---|---|---
$\ell$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  |
$2$ |  $19$ |  $20.5$ |  $28$ |  $50.6$ |  $34$ |  $67.7$ |  $37$ |  $78.4$ |  $37$ |  $87.2$ |  $37$ |  $95.1$
$3$ |  $22$ |  $33.5$ |  $31$ |  $66.8$ |  $36$ |  $84.0$ |  $39$ |  $95.7$ |  $41$ |  $104.8$ |  $41$ |  $110.6$
$4$ |  $26$ |  $48.7$ |  $34$ |  $83.9$ |  $39$ |  $101.0$ |  $42$ |  $111.8$ |  $43$ |  $119.8$ |  $45$ |  $129.2$
$5$ |  $30$ |  $64.5$ |  $37$ |  $100.8$ |  $42$ |  $119.1$ |  $45$ |  $130.6$ |  $46$ |  $141.8$ |  $47$ |  $147.7$
Table 6: Iteration count (it) and condition numbers ($\kappa$); wrench; $\mu=1,\lambda=100$.  | $\mathpzc{p}=1$ | $\mathpzc{p}=2$ | $\mathpzc{p}=3$ | $\mathpzc{p}=4$ | $\mathpzc{p}=5$ | $\mathpzc{p}=6$
---|---|---|---|---|---|---
$\ell$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  |
$2$ |  $10$ |  $2.9$ |  $16$ |  $6.8$ |  $20$ |  $10.6$ |  $22$ |  $14.0$ |  $23$ |  $16.3$ |  $23$ |  $17.9$
$3$ |  $12$ |  $3.6$ |  $17$ |  $8.1$ |  $21$ |  $12.3$ |  $23$ |  $15.3$ |  $24$ |  $17.6$ |  $24$ |  $19.4$
$4$ |  $13$ |  $4.4$ |  $19$ |  $9.1$ |  $22$ |  $13.4$ |  $24$ |  $16.5$ |  $25$ |  $18.7$ |  $26$ |  $20.4$
$5$ |  $14$ |  $5.0$ |  $20$ |  $10.0$ |  $23$ |  $14.2$ |  $25$ |  $17.3$ |  $26$ |  $19.6$ |  $27$ |  $21.2$
Table 7: Iteration count (it) and condition numbers ($\kappa$); wrench;
$\mu=1,\lambda=100$; including edge averages of normal component of
displacement (33) as primal degrees of freedom.

Next, we demonstrate that the IETI-DP solver also works well if $\lambda/\mu$
is increased, which corresponds to the case that the material is less
compressible. The results are presented in Tables 6 and 7. We observe that the
iteration counts only mildly grow compared to the case $\lambda=\mu=1$. The
condition numbers seem to be much more susceptible to this change. All effects
are reduced if, additionally to the corner values, also the normal components
of the edge averages are chosen as primal degrees of freedom.

In Table 8, we consider the dependence of the performance of the IETI-DP
solver on the length of the wrench. The wrench depicted in Figure 4 has handle
with a length of $10$ (with $x$-coordinate between $2$ and $12$). The other
choices have half, double or quadruple the length. For the first three rows,
we estimated Korn’s constant using a singular value decomposition for
$\mathpzc{p}=3$ and $\ell=2$ ($\ell=1$ for $K=112$). For $K=208$, the (global)
singular value decomposition was not feasible. Again, the homogeneous Neumann
boundary conditions are imposed on the top sides of the 4 patches on the top-
right side of the handle. Although the global Korn constant deteriorates with
the length of the handle, the iteration counts and condition number estimates
show that the solver does not suffer.

Finally, we also aim to illustrate on the discretization error. Namely, in
Figure 5, we present the total bending of the wrench (with $K=64$), i.e., we
consider the value of the $y$-component at the lower right corner of the
wrench. Here, the linear system is solved up to an accuracy of
$\varepsilon=10^{-12}$. In the figures, we present the number of refinement
levels on the $x$-axis, the total bending on the $y$-axes. In all
configurations the total bending seems to converge for $\ell\to\infty$ (which
is to be expected since we have shown a discretization error estimate
$\|\mathbf{u}-\mathbf{u}_{h}\|_{H^{1}(\Omega)}\to 0$ for $h\to 0$). We observe
that the total bending for small $\ell$ and $\mathpzc{p}$ is too small, which
means that the material appears as too stiff. This is what one expects in case
of geometry locking. In all cases, we observe convergence for $\ell\to\infty$.
The convergence is comparably slow for $\mathpzc{p}=1$. Here, the relative
error (compared to the total bending as computed for $\ell=5,\mathpzc{p}=6$)
is around $5\,\%$ for $\lambda=1$ and $62\,\%$ for $\lambda=100$. For
$\mathpzc{p}>1$, already the coarsest grid solution ($\ell=1$) gives rather
good approximations for the total bending; the error is below $1\,\%$ for
$\lambda=1$ and below $3\,\%$ for $\lambda=100$. This reflects the
observations from Section 3.

|  |  | $\mathpzc{p}=1$ | $\mathpzc{p}=2$ | $\mathpzc{p}=3$ | $\mathpzc{p}=4$ | $\mathpzc{p}=5$ | $\mathpzc{p}=6$
---|---|---|---|---|---|---|---|---
length | K | $\alpha_{\Omega,\Gamma_{D}}^{-2}$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  it |  $\kappa$ |  |
$5$ | $40$ | $\approx 402$ |  $16$ |  $5.4$ |  $19$ |  $7.9$ |  $21$ |  $9.4$ |  $22$ |  $10.5$ |  $22$ |  $11.5$ |  $23$ |  $12.4$
$10$ | $64$ | $\approx 1500$ |  $16$ |  $5.4$ |  $19$ |  $7.9$ |  $21$ |  $9.4$ |  $22$ |  $10.6$ |  $23$ |  $11.5$ |  $23$ |  $12.4$
$20$ | $112$ | $\approx 5794$ |  $16$ |  $5.4$ |  $20$ |  $7.9$ |  $21$ |  $9.4$ |  $22$ |  $10.6$ |  $23$ |  $11.5$ |  $23$ |  $12.4$
$40$ | $208$ | $-$ |  $16$ |  $5.4$ |  $19$ |  $7.9$ |  $21$ |  $9.4$ |  $23$ |  $10.6$ |  $23$ |  $11.6$ |  $23$ |  $12.5$
Table 8: Iteration count (it) and condition numbers ($\kappa$) for the wrench;
$\mu=1,\lambda=1,\ell=3$

$1$$2$$3$$4$$5$$-4.3$$-4.2$$-4.1$$-4$$\cdot
10^{-2}$$\ell$$\mathpzc{p}=1$$\mathpzc{p}=2$$\mathpzc{p}=3$$\mathpzc{p}=6$

$1$$2$$3$$4$$5$$-1$$-1.5$$-2$$-2.5$$-3$$\cdot
10^{-2}$$\ell$$\mathpzc{p}=1$$\mathpzc{p}=2$$\mathpzc{p}=3$$\mathpzc{p}=6$

Figure 5: Total bending for the wrench; $\mu=1,\lambda=1$ (left) and
$\lambda=100$ (right).

###  6.3 Concluding remarks

In this paper, a IETI-DP solution strategy was proposed and a bound for the
condition number of the preconditioned system was proven. This bound is
explicit with respect to the grid size, the spline degree, and the local Korn
constants. This means that if the patching of the domain is done well, the
solution strategy will not suffer from geometry locking. The numerical
experiments agree with our theoretical findings. Material locking, which
occurs when $\lambda\gg\mu$, was not addressed in this paper. This will be
addressed in a forthcoming paper.

## Acknowledgments

This work was supported by the Austrian Science Fund (FWF): P31048. The first
author acknowledges support from the Research Council of Norway, grant 301013.
This support is gratefully acknowledged.

## References

  * [1] G. Acosta and R. G. Durán.  Divergence operator and related inequalities.  Springer, 2017.
  * [2] D. Arnold, R. Falk, and R. Winther.  Mixed finite element methods for linear elasticity with weakly imposed symmetry.  Mathematics of Computation, 76(260):1699 – 1723, 2007.
  * [3] F. Auricchio, L. Beirão da Veiga, A. Buffa, C. Lovadina, A. Reali, and G. Sangalli.  A fully “locking-free” isogeometric approach for plane linear elasticity problems: A stream function formulation.  Computer methods in applied mechanics and engineering, 197(1-4):160 – 172, 2007.
  * [4] L. Beirão da Veiga, A. Buffa, G. Sangalli, and R. Vázquez.  Mathematical analysis of variational isogeometric methods.  Acta Numerica, 23:157 – 287, 2014.
  * [5] D. Boffi, F. Brezzi, M. Fortin, et al.  Mixed finite element methods and applications, volume 44.  Springer, 2013.
  * [6] S. Brenner and L. R. Scott.  The Mathematical Theory of Finite Element Methods (Third edition).  Texts in Applied Mathematics. Springer, 2008.
  * [7] J. A. Cottrell, T. J. R. Hughes, and Y. Bazilevs.  Isogeometric Analysis – Toward Integration of CAD and FEA.  John Wiley & Sons, 2009.
  * [8] C. Farhat, M. Lesoinne, P. LeTallec, K. Pierson, and D. Rixen.  FETI-DP: a dual–primal unified FETI method—part I: A faster alternative to the two-level FETI method.  International journal for numerical methods in engineering, 50(7):1523 – 1544, 2001.
  * [9] C. Farhat, M. Lesoinne, and K. Pierson.  A scalable dual-primal domain decomposition method.  Numerical linear algebra with applications, 7(7-8):687 – 714, 2000\.
  * [10] P. L. Gould and Y. Feng.  Introduction to linear elasticity.  Springer, 1994.
  * [11] C. Hofer and U. Langer.  Dual-primal isogeometric tearing and interconnecting solvers for multipatch dG-IgA equations.  Computer Methods in Applied Mechanics and Engineering, 316:2 – 21, 2017.
  * [12] T. J. Hughes, J. A. Cottrell, and Y. Bazilevs.  Isogeometric analysis: CAD, finite elements, NURBS, exact geometry and mesh refinement.  Computer methods in applied mechanics and engineering, 194(39-41):4135 – 4195, 2005.
  * [13] A. Klawonn and O. Rheinbach.  Robust FETI-DP methods for heterogeneous three dimensional elasticity problems.  Computer Methods in Applied Mechanics and Engineering, 196(8):1400 – 1414, 2007.  Domain Decomposition Methods: recent advances and new challenges in engineering.
  * [14] A. Klawonn and O. Widlund.  Dual‐primal FETI methods for linear elasticity.  Communications on Pure and Applied Mathematics, 59:1523 – 1572, 11 2006.
  * [15] S. K. Kleiss, C. Pechstein, B. Jüttler, and S. Tomar.  IETI–isogeometric tearing and interconnecting.  Computer Methods in Applied Mechanics and Engineering, 247:201 – 215, 2012.
  * [16] J. Mandel, C. R. Dohrmann, and R. Tezaur.  An algebraic theory for primal and dual substructuring methods by constraints.  Appl. Numer. Math., 54(2):167 – 193, 2005.
  * [17] M. Montardini, G. Sangalli, R. Schneckenleitner, S. Takacs, and M. Tani.  A IETI-DP method for discontinuous Galerkin discretizations in Isogeometric Analysis with inexact local solvers.  Technical report, 2022.  (Available as arXiv preprint arXiv:2206.08416).
  * [18] L. F. Pavarino and S. Scacchi.  Isogeometric block FETI-DP preconditioners for the Stokes and mixed linear elasticity systems.  Computer Methods in Applied Mechanics and Engineering, 310:694 – 710, 2016.
  * [19] E. Sande, C. Manni, and H. Speleers.  Sharp error estimates for spline approximation: explicit constants, $n$-widths, and eigenfunction convergence.  Math. Models Methods Appl. Sci., 29(6):1175 –– 1205, 2019.
  * [20] R. Schneckenleitner and S. Takacs.  Condition number bounds for IETI-DP methods that are explicit in $h$ and $p$.  Mathematical Models and Methods in Applied Sciences, 30(11):2067 – 2103, 2020.
  * [21] R. Schneckenleitner and S. Takacs.  IETI-DP methods for discontinuous Galerkin multi-patch Isogeometric Analysis with T-junctions.  Computer Methods in Applied Mechanics and Engineering, 393:114694, 2022.
  * [22] R. Schneckenleitner and S. Takacs.  Inexact IETI-DP for conforming isogeometric multi-patch discretizations.  In Spectral and High Order Methods for Partial Differential Equations ICOSAHOM 2021, 2023.  To appear. (Also available as arXiv preprint: arXiv:2110.06087).
  * [23] J. Sogn and S. Takacs.  Stable discretizations and IETI-DP solvers for the Stokes system in multi-patch Isogeometric Analysis.  ESAIM: Mathematical Modelling and Numerical Analysis, 57(2):921 –– 925, 2022.
  * [24] R. Stenberg.  On the construction of optimal mixed finite element methods for the linear elasticity problem.  Numerische Mathematik, 48(4):447 – 462, 1986.
  * [25] S. Takacs.  Robust approximation error estimates and multigrid solvers for isogeometric multi-patch discretizations.  Mathematical Models and Methods in Applied Sciences, 28(10):1899 – 1928, 2018.
  * [26] S. Takacs and T. Takacs.  Approximation error estimates and inverse inequalities for B-splines of maximum smoothness.  Math. Models Methods Appl. Sci., 26(7):1411 – 1445, 2016.
  * [27] O. Widlund, S. Zampini, S. Scacchi, and L. F. Pavarino.  Block FETI–DP/BDDC preconditioners for mixed isogeometric discretizations of three-dimensional almost incompressible elasticity.  Mathematics of Computation, 90(330):1773 – 1797, 2021.