Nicole Spillane
email: spillane@ann.jussieu.fr

* -Daniel J Rixen

) = Ker

(p ⊤ * ,N F P * ,N )g

: = I -G

P * ,N P

P ⊤ * ,N And P * ,N P

P ⊤ * ,N + P * ,N G

N = P ⊤ * ,N P ⊤

P ⊤ * ,N (iii) P ⊤

P ⊤ * ,N = P ⊤ * ,N (i -F P

P ⊤ * ,N = P ⊤ * ,N P ⊤

P ⊤ * ,N = P ⊤

P ⊤ * ,N P ⊤

N P

Then

P ⊤ * ,N Μ = P * ,N P

P ⊤ * ,N ) -1 P

Automatic spectral coarse spaces for robust FETI and BDD algorithms

Keywords: Domain decomposition, FETI, BDD, robustness, scalability, varying coefficients, irregular partitions, interface heterogeneity

We introduce spectral coarse spaces for the BDD (Balanced Domain Decomposition) and FETI (Finite Element Tearing and Interconnecting) methods. These coarse spaces are specifically designed for the twolevel methods to be scalable and robust with respect to the coefficients in the equation and the choice of the decomposition. We achieve this by solving generalized eigenvalue problems on the interfaces between subdomains to identify the modes which slow down convergence. Theoretical bounds for the condition numbers of the preconditioned operators which depend only on a chosen threshold and the maximal number of neighbours of a subdomain are presented and proved. For FETI there are two versions of the two-level method: one based on the full Dirichlet preconditioner and the other on the, cheaper, lumped preconditioner. Some numerical tests confirm these results.

INTRODUCTION

In domain decomposition it is a real challenge to solve problems with a decomposition given by an automatic partitioner [START_REF] Karypis | METIS: A software package for partitioning unstructured graphs, partitioning meshes, and computing fill-reducing orderings of sparse matrices[END_REF][START_REF]SCOTCH: Static mapping, graph partitioning, and sparse matrix block ordering package[END_REF] which does not take into account all the difficulties in the problem for the simple reason that there are too many. One well known challenge for elliptic problems is when the coefficients in the equation are highly heterogeneous. This is often the case in practical applications. Classical coarse spaces are known to give good results when the jumps in the coefficients are across subdomain interfaces (see e.g. [START_REF] Dryja | Multilevel Schwarz methods for elliptic problems with discontinuous coefficients in three dimensions[END_REF][START_REF] Mandel | Balancing domain decomposition for problems with large jumps in coefficients[END_REF][START_REF] Dohrmann | An overlapping Schwarz algorithm for almost incompressible elasticity[END_REF][START_REF] Dohrmann | Hybrid domain decomposition algorithms for compressible and almost incompressible elasticity[END_REF]) or inside the subdomains and not near their boundaries (cf. [START_REF] Pechstein | Scaling up through domain decomposition[END_REF][START_REF] Pechstein | Analysis of FETI methods for multiscale PDEs[END_REF]). However, when the discontinuities are along subdomain interfaces, classical results break down, and one observes very bad convergence of the iterative solvers for the interface problem (see e.g. [START_REF] Bhardwaj | Application of the FETI method to ASCI problems: Scalability results on a thousand-processor and discussion of highly heterogeneous problems[END_REF][START_REF] Klawonn | Robust FETI-DP methods for heterogeneous three dimensional elasticity problems[END_REF]). It is also well known that non-smooth decompositions (where the interfaces are jagged) [START_REF] Klawonn | An analysis of a feti-dp algorithm on irregular subdomains in the plane[END_REF] or bad aspect ratios of the domains [START_REF] Farhat | Mesh partitioning for implicit computations via iterative domain decomposition: impact and optimization of the subdomain aspect ratio[END_REF] can also lead to poor convergence.This is what we work to improve: we aim to design a method for which the convergence rate does not depend on the choice of the decomposition into subdomains or on any of the coefficients in the equations.

In order to achieve this we will use the strategy introduced in the additive Schwarz framework by [START_REF] Spillane | A robust two level domain decomposition preconditioner for systems of PDEs[END_REF][START_REF] Spillane | Abstract robust coarse spaces for systems of PDEs via generalized eigenproblems in the overlaps[END_REF] and [START_REF] Efendiev | Robust domain decomposition preconditioners for abstract symmetric positive definite bilinear forms[END_REF]. This strategy is based on the abstract theory of the two-level additive Schwarz method [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF]. The strategy is to write the Schwarz theory up to the point where it depends on the set of equations we are dealing with and where assumptions on the coefficient distribution with respect to the decomposition into subdomains are needed to write estimates which do not depend on the parameters. For the Darcy equation (-∇ • ∇(αu) = b) with the minimal coarse space (the constant functions) the Poincaré inequality and trace theorem are needed to complete the proof and they require quite strong assumptions. Instead, the authors in [START_REF] Efendiev | Robust domain decomposition preconditioners for abstract symmetric positive definite bilinear forms[END_REF][START_REF] Spillane | Abstract robust coarse spaces for systems of PDEs via generalized eigenproblems in the overlaps[END_REF][START_REF] Spillane | A robust two level domain decomposition preconditioner for systems of PDEs[END_REF] propose to solve a generalized eigenvalue problem in each subdomain which selects what modes of the solution satisfy the required estimates for a chosen constant. The other modes, which do not satisfy the estimate, are used to build the coarse space and are basically taken care of with a direct solve in the coarse space. This is what we will refer to as the Schwarz-GenEO coarse space (Generalized Eigenvalues in the Overlaps).

It leads to a two-level method with a convergence rate chosen a priori for problems described by symmetric positive definite matrix. The idea to use eigenvalue problems to build a coarse space is not new, it was first explored in the algebraic multigrid community. In [START_REF] Brezina | An iterative method with convergence rate chosen a priori[END_REF], a strategy to build a coarse space based on spectral information is presented that allows to achieve any a priori chosen target convergence rate. This idea was further developed and implemented in the spectral AMGe method in [START_REF] Chartier | Spectral AMGe (ρAMGe)[END_REF]. More recently, in the framework of two-level overlapping Schwarz, [START_REF] Galvis | Domain decomposition preconditioners for multiscale flows in highcontrast media[END_REF][START_REF] Galvis | Domain decomposition preconditioners for multiscale flows in high contrast media: Reduced dimension coarse spaces[END_REF][START_REF] Nataf | A two level domain decomposition preconditioner based on local Dirichlet-to-Neumann maps[END_REF][START_REF] Dolean | A coarse space construction based on local Dirichlet to Neumann maps[END_REF][START_REF] Efendiev | Robust domain decomposition preconditioners for abstract symmetric positive definite bilinear forms[END_REF][START_REF] Spillane | A robust two level domain decomposition preconditioner for systems of PDEs[END_REF][START_REF] Spillane | Abstract robust coarse spaces for systems of PDEs via generalized eigenproblems in the overlaps[END_REF]] also build coarse spaces for problems with highly heterogeneous coefficients by solving local eigenproblems. However, compared to the earlier works in the AMG context all of these approaches focus on generalized eigenvalue problems. We can distinguish three sets of methods that differ by the choice of the bilinear form on one side of the generalized eigenproblem. First, in the work of [START_REF] Galvis | Domain decomposition preconditioners for multiscale flows in highcontrast media[END_REF][START_REF] Galvis | Domain decomposition preconditioners for multiscale flows in high contrast media: Reduced dimension coarse spaces[END_REF] for the Darcy equation it is the local mass matrix, or a 'homogenised' version obtained by using a multiscale partition of unity. In [START_REF] Nataf | A two level domain decomposition preconditioner based on local Dirichlet-to-Neumann maps[END_REF][START_REF] Dolean | A coarse space construction based on local Dirichlet to Neumann maps[END_REF] it corresponds to an L 2 -product on the subdomain boundary, so that the problem can be reduced to a generalized eigenproblem for the Dirichlet to Neumann operator. This method was analysed in [START_REF] Dolean | Analysis of a two-level Schwarz method with coarse spaces based on local Dirichlet-to-Neumann maps[END_REF]. The latest set of papers, of which this one is inspired, [START_REF] Efendiev | Robust domain decomposition preconditioners for abstract symmetric positive definite bilinear forms[END_REF][START_REF] Spillane | A robust two level domain decomposition preconditioner for systems of PDEs[END_REF][START_REF] Spillane | Abstract robust coarse spaces for systems of PDEs via generalized eigenproblems in the overlaps[END_REF], uses yet another type of bilinear form inspired by the theory. There have also been some recent multilevel extensions of some of the above approaches [START_REF] Efendiev | Spectral element agglomerate algebraic multigrid methods for elliptic problems with high contrast coefficients[END_REF][START_REF] Efendiev | Multiscale spectral AMGe solvers for high-contrast flow problems[END_REF][START_REF] Willems | Robust multilevel methods for general symmetric positive definite operators[END_REF]. The approach in [START_REF] Napov | Algebraic analysis of aggregation-based multigrid[END_REF][START_REF] Napov | An algebraic multigrid method with guaranteed convergence rate[END_REF], in the multigrid framework is also comparable.

The purpose of this paper is to extend the GenEO strategy [START_REF] Efendiev | Robust domain decomposition preconditioners for abstract symmetric positive definite bilinear forms[END_REF][START_REF] Spillane | A robust two level domain decomposition preconditioner for systems of PDEs[END_REF][START_REF] Spillane | Abstract robust coarse spaces for systems of PDEs via generalized eigenproblems in the overlaps[END_REF] to the BDD (Balancing Domain Decomposition) algorithm and the FETI (Finite Element Tearing and Interconnecting) algorithm. These are two well known non overlapping domain decomposition methods. Up until now the GenEO strategy has been applied in the context of overlapping Schwarz which was first introduced in [START_REF] Matsokin | A Schwarz alternating method in a subspace[END_REF]. The idea of a coarse space correction goes back to [START_REF] Nicolaides | Deflation of conjugate gradients with applications to boundary value problems[END_REF][START_REF] Bramble | The construction of preconditioners for elliptic problems by substructuring[END_REF] and the two-level overlapping Schwarz preconditioner is due to [START_REF] Dryja | Some domain decomposition algorithms for elliptic problems. Iterative methods for large linear systems[END_REF]. As for the Balancing Domain Decomposition (BDD) method, it is the work of [START_REF] Mandel | Balancing domain decomposition[END_REF] who added a coarse space to the preexisting Neumann Neumann method [START_REF] De Roeck | Analysis and test of a local domain decomposition preconditioner[END_REF] to deal with singularities in the local problems. We will refer to the analysis of BDD in [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF] which is very closely related to the analysis of the two-level Schwarz preconditioner. Finally, the FETI algorithm was first introduced in [START_REF] Farhat | A method of finite element tearing and interconnecting and its parallel solution algorithm[END_REF] and the convergence proof is due to [START_REF] Mandel | Convergence of a substructuring method with lagrange multipliers[END_REF][START_REF] Tezaur | Analysis of lagrange multiplier based domain decomposition[END_REF]. It is generalized in [START_REF] Klawonn | Feti and neumann-neumann iterative substructuring methods: connections and new results[END_REF]. Coarse spaces for the FETI method are introduced first in [START_REF] Farhat | A scalable lagrange multiplier based domain decomposition method for time-dependent problems[END_REF] and further developed in [START_REF] Farhat | The two-level FETI method for static and dynamic plate problems -part I: An optimal iterative solver for biharmonic systems[END_REF][START_REF] Farhat | The two-level FETI method -part II: Extension to shell problems. parallel implementation and performance results[END_REF]. In [START_REF] Toselli | A feti domain decomposition method for edge element approximations in two dimensions with discontinuous coefficients[END_REF] a two level FETI method is also introduced for a particular problem and a convergence result is proved. However we will follow a very different approach here both for choosing the coarse space and also for writing the proof. In both cases (BDD and FETI) the generalized eigenvalue problem which we solve is used to prove a bound for the largest eigenvalue of the preconditioned operator. As usual the lower bound for the eigenvalues of the preconditioned operator is 1 regardless of the coarse space.

The rest of the article is organized as follows. In Section 1 we introduce the notation which will be needed for both algorithms. In Section 2 we introduce the two-level GenEO preconditioner for the BDD algorithm. And in Section 3 we introduce the two-level preconditioner for the FETI algorithm. The definitions of each of the coarse spaces with the corresponding generalized eigenvalue problems can be found in Definitions 2.3 and 3.7 respectively. These generalized eigenvalue problems are chosen specifically to ensure that the so called stable splitting properties in Lemmas 2.8 and 3.12 are satisfied. As for the convergence results they are stated (and proved) in Theorems 2.11 and 3.14. Finally in section 4 we give a few numerical results.

NOTATION FOR FETI AND BDD

For a given domain Ω ∈ R d and a finite dimensional Hilbert space W h (Ω), given a symmetric, positive definite bilinear form,

â(•, •) : W h (Ω) × W h (Ω) → R, (1.1)
and an element ĝ ∈ W h (Ω) ′ , we consider the problem of finding u ∈ W h (Ω), such that

â(u, v) = ĝ(v), ∀ v ∈ W h (Ω). (1.2)
In order to introduce the BDD and FETI algorithms we will need to introduce notation for discrete operators at the global and local (on each subdomain) levels.

Problem setting

We begin by rewriting Problem (1.2) in an algebraic framework. As usual in the finite element setting, we start with a triangulation T h of Ω: Ω = τ ∈T h τ and a basis {φ k } 1≤k≤N for the finite element space W h (Ω).

Assumption 1.1

Given any element τ of the mesh T h , let W h (τ ) := {u| τ : u ∈ W h (Ω)}. We assume that for each element τ ∈ T h , there exists a symmetric positive semi-definite (spsd) bilinear form a τ :

W h (τ ) × W h (τ ) → R, such that â(u, v) = τ ∈T h a τ (u| τ , v| τ ), ∀ u, v ∈ W h (Ω),
and an element

g τ ∈ W h (τ ) ′ such that ĝ(v) = τ ∈T h g τ (v| τ ), ∀ v ∈ W h (Ω).
The stiffness matrix is assembled with the following entries

( K) kl := â(φ k , φ l ) = τ ∈T h a τ (φ k|τ , φ l|τ ) , ∀ k, l = 1, . . . , n, (1.3)
and the discrete right hand side f ∈ R n is defined by the entries

( f ) k := ĝ(φ k ) = τ ∈T h g τ (φ k | τ ) , ∀ k = 1, . . . , n.
As is quite customary we identify vectors of degrees of freedom, which are in some spaces R m , with the associated finite element functions. Operators between the spaces are represented as matrices, and we frequently commit an abuse of notation by using matrices and operators interchangeably. With this abuse of notation the original problem (1.2) is equivalent to the linear system:

find u ∈ W h (Ω) such that Ku = f , (1.4)
with K symmetric, positive definite (spd).

Local setting and notation

Local Setting We introduce a partition of the global domain Ω into N non-overlapping subdomains Ω i which are resolved by the mesh

Ω = N i=1
Ωi and

Ω i ∩ Ω j = ∅, i = j,
and the resulting set of boundaries between subdomains

Γ := i =i ′ Ωi ∩ Ωi ′ .
The reason why we have required the information on the non-assembled stiffness matrices is that we want to have access to local matrices for any choice of the partition into subdomains. In order to do this we also need to define local finite element spaces and local bilinear forms.

Assumption 1.2

The basis functions φ k are continuous on Ω. In particular for any subset D ⊂ Ω the restriction φ k|D of φ k to D is well defined.

where 0 |D : D → R is identically zero. Then the finite element space on D is defined as

W h (D) := {u |D ; u ∈ W h (Ω)} = span{φ k|D ; k ∈ dof (D)}. (1.6)
The second equality in the definition of W h (D) is an immediate consequence.

Definition 1.4 (Local bilinear forms and local right hand sides)

For any open subset D ⊂ Ω which is resolved by the mesh T h , let the local bilinear form on D be

a D : W h (D) × W h (D) → R; a D (v, w) := τ ⊂D a τ (v |τ , w |τ ), (1.7)
and the local right hand side be the element

g D ∈ W ′ h (D); g D (v) := τ ⊂D g τ (v |τ ). (1.8)
For any i = 1, . . . , N , the space of finite element functions on each Ω i follows from (1.6) with

D = Ω i : W h (Ω i ) = {u |Ωi ; u ∈ W h (Ω)},
as well as the trace spaces for D = ∂Ω i ∩ Γ:

W i := W h (Γ ∩ ∂Ω i ) = {u |Γ∩∂Ωi ; u ∈ W h (Ω)}.
Finally, we define the product space

W := N i=1 W i .
We know from (1.6) that W i = span{φ k|∂Ω i∩Γ ; k ∈ dof (∂Ω i ∩ Γ)}, we make the further assumption that this set of functions is a basis of W i .

Assumption 1.5

The set

{φ k|∂Ω i∩Γ ; k ∈ dof (∂Ω i ∩ Γ)} is a basis of W i .
Throughout the analysis, we will consider elements in the product space W . Each component u i ∈ W i is defined on a part Γ ∩ ∂Ω i of the boundary and two contributions from two neighbouring subdomains do not necessarily match on the shared interface. This is a result of the partition of Ω into subdomains. Our finite element approximation of the elliptic problem is, however, based on functions in W h (Ω) which are defined on the whole of Ω with one value per degree of freedom. We denote the space of restrictions of these functions to the set of internal boundaries Γ by Ŵ :

Ŵ := W h (Γ) = {u |Γ ; u ∈ W h (Ω)} = span{φ k|Γ ; k ∈ dof (Γ)} .
(1.9)

Next we introduce interpolation (prolongation) operators R ⊤ i : W i → Ŵ for i = 1, . . . , N :

∀ u i = k∈dof (Γ∩∂Ωi) α k i φ k|Γ∩∂Ω i (α k i ∈ R); R ⊤ i u i := k∈dof (Γ∩∂Ωi) α k i φ k|Γ .
These are the natural interpolation operators represented by boolean matrices: the continuous global function R ⊤ i u i ∈ Ŵ shares the same values as u i for degrees of freedom in dof (Γ ∩ ∂Ω i ) and has no contributions from any other degrees of freedom. The corresponding restriction operator R i : Ŵ → W i is defined as

∀ u = k∈dof (Γ) α k φ k|Γ (α k ∈ R); R i u := k∈dof (Γ∩∂Ωi) α k φ k|Γ∩∂Ω i . We note that Ŵ ⊂ W and Ŵ = N i=1 R ⊤ i W i .
It is obvious from the definition of R ⊤ i and Assumption 1.5 that for i = 1, . . . , N and u i ∈ W i :

u i = 0 |Γ∩∂Ωi ⇔ R ⊤ i u i = 0 |Γ .
(1.10)

Stiffness matrices The local stiffness matrix K

i : W h (Ω i ) → W h (Ω i )
is the matrix associated with bilinear form a Ωi defined by (1.7) for D = Ω i . From these, the stiffness matrix on the product space is defined as

K : W h (Ω 1 ) × . . . W h (Ω N ) → W h (Ω 1 ) × . . . W h (Ω N ); K :=    K 1 0 . . . 0 0 K 2 . . . 0 . . . . . . . . . . . . 0 0 . . . K N    (1.11) so that Ku = (K 1 u 1 , . . . , K N u N ) ⊤ , ∀ u = (u 1 , . . . , u N ) ⊤ ∈ W h (Ω 1 ) × . . . W h (Ω N ).
(1.12)

Schur complement matrices

The degrees of freedom dof (Ω i ) in W h (Ω i ) can be split into the set b i := dof (Γ ∩ ∂Ω i ) of degrees of freedom that are also in the trace space W i and the remainder

I i := dof (Ω i ) \ dof (Γ ∩ ∂Ω i ).
This way we can rewrite the local stiffness matrix in block formulation

K i = K bibi i K biIi i K Iibi i K IiIi i .
The interior variables of any subdomain are then eliminated in work that can be parallelized across the subdomains. The resulting matrices are the local Schur complements

S i : W i → W i ; S i := K bibi i -K biIi i (K IiIi i ) -1 K Iibi i , i = 1, . . . , N, (1.13)
and the Schur complement on the product space is

S : W 1 × . . . W N W → W 1 × . . . W N W ; S :=    S 1 0 . . . 0 0 S 2 . . . 0 . . . . . . . . . . . . 0 0 . . . S N    (1.14)
so that Su = (S 1 u 1 , . . . , S N u N ) ⊤ , ∀ u = (u 1 , . . . , u N ) ⊤ ∈ W.

(1.15)

The Schur complement S on the product space W admits the following counterpart Ŝ for functions in Ŵ :

Ŝ : Ŵ → Ŵ ; Ŝu := N i=1 R ⊤ i S i R i u. (1.16)
We notice that this is the usual Schur complement for the global problem reduced to the set Γ of internal boundaries:

Ŝ = Kbb -KbI ( KII ) -1 KIb , (1.17)
where Kbb , KbI , KII and KIb are the components in the bloc formulation of K K = Kbb KbI KIb KII , b := dof (Γ) and I := dof (Ω) \ dof (Γ).

(1.18)

Lumped matrices

In the FETI literature the lumped version of the stiffness matrix is the extraction of the entries in the stiffness matrix which correspond to boundary degrees of freedom.

We have already introduced Kbb and K bibi i , let K bb be the counterpart on the product space W :

K bb : W 1 × . . . W N W → W 1 × . . . W N W ; K bb :=     K b1b1 1 0 . . . 0 0 K b2b2 2 . . . 0 . . . . . . . . . . . . 0 0 . . . K bN bN N     . (1.19)
We notice that Kbb = N i=1 R ⊤ i K bibi i R i and the next Lemma gives an important relation between lumped matrices and Schur complement matrices.

Lemma 1.6

For any û ∈ Ŵ and any u ∈ W the following inequalities hold Ŝ û, û ≤ Kbb û, û and Su, u ≤ K bb u, u .

Proof

Let û ∈ Ŵ . Then by definition of Ŝ Ŝ û, û = ( Kbb -KbI ( KII ) -1 KIb )û, û = Kbb û, û -( KII ) -1 KIb û, KIb û .
The first inequality follows by noticing that ( KII ) -1 KIb û, KIb û ≥ 0 because ( KII ) -1 is spd. For the second, let u ∈ W . Then by definition of S Su, u

= N i=1 S i u i , u i = N i=1 (K i bibi -K i biIi (K i IiIi ) -1 K i Iibi )u i , u i = K bb u, u - N i=1 (K i IiIi ) -1 K i Iibi u i , K i Iibi u i .
And the second inequality follows by noticing that

(K i IiIi ) -1 K i Iibi u i , K i Iibi u i ≥ 0 for any i = 1, . . . , N because (K i IiIi ) -1 is spd.
Right hand sides In order to reduce the problem to the set of interfaces between subdomains, we define the following right hand side

fΓ := f b -KbI ( KII ) -1 f I , (1.20)
which is the right hand side of the original problem (1.4) condensed onto the degrees of freedom in Ŵ . As for the right hand side on the product space W , for each subdomain i = 1, . . . , N : first let f i be the local right hand side given by (1.8) with D = Ω i . Then condense it onto the interfaces following:

f Γ,i := f bi i -K biIi i (K IiIi i ) -1 f Ii i . (
We have used the identification between the finite element representation of f i and its vector representation.) Finally, the right hand side for the problem condensed onto the space W is

f Γ =   f Γ,1 . . . f Γ,N   . (1.21)
Most of this notation is summed up in Table I at the beginning of the article. Some comments are given in subsection 1.4, along with an important lemma on which of these matrices are positive definite.

Remark 1.7 Assumption 1.1 is actually stronger than what we really need but enables the use of any partition into subdomains and allowed us to define each component of the algorithm thoroughly. For a given non overlapping partition into subdomains it is enough to have access to the local matrices K i on each subdomain, the local right hand sides f i , the local-global interpolation operators R ⊤ i and the information on the boundary of each subdomain Γ ∩ ∂Ω i .

Partition of unity and weighted operators

An important role in the description of the BDD algorithms is played by a weighting (counting) function on W . As in the original GenEO algorithm [START_REF] Spillane | A robust two level domain decomposition preconditioner for systems of PDEs[END_REF][START_REF] Spillane | Abstract robust coarse spaces for systems of PDEs via generalized eigenproblems in the overlaps[END_REF] this induces partition of unity operators Ξ i which act directly on the degrees of freedom of the finite element functions. Then for any function u i ∈ W i written as

u i = k∈dof (Γ∩∂Ωi) α k i φ k|Γ∩∂Ω i , α k i ∈ R,
the local partition of unity operator Ξ i : W i → W i is defined by:

Ξ i (u i ) := k µ k i α k i φ k|Γ∩∂Ω i ,
where

µ k i is the k-th entry in µ i . The inverse Ξ -1 i : W i → W i is defined by: Ξ -1 i (u i ) := k 1 µ k i α k i φ k|Γ∩∂Ω i .
It is clear that the Ξ i define a partition of unity from Ŵ onto the product space

W = W 1 × • • • × W N in the sense that u = N i=1 R ⊤ i Ξ i (R i u) ∈Wi , ∀ u ∈ Ŵ .
It is also clear that Ξ -1 i is the inverse of Ξ i since any

u i ∈ W i satisfies Ξ -1 i (Ξ i (u i )) = Ξ i (Ξ -1 i (u i )) = u i .
Remark 1.9 Two common choices for µ are the multiplicity scaling where µ k i is chosen as

(#{i = 1, . . . , N ; k ∈ dof (Γ ∩ ∂Ω i )})
-1 and the K-scaling where µ depends on the diagonal entries of the stiffness matrices [START_REF] Rixen | A simple and efficient extension of a class of substructure based preconditioners to heterogeneous structural mechanics problems[END_REF][START_REF] Klawonn | Feti and neumann-neumann iterative substructuring methods: connections and new results[END_REF]. In the numerical result section we mostly use K-scaling.

We introduce the local bilinear forms which correspond to the local Schur complements S i as follows. For i = 1, . . . , N define

s i : W i × W i → R, s i (u i , v i ) := S i u i , v i ; ∀u i , v i ∈ W i . (1.22)
Next we use the partition of unity operators to define weighted versions of the Schur complements which will be instrumental in defining the BDD algorithm.

Definition 1.10 (Weighted Schur complements)

For any i = 1, . . . , N , let si :

W i × W i → R be the bilinear form defined by si (u i , v i ) := s i (Ξ -1 i (u i ), Ξ -1 i (v i )); ∀ u i , v i ∈ W i , (1.23)
where s i is the local Schur complement, and Ξ -1 i is the inverse partition of unity operator introduced in Definition 1.8. Next, let the matrix Si : W i → W i be the matrix counterpart of si :

Summary of the notation and complements

We have introduced quite a lot of notation. Table I at the beginning of the article sums up most of the notation which will appear in the description of the algorithms and the reference to where it is first introduced. Some of the operators are introduced for the first time (â bb , a bb , ŝ and s) as the bilinear forms associated with a matrix. More precisely, let âbb and ŝ be defined as

âbb : Ŵ × Ŵ → R; âbb (û, v) := Kbb û, v and ŝ : Ŵ × Ŵ → R; ŝ(û, v) := Ŝ û, v ,
for any û and v ∈ Ŵ , and let a bb and s be defined as

a bb : W × W → R; a bb (u, v) := K bb u, v and s : W × W → R; s(u, v) := Su, v ,
for any u and v ∈ W .

The operators with a • always correspond to functions defined either on the whole of Ω or the whole of Γ. The subscript i always refers to a local operator defined on a subdomain Ω i or its boundary. Operators without a • or a subscript i are defined on the product spaces. Finally operators Si are weighted by the inverse partition of unity operators.

In many cases the local stiffness matrices K i are not spd on all floating subdomains. (A floating subdomain is a subdomain which does not touch the Dirichlet part of the boundary). For example, in the case of the Darcy equation, the kernel of K i for a floating subdomain is the set of constant functions. In the case of linear elasticity, the kernel of K i is the set of rigid body motions. It is easy to see that these kernels induce kernels for the corresponding Schur complements S i as well as their weighted counterparts Si and, possibly, the lumped matrices K bibi i .The next lemma makes precise which matrices are positive definite. They are all symmetric positive semi definite.

Lemma 1.11

The stiffness matrix K, lumped stiffness matrix K bb and Schur complement S, which correspond to the product spaces, can be singular. Their respective counterparts, K, Kbb and Ŝ, on the original spaces of functions W h (Ω) and Ŵ are symmetric positive definite. Finally, under Assumption 1.5 each of the local matrices R i Kbb R ⊤ i and R i ŜR ⊤ i is also symmetric positive definite. Proof The fact that K and Ŝ are positive definite is clear because the original problem is well posed. The positive definiteness of Kbb follows from Lemma 1.6 and the positive definiteness of Ŝ:

let u ∈ Ŵ Kbb u, u = 0 ⇒ Ŝu, u = 0 ⇒ u = 0.
The positive definiteness of R i ŜR ⊤ i and R i Kbb R ⊤ i is obvious from the positive definiteness of K and Ŝ and (1.10) which is a direct consequence of Assumption 1.5.

Remark 1.12

Note that in nearly all practical cases K bb is also symmetric positive definite.

We are now ready to introduce the BDD preconditioner.

BALANCING DOMAIN DECOMPOSITION

The problem which we solve is the original problem (1.4) reduced to the set Γ of interfaces between subdomains: find u ∈ Ŵ such that Ŝu = fΓ .

(2.1)

One level BDD preconditioner in the abstract Schwarz framework [16]

The only thing that is needed in order to define the one-level preconditioner is a solver on each subdomain. Then we will precondition the global problem (2.1) with a sum of these local solves.

The usual BDD strategy is to use the weighted Schur complements Si introduced in Definition 1.10 to build local problems. Then each local solve is the solution of a Neumann problem: S † i .

Definition 2.1 (One level preconditioner)

For each i = 1, . . . , N , let Pi and P i be defined as

Pi := S † i R i Ŝ and P i := R ⊤ i Pi , (2.2)
where S † i is a pseudo inverse of Si . Equivalently for any u ∈ Ŵ , Pi u is the unique vector in range( S † i ) which satisfies

si ( Pi u, v i ) = ŝ(u, R ⊤ i v i ), ∀ v i ∈ W i . (2.
3)

The one-level preconditioner is the sum of local solves

N i=1 R ⊤ i S † i R i so the one-level preconditioned operator is N i=1 P i .
The next lemma gives a lower bound on the eigenvalues of the one-level preconditioned operator. It does not depend on the specific choice of the pseudo inverse or on any coarse space.

Essentially what we do is check that a stable splitting assumption (Assumption 2.2 in [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF]) holds on the whole of Ŵ . Then we give the result of Lemma 2.5 in [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF] which is that this implies a lower bound for the condition number of the one-level preconditioned operator. One of the assumptions in [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF] is that the local bilinear forms ( Si in this case) be positive definite. Here they are only positive semi definite but the proof goes through in the exact same way so we don't give it again. Lemma 2.2 (Stable splitting -Lower bound for the eigenvalues of the preconditioned operator) For any u ∈ Ŵ there exists a stable splitting

(v 1 , . . . , v N ) of u onto W = W 1 × • • • × W N : u = R ⊤ 1 v 1 + • • • + R ⊤ N v N ; v i ∈ W i and N i=1 si (v i , v i ) ≤ ŝ(u, u). (2.4)
This implies that the one-level preconditioned operator satisfies

ŝ(u, u) ≤ ŝ N i=1 P i u, u for any u ∈ Ŵ .
(2.5)

Proof

Let u ∈ Ŵ . The fact that, by definition, the operators Ξ i define a partition of unity allows us to write an obvious splitting of u onto W :

(v i := Ξ i (R i u), ∀ i = 1, . . . , N ) ⇒ u = N i=1 R ⊤ i v i .
We prove (2.4) for this splitting using only the definitions of si and ŝ:

N i=1 si (v i , v i ) = N i=1 s i (Ξ -1 i (Ξ i (R i u)), Ξ -1 i (Ξ i (R i u)) = N i=1 s i (R i u, R i u) = ŝ(u, u).
The second part of the lemma is the result of Lemma 2.5 in [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF], we refer the reader to there for the proof.

The fact that (2.5) provides a lower bound for the eigenvalues of the preconditioned operator N i=1 P i is easy to see: suppose u is an eigenvector associated with eigenvalue λ, then

N i=1 P i u = λu ⇒ Ŝ N i=1 P i u = λ Ŝu ⇒ ŝ( N i=1 P i u, u) = λŝ(u, u),
and (2.5) implies that λ ≥ 1.

In other words the lower bound for the eigenvalues of the preconditioned operator does not depend on the choice of the coarse space. This is a big difference with the Additive Schwarz method where the proof of a lower bound depends very strongly on the choice of the coarse space and on restrictive assumptions on the coefficient distribution. This is why the Schwarz-GenEO strategy in [START_REF] Spillane | Abstract robust coarse spaces for systems of PDEs via generalized eigenproblems in the overlaps[END_REF] is precisely to build an enriched coarse space for which the stable splitting property and thus a lower bound for the spectrum of the preconditioned operator hold regardless of the partition into subdomains and the coefficient distribution. Luckily, the upper bound for the eigenvalues of the Additive Schwarz operator depends only on the number of neighbours of each subdomain enabling the proof of a bound for the condition number of the preconditioned operator.

Here the situation is reversed: Lemma 2.2 gives a lower bound for the eigenvalues of the preconditioned operator which does not depend on the choice of the coarse space thanks to the adequate weighting of the local solvers. However the upper bound requires more work and with the usual coarse space it can only be independent of the coefficients in the equation if some assumptions on the coefficient distribution are satisfied. The GenEO strategy will enable us to waive all of these assumptions.

GenEO coarse space for BDD

The abstract Schwarz theory tells us that the upper bound for the eigenvalues of the preconditioned operator is implied by the stability of the local solvers si on the local subspaces once the coarse components have been removed (this is made explicit in Lemma 2.8). This is where the GenEO strategy comes in. We solve a generalized eigenvalue problem which identifies the 'bad' modes: in this case those for which we cannot ensure that the local solver is stable for a constant independent of the coefficients in the equations. These 'bad' modes are then used to span the coarse space, and the local solvers are stable on all remaining local components (the 'good' components). More precisely, the next two definitions introduce the generalized eigenvalue problem, the coarse space and the corresponding two-level BDD-GenEO preconditioners.

Definition 2.3 (GenEO coarse space for BDD)

For each subdomain i = 1, . . . , N , find the eigenpairs

(p k i , λ k i ) ∈ W i × R + of the generalized eigenvalue problem: si (p k i , v i ) = λ k i âbb (R ⊤ i p k i , R ⊤ i v i ) for any v i ∈ W i . (2.6)
Next, given a threshold K i > 0 for each subdomain, define the coarse space as

W 0 = span{R ⊤ i p k i ; λ k i < K i , i = 1, . . . , N } ⊂ Ŵ . (2.7)
Let the interpolation operator R ⊤ 0 be the matrix whose columns are the coarse basis functions

{R ⊤ i p k i ; λ k i < K i , i = 1, . . . , N }.
Finally, let the coarse solver be the exact solver on W 0 :

S 0 := R 0 ŜR ⊤ 0 ,
and P 0 be the Ŝ-orthogonal projection operator defined by

P 0 := R ⊤ 0 S † 0 R 0 Ŝ. (2.8)
This definition gives rise to a few immediate remarks.

Remark 2.4

(i) The operator R ⊤ 0 is a mapping between the coordinates of a vector from W 0 in the set of coarse basis functions and its representation in Ŵ (range(R ⊤ 0 ) ⊂ Ŵ ). Its transpose R 0 is a restriction operator which maps an element in Ŵ to the coordinates of its l 2 projection onto W 0 in the set of coarse basis functions.

(ii) Eigenvalue 0 for eigenproblem (2.6) is associated with the kernel of si so in some sense the coarse space will take care of the fact that si is not necessarily coercive. Note that if the coarse space would include only the kernel of si , one would obtain the usual coarse grid of the BDD.

(iii) In the definition of P 0 we used a pseudo inverse S † 0 because the columns of R ⊤ 0 are not necessarily linearly independent. The pseudo inverse is defined up to an element in Ker(R ⊤ 0 ) and the specific choice of the pseudo inverse makes no difference because the application of S † 0 is followed by an application of R ⊤ 0 .

(iv) The fact that P 0 is an Ŝ-orthogonal projection can be proved easily using the definitions of P 0 and S 0 and it is equivalent to the fact that P 0 is self adjoint with respect to S 0 .

We are now ready to introduce the BDD-GenEO preconditioner. There are two ways to add the second level once that we have chosen the coarse space: either we use a deflation based preconditioner (2.10) with the preconditioned conjugate gradient (PCG) algorithm or we use the projected preconditioned conjugate gradient (PPCG) algorithm in the space range(I -P 0 ) with the projected preconditioner (2.9). Both alternatives will lead to essentially identical convergence bounds.

Definition 2.5 (Two-level preconditioners)

Recall that, according to (2.2) and (2.8), we have defined

P i = R ⊤ i S † i R i Ŝ for any i = 1, . . . , N and P 0 = R ⊤ 0 S † 0 R 0 Ŝ.
Then define the projected preconditioned operator as

P proj := N i=1 (I -P 0 ) ⊤ P i (I -P 0 ), (2.9)
and the deflation based preconditioned operator as

P def := P 0 + N i=1 (I -P 0 ) ⊤ P i (I -P 0 ). (2.10)
In the remainder of this subsection we show that the BDD-GenEO coarse space leads to an upper bound for the eigenvalues of the preconditioned operators which does not depend on the number of subdomains or the coefficients in the equations. Instead it depends on the thresholds K i which were introduced to select the coarse basis functions. First we give some properties of the family of generalized eigenvectors (Lemma 2.6). Then we use these properties to show that the local bilinear forms are stable on the deflated local subspaces (Lemma 2.8) and the upper bound follows from there (Lemma 2.10).

Lemma 2.6

For a given subdomain i = 1, . . . , N , the eigenpairs (p k i , λ k i ) of generalized eigenproblem (2.6) can be chosen so that the set

{p k i } k of eigenvectors is an orthonormal basis of W i with respect to the inner product induced âbb (R ⊤ i •, R ⊤ i •). This writes âbb (R ⊤ i p k i , R ⊤ i p k i ) = 1; and âbb (R ⊤ i p k i , R ⊤ i p k ′ i ) = 0, k = k ′ .
An orthogonality type property with respect to si (which is not necessarily coercive) also holds:

si (p k i , p k ′ i ) = 0, k = k ′ . Proof Lemma 1.11 tells us that R i Kbb R ⊤ i is positive definite on W i so we may indeed speak of a âbb (R ⊤ i • , R ⊤ i • ) orthonormal basis of W i .
Then for the proof see e.g. [START_REF] Leborgne | Valeurs et vecteurs propres: définition[END_REF].

Remark 2.7

The fact that the generalized eigenproblem (2.6) is equivalent to a non-generalized eigenproblem implies that all eigenvalues are finite. Because both matrices are symmetric positive semi definite, the eigenvalues are also non negative: for any k, 0 ≤ λ k i < +∞.

The next lemma states that the local solvers are stable and strongly relies on the definition of the GenEO coarse space. In fact the purpose of the GenEO strategy is specifically to ensure that Lemma 2.8 holds. This corresponds to Assumption 2.4 in [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF].

Lemma 2.8 (Stability of the local solvers) Suppose the pseudo inverse S †

i in Definition 2.1 is chosen such that range( S † i ) = span{p k i ; λ k i > 0}. Then for any i = 1, . . . , N , the local solvers are stable in the sense

ŝ(R ⊤ i u i , R ⊤ i u i ) ≤ 1 K i si (u i , u i ), ∀ u i ∈ range( Pi (I -P 0 )),
where the K i are the thresholds that were used to select eigenvectors for the coarse space in Definition 2.3.

Proof

We may indeed choose range(

S † i ) = span{p k i ; λ k i > 0}
because the pseudo inverse of an operator is defined up to an element in the kernel of this operator. Precisely there are an infinity of pseudo inverse and we may choose the range of S † i among all the spaces which satisfy range(

S † i ) ⊕ Ker( Si ) = W i . Here, Ker( Si ) = span{p k i ; λ k i = 0}
and the set of all p k i is a basis of W i so our choice fits this limitation.

Next we prove that range( Pi (I -

P 0 )) = range( S † i R i Ŝ(I -P 0 )) ⊂ span {p k i } K .
where we have introduced the notation {p k i } K for the set of good eigenvectors

{p k i } K = {p k i ; λ k i ≥ K i }.
We will use the following linear algebra identity:

Ker((I -P 0 ) ⊤ ŜR ⊤ i ) ⊕ ⊥ range(R i Ŝ(I -P 0 )) = W i , (2.11)
where the symbol ⊥ refers to the l 2 orthogonality between both spaces and ⊕ means that the sum is direct. By definition (2.8) of P 0 , (I -

P 0 ) ⊤ = I -ŜR ⊤ 0 S † 0 R 0 so range( ŜR ⊤ 0 ) ⊂ Ker((I -P 0 ) ⊤ ).
In particular, for a given i = 1, . . . , N : span{

ŜR ⊤ i p k i ; λ k i < K i } ⊂ Ker((I -P 0 ) ⊤ ), which implies span {p k i } K ⊂ Ker((I -P 0 ) ⊤ ŜR ⊤ i ).
(2.12)

Next we use another linear algebra identity:

W i is finite dimensional so span {p k i } K ⊕ ⊥ span{p k i ; λ k i < K i } ⊥ = W i . (2.13)
According to Lemma 2.6 the

{p k i } K form a R i Kbb R ⊤ i -orthonormal basis of W i so p k i , R i Kbb R ⊤ i p k ′ i = 0, ∀k = k ′ .
This implies that span{R

i Kbb R ⊤ i p k i ; λ k i ≥ K i } ⊂ span {p k i } K ⊥ .
The equality between these subsets follows by a dimensional argument: the set

{p k i } K forms a basis of W i and R i Kbb R ⊤ i is spd so rank{R i Kbb R ⊤ i p k i ; λ k i ≥ K i } = rank{p k i ; λ k i ≥ K i } = rank {p k i } K ⊥ ,
and in turn the inclusion becomes an equality:

span{R i Kbb R ⊤ i p k i ; λ k i ≥ K i } = span {p k i } K ⊥ .
Injecting this into (2.13) implies

span {p k i } K ⊕ ⊥ span{R i Kbb R ⊤ i p k i ; λ k i ≥ K i } = W i . ( 2
(R i Ŝ(I -P 0 )) ⊂ span{R i Kbb R ⊤ i p k i ; λ k i ≥ K i },
where the argument is:

(E 1 ⊕ ⊥ E 2 = E 3 ⊕ ⊥ E 4 and E 1 ⊂ E 3 ) ⇒ E 4 ⊂ E 2 ,
for any vector spaces E 1 , . . . , E 4 . By definition of eigenproblem (2.6),

λ k i R i Kbb R ⊤ i p k i = Si p k i so range(R i Ŝ(I -P 0 )) ⊂ span{ Si p k i ; λ k i ≥ K i }.
Finally, for the specific choice of the pseudo inverse S † i it follows that range(

S † i R i Ŝ(I -P 0 )) = range( Pi (I -P 0 )) ⊂ span {p k i } K .
Now we prove the inequality in the lemma. Any u i ∈ range( Pi (I -P 0 ))

writes u i = {k;λ k i ≥Ki} α k i p k i for some coefficients α k i ∈ R. From Lemma 1.6, it is obvious that ŝ(R ⊤ i u i , R ⊤ i u i ) ≤ âbb (R ⊤ i u i , R ⊤ i u i ) = âbb   R ⊤ i {k;λ k i ≥Ki} α k i p k i , R ⊤ i {k;λ k i ≥Ki} α k i p k i   .
Using successively the first orthogonality property in Lemma 2.6, the definition of the eigenproblem and the second orthogonality property in Lemma 2.6 we get

ŝ(R ⊤ i u i , R ⊤ i u i ) ≤ {k;λ k i ≥Ki} α k i 2 âbb (R ⊤ i p k i , R ⊤ i p k i ) = {k;λ k i ≥Ki} 1 λ k i α k i 2 si (p k i , p k i ) ≤ 1 K i {k;λ k i ≥Ki} α k i 2 si (p k i , p k i ) = 1 K i si (u i , u i ).
Remark 2.9 (Local stability, Exact solvers, and Choice of the eigenproblem)

The bilinear form on the left hand side of the inequality in the lemma is

ŝ(R ⊤ i •, R ⊤ i •)
. This is the so called exact solver on subdomain i for the global problem given by Ŝ. The exact solvers are by definition the solvers which are used to build the Additive Schwarz preconditioner. For the problem Ŝu = fΓ the Additive Schwarz preconditioner would be

N i=1 R ⊤ i ŜR i .
If these exact solvers were used instead of Si the upper bound for the eigenvalues of the preconditioned operator would depend only on a constant related to the number of neighbours (introduced in the next lemma). The nice bound that we have for the lowest eigenvalue of the preconditioned operator would no longer hold though. The most straightforward generalized eigenproblem which arises from the theory is

si (p k i , v i ) = λ k i ŝ(R ⊤ i p k i , R ⊤ i v i ) for any v i ∈ W i , (2.15)
so the eigensolve operates some sort of spectral comparison between the exact solver (on the right) and the one which we actually use (on the left). We then isolate the modes for which the chosen preconditioner is not a good enough approximation in the coarse space and use a direct solve on these modes. It is however expensive to assemble and to solve (2.15). This is is why in this article we have chosen to go through only with eigenproblem (2.6) where ŝ is replaced by âbb . For a coarse space based on Eigenproblem (2.15) the theory goes through to the exact same final estimate simply by replacing âbb by ŝ in the proofs.

The following lemma gives a consequence of the stability of the local solvers. It is very narrowly related to Lemma 2.6 in [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF].

Lemma 2.10 (Upper bound for the eigenvalues of the preconditioned operator) The stability of each of the local solvers which was proved in Lemma 2.8 implies

ŝ N i=1 P i u, u ≤ N max 1≤i≤N 1 K i ŝ(u, u) ∀ u ∈ range(I -P 0 ),
where N is the maximal number of neighbours of a subdomain (including itself) in the sense:

N := max 1≤i≤N #{j; R j R ⊤ i = 0} .

Proof

This is basically the proof of Lemma 2.6 in [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF] but where we have chosen not to rely on strengthened Cauchy Schwarz inequalities. Instead we make the number of neighbours of a subdomain appear explicitly. Let u ∈ range(I -P 0 ), then

ŝ(P i u, P i u) = ŝ(R ⊤ i Pi u, R ⊤ i Pi u) ≤ 1 K i si ( Pi u, Pi u) (Lemma 2.8) = 1 K i ŝ(u, R ⊤ i Pi u) (definition of Pi (2.3)) = 1 K i ŝ(u, P i u).
We use the fact that P i = R ⊤ i Pi and the definition of ŝ to write

ŝ(P i u, u) = N j=1 s j (R j R ⊤ i Pi , R j u) = {j;Rj R ⊤ i =0} s j (R j R ⊤ i Pi , R j u).
We apply the Cauchy Schwarz inequality first for s j then for the Euclidean inner product to this and inject the previous result (in the last step)

ŝ(P i u, u) ≤ {j;Rj R ⊤ i =0} s j (R j R ⊤ i Pi , R j R ⊤ i Pi ) 1/2 s j (R j u, R j u) 1/2 ≤   {j;Rj R ⊤ i =0} s j (R j R ⊤ i Pi , R j R ⊤ i Pi )   1/2   {j;Rj R ⊤ i =0} s j (R j u, R j u)   1/2 = ŝ(P i u, P i u) 1/2   {j;Rj R ⊤ i =0} s j (R j u, R j u)   1/2 ≤ 1 K i ŝ(u, P i u) 1/2   {j;Rj R ⊤ i =0} s j (R j u, R j u)   1/2 .
Raising to the square and simplifying by ŝ(P i u, u) yields

ŝ(P i u, u) ≤ 1 K i   {j;Rj R ⊤ i =0} s j (R j u, R j u)   .
Finally summing these inequalities over i gives the result.

Main theorem: convergence bound for BDD with the GenEO coarse space

We are now ready to give the estimates for the condition number of BDD with the GenEO coarse space.

Theorem 2.11 (Main theorem for BDD with the GenEO coarse space)

The condition number for BDD solved in range(I -P 0 ) with the projected additive operator (2.9) satisfies

κ (P proj ) ≤ N max 1≤i≤N 1 K i . ( 2

.16)

As for the condition number of the deflated operator (2.10) with the GenEO coarse space, it satisfies

κ (P def ) ≤ max 1, N max 1≤i≤N 1 K i .
(2.17)

These bounds depend only on the chosen thresholds K i which we use to select eigenvectors for the coarse space in Definition 2.3 and on the maximal number N of neighbours of a subdomain:

N = max 1≤i≤N #{j; R j R ⊤ i = 0} .

Proof

The proof of this theorem is the proof of Theorem 2.13 in [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF]. The fact that the local solvers ( S † i here) are not spd does not play a role in the proof. The idea is to prove the following bounds:

ŝ(u, u) ≤ ŝ(P proj u, u) ≤ N max 1≤i≤N 1 K i ŝ(u, u); u ∈ range(I -P 0 ), (2.18)
and

ŝ(u, u) ≤ ŝ(P def u, u) ≤ max 1, N max 1≤i≤N 1 K i ŝ(u, u); u ∈ Ŵ . (2.19)
Following Lemma C.1 in the appendix of [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF] these bounds imply the bounds for the condition numbers. They are proved using Lemma 2.2 and Lemma 2.10 combined with the fact that P 0 is an ŝ-orthogonal projection.

Remark 2.12

The fact that K i can be chosen such that N max 1≤i≤N or with notation inspired by the usual FETI notation:

G ⊤ N λ = R N f Γ , G N := BR ⊤ N . (3.5)
Injecting (3.3) into the second equation in (3.2) we get

BS † B ⊤ λ -G N α = BS † f Γ , for some α ∈ range(R N ).
We may again rewrite the problem using a saddle point formulation as

F -G N G ⊤ N 0 λ α = d e , (3.6)
where

F := BS † B ⊤ , d := BS † f Γ , e := R N f Γ , and again G N = BR ⊤ N . (3.7)
In order to homogenize the second equation and bring the problem down to a single equation we decompose λ into λ = λ + λ N where G ⊤ N λ = 0 and G ⊤ N λ N = e. Then we introduce a projection operator P N as follows: let Q : U → U be a self-adjoint matrix which is positive definite on range(G N ), then define

P N : U → U ; P N := I -QG N (G ⊤ N QG N ) -1 G ⊤ N . (3.8) Remark 3.1
It is straightforward to prove that P N is a projection operator from U onto Ker(G ⊤ N ) and that its transpose

P ⊤ N = I -G N (G ⊤ N QG N ) -1 G ⊤ N Q is a Q-orthogonal projection.
It is however less obvious to prove that the inverse (G ⊤ N QG N ) -1 is well defined. This can be derived from the fact that

Q is positive definite on range(G N ) so G ⊤ N QG N β = 0 implies G N β = 0 ⇔ BR ⊤ N β = 0.
In other words R ⊤ N β ∈ Ker(S) ∩ Ker(B) and this intersection is zero because the problem is well posed. † Finally β = 0 and (G ⊤ N QG N ) -1 is well defined. The system which we solve is the projected system into the space

V N := Ker(G ⊤ N ) = range(P N ). (3.9)
For the choice

λ N := QG N (G ⊤ N QG N ) -1 R N f Γ (which fulfills the condition G ⊤ N λ N = e) the problem is: find λ ∈ V N and α ∈ range(R N ) such that F λ -G N α = d -F λ N .
(3.10)

Testing this against elements in V N yields the final form of the problem before preconditioning

P ⊤ N F λ = P ⊤ N (d -F λ N ), (3.11)
whereas testing against function in range(I -P N ) allows us to define the component α of the solution completely with respect to λ:

(I -P ⊤ N )G N α = (I -P ⊤ N )(F λ -d + F λ N ) ⇔ α = (G ⊤ N QG N ) -1 G ⊤ N Q(F λ -d),
where we simply used a multiplication by (G ⊤ N QG N ) -1 G ⊤ N Q to write the equivalence. Next we introduce the two usual FETI preconditioners.

Usual preconditioners for FETI

We first need to introduce diagonal scaling matrices D i : W i → W i for each i = 1, . . . , N . These are the matrix counterparts of the partition of unity operators Ξ i used in the BDD section. Then space. We will consider two different preconditioners for (3.11): the Dirichlet preconditioner with the subscript D and the lumped preconditioner with the subscript L [START_REF] Farhat | A method of finite element tearing and interconnecting and its parallel solution algorithm[END_REF]. When scaled, those preconditioners can be written as the following operators on U [START_REF] Klawonn | Feti and neumann-neumann iterative substructuring methods: connections and new results[END_REF]:

M -1 D = D -1 B ⊤ (BD -1 B ⊤ ) † ⊤ S D -1 B ⊤ (BD -1 B ⊤ ) † (3.12) M -1 L = D -1 B ⊤ (BD -1 B ⊤ ) † ⊤ K bb D -1 B ⊤ (BD -1 B ⊤ ) † . (3.13)
We use the subscript * to refer to either of these preconditioners generically: if * denotes D then

M -1 * = M -1 D is the Dirichlet preconditioner and if * denotes L then M -1 * = M -1
L is the Lumped preconditioner. When the diagonal scaling matrix D is chosen to be the diagonal of the local operator matrix K, the scaling in the preconditioners (3.12,3.13) are equivalent to so-called super-lumped scaling (or K-scaling) originally proposed in [START_REF] Rixen | A simple and efficient extension of a class of substructure based preconditioners to heterogeneous structural mechanics problems[END_REF].

Remark 3.2

In (3.12,3.13) we have used a pseudo inverse where the usual FETI theory uses an inverse. This has no impact on what follows. Indeed, (BD -1 B ⊤ ) † is defined up to an additive element in Ker(BD -1 B ⊤ ) and we have the inclusion Ker(BD

-1 B ⊤ ) ⊂ Ker(B ⊤ ) since λ ∈ Ker(BD -1 B ⊤ ) ⇒ D -1 B ⊤ λ ∈ Ker(B) ⇒ B ⊤ λ = Dv for some v ∈ Ker(B), and Ker(B) = (range(B ⊤ )) ⊥ so v ⊤ B ⊤ λ = v ⊤ Dv = 0 ⇒ v = 0 ⇒ λ ∈ Ker(B ⊤ ).
The operator (BD -1 B ⊤ ) † is applied to elements in range(B) = Ker(B ⊤ ) ⊥ so this application is well defined. Moreover the application of (BD -1 B ⊤ ) † is followed by an application of B ⊤ so D -1 B ⊤ (BD -1 B ⊤ ) † is uniquely defined independently of the choice of the pseudo inverse. This pseudo inverse can be avoided by defining scaling matrices directly on the space of Lagrange multipliers which is done for instance in the redundant Lagrange multiplier section of [START_REF] Klawonn | Feti and neumann-neumann iterative substructuring methods: connections and new results[END_REF]. For sensible choices both approaches can lead to identical preconditioners and in practical implementations the scaling matrices are actually never computed explicitly as is explained in [START_REF] Rixen | A simple and efficient extension of a class of substructure based preconditioners to heterogeneous structural mechanics problems[END_REF].

Using the subscript * for either D or L, the preconditioned operator is M -1 * P ⊤ N F . Because we solve the system using a projected conjugate gradient method we require that the search directions remain in V N . Therefore we actually solve: find λ ∈ V N such that

P N M -1 * P ⊤ N F λ = P N M -1 * P ⊤ N (d -F λ N ). (3.14)
Because of the projection step (3.11) and the choice

λ N := QG N (G ⊤ N QG N ) -1 R N f
this is already a two-level preconditioner where the coarse space is Ker(P N ) = range(QG N ) = range(QBR ⊤ N ). The PPCG solver is initialized with λ N and the entire solution space is λ N + V N . We will refer to P N as the natural coarse space projector.

The theoretical study of the preconditioner is related to operator

P D : W → W ; P D := D -1 B ⊤ (BD -1 B ⊤ ) † B, (3.15)
where D : W → W is the diagonal scaling matrix already introduced. This is a projection that is orthogonal in the scaled l 2 inner product x ⊤ Dy (x, y ∈ W ). The next two lemmas follow essentially by noticing that BP D u = Bu. They are Lemmas 4.1 and 4.3 in [START_REF] Klawonn | Feti and neumann-neumann iterative substructuring methods: connections and new results[END_REF]. We give the proofs for sake of completeness because they are short.

Lemma 3.3

For any µ ∈ U there exists ũ ∈ range(P D ) such that µ = B ũ.

Proof

By definition of U there exists u ∈ W such that µ = Bu. Now take ũ = P D u, B ũ = Bu = µ.

Lemma 3.4

Let u ∈ W , then

P D u = u -E D u, (3.16)
where E D u : W → W is an averaging operator defined by its components as:

(E D u) i = R i N j=1 R ⊤ j D j u j .

Proof

We start by noticing that B(u -P D u) = 0. This means that u -P D u matches at the interfaces and thus its weighted average satisfies E D (u -P D u) = u -P D u. A sufficient condition to ensure that the result holds is now

E D P D u = 0.
By definition of E D , E D P D u is a D-weighted average of the values of P D u which correspond to the same global dof. One way to compute the averaged value for global dof k is to first compute DP D u = B ⊤ (BD -1 B ⊤ ) † Bu and then sum the contributions from the different subdomains for which k is a degree of freedom. This is the same as computing an l 2 scalar product between B ⊤ (BD -1 B ⊤ ) † Bu and the function e x ∈ W which is zero everywhere except at the degrees of freedom which correspond to global dof k. By definition Be x = 0. The orthogonality of Ker(B) and range(B ⊤ ) allows us to conclude that Be x , B ⊤ (BD -1 B ⊤ ) † Bu = 0 and thus

E D P D u = 0.
The operator G * ,0 is a mapping between the coordinates of a vector from U * ,0 in the set of coarse basis functions and its representation in U . Its transpose G ⊤ * ,0 is a restriction operator which maps an element in W to the coordinates of its l 2 projection onto W * ,0 in the set of coarse basis functions. The main difference with the coarse space for BDD is that we have left out the zero eigenvalues which correspond to the kernel of S because they are already taken care of by the natural coarse space through P N .

Remark 3.8

One common point with the BDD GenEO eigenvalue problem is that one of the operators (S i ) is a non assembled operator on the local space W i whereas the other (B ⊤ i M -1 * B i ) is an assembled operator restricted to the local space W i . This time the words assembled and restricted are to be understood in the FETI context and rely on the mappings B i between the degrees of freedom in W i and the Lagrange multipliers in U . In the same way as for BDD, the role of the GenEO eigenvalue problem for FETI can be interpreted as finding the modes necessary for describing the discrepancy between the interface behavior as seen from a single domain (left hand side of (3.19)), and the assembled interface operator F -1 , approximated by M -1 * (right hand side of (3.19)). The idea is then to introduce those differences, which will not be well accounted for by the preconditioner, into the coarse space.

Once again in proving our estimate for the condition number we will take advantage of the orthogonality type properties which result from the generalized eigenvalue problem. Lemma 3.9 Let * denote either D or L. For a given subdomain i = 1, . . . , N , the eigenpairs (q k i , Λ k i ) of the generalized eigenproblem (3.19) can be chosen so that the set {q k i } k of eigenvectors is an orthonormal basis of W i with respect to the inner product induced by B ⊤ i M -1 * B i . This writes

M -1 * B i q k i , B i q k i = 1; and M -1 * B i q k i , B i q k ′ i = 0, k = k ′ .
An orthogonality type property with respect to S i (which is not necessarily coercive) also holds:

S i q k i , q k ′ i = 0, k = k ′ .

Proof

We proved in Lemma 3.5 that M -1 * is spd on range(G N ) = Ker(P ⊤ N ). We also proved in Lemma 3.6 that M -1 * is spd on

V ′ N = range(P ⊤ N ). So M -1 * is spd on Ker(P ⊤ N ) ⊕ range(P ⊤ N ) = U . Finally by definition of B i , B i u i = 0 implies u i = 0 so B ⊤ i M -1
For two spd matrices M 1 and M 2 of same size, the spectrum of M 1 M 2 is identical to the spectrum of M 2 M 1 . Following this idea we decide to look at the problem in reverse: Is F a good preconditioner for M -1 * ? The reason why we do this is that then we recognize an abstract Schwarz type preconditioner

F = N i=1 B i S † i B ⊤ i .
In this framework, the local subspaces are the W i and the local solvers are the pseudo inverses S † i of the local bilinear forms S i . The prolongation operators are the B i : W i → U and the restriction operators are the B ⊤ i : U → W i . Taking advantage of the abstract Schwarz framework, in Lemmas 3.11 and 3.13 we will prove the same estimates as in the BDD subsection for F viewed as the preconditioner and M -1 * viewed as the matrix problem. In the proof of our final theorem it will become apparent that these estimates allow to prove the condition number of FETI with the two-level preconditioners given by (3.22). In the next Lemma, applying the exact same strategy as in Lemma 2.2 we give an estimate related to a lower bound for the eigenvalues of the preconditioned operator F P * ,N P * ,0 M -1 * . This bound does not depend on the choice of the coarse space. Lemma 3.11 (Stable splitting -Lower bound for the eigenvalues of the preconditioned operator) For any µ ∈ V ′ * ,N there exists a stable splitting

(v 1 , . . . , v N ) ∈ W 1 × • • • × W N of µ : µ = B 1 v 1 + . . . B N v N ; v i ∈ W i and N i=1 Proof Let µ ∈ V ′ * ,N and let v i = D -1 i B ⊤ i (BD -1 B ⊤ ) † µ for each i = 1, . . . , N
. This provides a splitting of µ:

N i=1 B i v i = N i=1 B i D -1 i B ⊤ i (BD -1 B ⊤ ) † µ = (BD -1 B ⊤ )(BD -1 B ⊤ ) † µ = µ, since µ ∈ range(BD -1 B ⊤ ) = range(B) = U .
Moreover, the splitting is stable:

N i=1 S i v i , v i = N i=1 S i D -1 i B ⊤ i (BD -1 B ⊤ ) † µ, D -1 i B ⊤ i (BD -1 B ⊤ ) † µ = SD -1 B ⊤ (BD -1 B ⊤ ) † µ, D -1 B ⊤ (BD -1 B ⊤ ) † µ = M -1 D µ, µ , ≤ M -1 * µ, µ ,
by Lemma 1.6. This is exactly (3.23). Now let µ ∈ range(P ⊤ * ,0 P ⊤ * ,N ), then M -1 * µ, µ = P * ,N P * ,0 M -1 * µ, µ . Moreover, the fact that the v i provide a splitting implies

M -1 * µ, µ = P * ,N P * ,0 M -1 * µ, N i=1 B i v i = N i=1 P * ,N P * ,0 M -1 * µ, B i (S † i S i )v i = N i=1 S i v i , S † i B ⊤ i P * ,N P * ,0 M -1 * µ .
Then we apply the Cauchy Schwarz inequality twice, first in the S i inner product and then in the l 2 inner product and finish by using (3.23)

M -1 * µ, µ ≤ N i=1 S i v i , v i 1/2 S i S † i B ⊤ i P * ,N P * ,0 M -1 * µ, S † i B ⊤ i P * ,N P * ,0 M -1 * µ 1/2 ≤ N i=1 S i v i , v i 1/2 N i=1 S i S † i B ⊤ i P * ,N P * ,0 M -1 * µ, S † i B ⊤ i P * ,N P * ,0 M -1 * µ 1/2 ≤ M -1 * µ, µ 1/2 P * ,N P * ,0 M -1 * µ, N i=1 B i S † i B ⊤ i P * ,N P * ,0 M -1 * µ 1/2 .
The result follows by raising to the square, simplifying by M -1 * µ, µ and recognizing

F = N i=1 B i S † i B ⊤ i .
The next lemma is the FETI counterpart of lemma 2.8 and the proof follows the exact same steps. We prove a crucial result which relies very strongly on the choice of the coarse space. In fact the coarse space was chosen specifically to ensure that this estimate holds. Lemma 3.12 (Stability of the local solvers) Let * denote either D or L. For each i = 1, . . . , N , let the pseudo inverse S † i be chosen such that range(S

† i ) = span{q k i ; Λ k i > 0}.
Then the following estimate for the local solver holds

Proof

First we prove that range(

S † i B ⊤ i M -1 * P ⊤ * ,0 P ⊤ * ,N ) ⊂ span{q k i ; Λ k i ≥ K i }.
We will use the following linear algebra identity

Ker(P * ,N P * ,0 M -1 * B i ) ⊕ ⊥ range(B ⊤ i M -1 * P ⊤ * ,0 P ⊤ * ,N ) = W i , (3.25)
where the symbol ⊥ refers to the l 2 orthogonality between both spaces and ⊕ means that the sum is direct. According to item (ii) in Lemma 3.10, P ⊤ * ,0 P ⊤ * ,N = P ⊤ * ,N P ⊤ * ,0 P ⊤ * ,N . This implies P * ,N P * ,0 = P * ,N P * ,0 P * ,N . So Ker(P * ,N ) ⊂ Ker(P * ,N P * ,0 ). It is also obvious that Ker(P * ,0 ) ⊂ Ker(P * ,N P * ,0 ). Using the definitions of these projections ((3.17 By definition of G * ,0 and G N , in particular, for each i = 1, . . . , N ,

span{M -1 * B i q k i ; Λ k i < K i } ⊂ Ker(P * ,N P * ,0 ), so span{q k i ; Λ k i < K i } ⊂ Ker(P * ,N P * ,0 M -1 * B i ). (3.26)
Following the same procedure as to prove (2.14) in Lemma 2.8, the first orthogonality property in Lemma 3.9 implies that

span{q k i ; Λ k i < K i } ⊕ ⊥ span{B ⊤ i M -1 * B i q k i ; Λ k i ≥ K i } = W i . (3.27) Putting (3.25), (3.26) and (3.27)
together tells us that range(

B ⊤ i M -1 * P ⊤ * ,0 P ⊤ * ,N ) ⊂ span{B ⊤ i M -1 * B i q k i ; Λ k i ≥ K i }.
Next the definition of eigenproblem (3.19)

, S i q k i = Λ k i (B ⊤ i M -1 * B i ) q k i , yields range(B ⊤ i M -1 * P ⊤ * ,0 P ⊤ * ,N ) ⊂ span{S i q k i ; Λ k i ≥ K i }.
Finally for the specific choice of the pseudo inverse

S † i it is obvious that range(S † i B ⊤ i M -1 * P ⊤ * ,0 P ⊤ * ,N ) ⊂ span{q k i ; Λ k i ≥ K i }.
Now it is easy to prove (3.24) using the orthogonality type properties in Lemma 3.9 and the definition of the eigenproblem. Any

u i ∈ range(S † i B ⊤ i M -1 * P ⊤ * ,0 P ⊤ * ,N ) writes u i = {k;Λ k i ≥Ki} α k i q k i for some coefficients α k i ∈ R, so: M -1 * B i u i , B i u i = {k;Λ k i ≥Ki} α k i 2 M -1 * B i q k i , B i q k i = {k;Λ k i ≥Ki} 1 Λ k i α k i 2 S i q k i , q k i ≤ 1 K i {k;Λ k i ≥Ki} α k i 2 S i q k i , q k i = 1 K i S i u i , u i
The next lemma is a direct consequence. It is the FETI counterpart of Lemma 2.10 and gives an estimate related to an upper bound for the eigenvalues of the preconditioned operator. The relationship will become apparent in the proof of the final theorem.

Lemma 3.13 (Upper bound for the eigenvalues of the preconditioned operator)

The following estimate holds

F M -1 * λ, M -1 * λ ≤ N max 1≤i≤N 1 K i M -1 * λ, λ for any λ ∈ range(P ⊤ * ,0 P ⊤ * ,N ), (3.28)
where N is the maximal number of neighbours of a subdomain (including itself) in the sense

N = max 1≤i≤N #{j; B ⊤ j B i = 0} .

Proof

In order to simplify notation lets write P * ,i := S † i B ⊤ i M -1 * and P * ,i := B i P * ,i . Let λ ∈ range(P ⊤ * ,0 P ⊤ * ,N ), then

M -1 * P * ,i λ, P * ,i λ = M -1 * B i P * ,i λ, B i P * ,i λ ≤ 1 K i S i P * ,i λ, P * ,i λ (Lemma 3.12) = 1 K i M -1 * λ, B i P * ,i λ (definition of P * ,i ) = 1 K i M -1 * λ, P * ,i λ (3.29)
Taking a close look at the definition of the preconditioners in (3.12) and (3.13) we notice that they can be written as a sum of local contributions:

M -1 * = N j=1 M -1 * ,j ; M -1 * ,j := D -1 j B ⊤ j (BD -1 B ⊤ ) † ⊤ S j D -1 j B ⊤ j (BD -1 B ⊤ ) † ,
and

M -1 * ,j B i u i , u i = 0 if and only if B ⊤ j B i = 0. A consequence of this is that M -1 * λ, P * ,i λ = M -1 * λ, B i P * ,i λ = {j;Bj B ⊤ i =0}
M -1 * ,j λ, B i P * ,i λ .

We apply the Cauchy Schwarz inequality for M -1 * ,j and then for the Euclidean inner product to this and inject the previous result

M -1 * λ, P * ,i λ ≤ {j;Bj B ⊤ i =0} M -1 * ,j λ, λ 1/2 M -1 * ,j P * ,i λ, P * ,i λ 1/2 ≤   {j;Bj B ⊤ i =0} M -1 * ,j λ, λ   1/2   {j;Bj B ⊤ i =0} M -1 * ,j P * ,i λ, P * ,i λ   1/2 =   {j;Bj B ⊤ i =0} M -1 * ,j λ, λ   1/2 M -1 * P * ,i λ, P * ,i λ 1/2 ≤   {j;Bj B ⊤ i =0} M -1 * ,j λ, λ   1/2 1 K i M -1 * λ, P * ,i λ 1/2
(from (3.29)).

Raising to the square and simplifying by M -1 * λ, P * ,i λ yields We give here a few numerical results to confirm the estimate for the condition number in the FETI case. The system of equations which we solve is related to two dimensional linear elasticity where the domain is clamped on the left hand side and subject to gravity. An important feature of the methods which we presented is that, given a FETI code, they do not demand a lot of implementation work: all the mathematical objects which are used to build the coarse space already appear in the algorithms.

All the results that follow were obtained using Freefem++ [START_REF] Hecht | Numerical Mathematics and Scientific Computation[END_REF] to build the problem matrices and visualize solutions and Matlab for the solving procedure. The test problems we present here are only small tests which we use to validate our theoretical results. Of course, a full validation of the efficiency of the method would require larger scale tests with an optimized code. Full reorthogonalization at each iteration is used in PPCG. The meshes are regular with quadrilateral elements and the finite element discretization of the two dimensional elasticity equation uses standard P 1 (linear) functions. There are two parameters in the linear elasticity system of equations: Young's modulus E and Poisson's ratio ν. Each time an iteration count is given, the stopping criterion is that the relative primal residual at the final iteration k reach 10 -4 :

N i=1 R ⊤ i S i D -1 i B ⊤ i (BD -1 B ⊤ ) † P ⊤ * ,0 P ⊤ * ,N (d -F λ k ) 2 fΓ 2 < 10 -4 .
The fact that this is indeed the primal residual is explained in [START_REF] Rixen | Extended preconditioners for FETI method applied to constrained problems[END_REF] and proved for instance in [START_REF] Mandel | An algebraic theory for primal and dual substructuring methods by constraints[END_REF].

Checkerboard coefficient distribution

We discretize a square of size 1 × 1 using 81 × 81 nodes. We use two different decompositions of this unit square: a regular decomposition into 8 × 8 regular subdomains (Figure 1 -left) and a decomposition into 64 subdomains obtained using Metis [START_REF] Karypis | METIS: A software package for partitioning unstructured graphs, partitioning meshes, and computing fill-reducing orderings of sparse matrices[END_REF] (Figure 1 -middle). Throughout this subsection, the scaling matrices are chosen to be the K-scaling matrices [START_REF] Rixen | A simple and efficient extension of a class of substructure based preconditioners to heterogeneous structural mechanics problems[END_REF][START_REF] Klawonn | Feti and neumann-neumann iterative substructuring methods: connections and new results[END_REF], meaning that in the definitions of the preconditioners (3.12) and (3.13) we set 4.1.1. The partition resolves the heterogeneities It is well known by now that in the case of a regular decomposition into subdomains which resolves the jumps in the coefficients and the Dirichlet preconditioner, the use of the K-scaling matrices (4.1) is sufficient to ensure good convergence. We check here that in these cases the (automatic) GenEO strategy is to do nothing special which is to say that no extra modes are selected to build the additional coarse space U 0 . Table II gives the results for the regular partition (Figure 1 -left) into subdomains and a constant coefficient distribution (E; ν) = (10 7 ; 0.4) as well as a checkerboard coefficient distribution (Figure 1 -right) where the coefficients take the values (E 1 ; ν 1 ) = (10 7 ; 0.4) and (E 2 ; ν 2 ) = (10 12 ; 0.3). We have solved each of these problems with the Dirichlet preconditioner and the Lumped preconditioner with and without the GenEO coarse space (we refer to these cases as FETI-GenEO and FETI-1 respectively). For each test we give the condition number κ of the preconditioned operator, the size of the GenEO coarse space #U 0 (if there is one) and the number it of iterations needed to reach convergence. The first thing that we notice is that in all four cases where the GenEO coarse space is used the estimate for the condition number is satisfied. In the Dirichlet preconditioner case, no modes where selected to build the coarse space which is what we expected since the K-scaling alone is known to be efficient. With the Lumped preconditioner case only few modes were selected (less than one per subdomain). This test indicates that the GenEO coarse grid circumvents the fact that the lumped preconditioner does not properly predict the corrections needed on the interface for checkerboard problems.

The partition does not resolve the heterogeneities

This time we use the automatic partition into 64 subdomains obtained using METIS [START_REF] Karypis | METIS: A software package for partitioning unstructured graphs, partitioning meshes, and computing fill-reducing orderings of sparse matrices[END_REF] (Figure 1 -middle). The coefficient distribution is still the checkerboard distribution shown on the right hand side of Figure 1 so the subdomain interfaces do not coincide with the jumps in the coefficients. The coefficients are a fixed (E 1 ; ν 1 ) = (10 7 ; 0.4) and a variable (E 2 ; ν 2 ) one. Table III gives the results for different values of (E 2 ; ν 2 ).

The middle line shows a case where the coefficients are constant throughout the subdomain ((E 2 ; ν 2 ) = (E 1 ; ν 1 )). Once again we observe that in all cases the condition number satisfies the estimate and that it hardly varies with the jumps in the coefficients. In the worse case the number of modes used to build the coarse space is 370 (less than 6 modes per subdomain on average). Because of bad numerical conditioning there are a few cases where the FETI-1 residual never reaches 10 -4 , instead it stagnates. In this case we report the iteration count before the plateau and the corresponding residual. Figure 2 shows a comparison between the convergence curves with and without the additional GenEO coarse space where this phenomenon can be observed. Figure 3 shows the spectrum of the preconditioned operators with and without the additional coarse space. The spectrum is represented in the complex plane but the imaginary part is always almost zero (imaginary parts result from numerical errors in the eigensolver). The zeros in the spectrum correspond to the coarse modes (either natural or GenEO) as well as the null space of B ⊤ . Whether the GenEO coarse space is used or not, the first non zero eigenvalue of the preconditioned operator is 1 which is what is expected.

Discontinuities along the interfaces

In this subsection we focus only on the GenEO coarse space for the Dirichlet preconditioner and we conduct a more extensive study. We use a partition into N regular subdomains of a rectangle of size N × b where b is the aspect ratio of each subdomain (see Figure 4). The discretization of each subdomain is n el × n el rectangular elements so that each element has the same aspect ratio as the (1) the relative residual reaches a plateau at 2 • 10 -4 after 142 iterations.

(2) the relative residual reaches a plateau at 3 • 10 -4 after 154 iterations.

(3) the relative residual reaches a plateau at 2 • 10 -3 after170 iterations. (4) the relative residual reaches a plateau at subdomain to which it belongs. The coefficient distribution consists of a constant value ν = 0.3 of Poisson's ratio and 7 layers of E (4 soft layers, 3 hard layers, see again Figure 4). Throughout this The parameters are: b = 1 (aspect ratio), n el = 21 (number of elements per direction per subdomain) and E 1 /E 2 = 10 -5 (jump in the coefficient). The spectrum is shown in Figure 5 along with the first 11 generalized eigenvectors and corresponding eigenvalues. We observe that there is a gap in the spectrum of the generalized eigenproblem after the 9-th generalized eigenvalue since λ 9 = 0.11 and λ 10 = 0.98. For this reason a judicious choice of the threshold for selecting eigenvectors which are put into the coarse space is for instance

K i = 0.15,
we will use this in all following numerical tests. With this criteria, the GenEO eigenproblem for a floating subdomain will provide 9 modes: the first three are rigid body modes included in the usual FETI natural coarse space, and 6 deformation modes that are included in the GenEO coarse grid. As can be seen in Figure 5 those deformation modes represents the behavior of the subdomain when the hard layers deform the soft ones. The 9 modes can be seen as a basis to describe the nearly rigid motion of the hard layers (3 modes for each of the 3 layers, amounting to 9 modes) and the basis spanned by those modes represent the behavior of the domain as if the hard layers were its backbones. In some sense the GenEO coarse space can be interpreted in this case as a skeleton of the overall problem describing the dominant behavior of the structure according to its hard layers.

Next we actually solve the problem for different numbers of subdomains, different aspect ratios and different discretizations. The results are shown in Table IV. The two level method with the GenEO coarse space is robust throughout all of these tests: the condition number varies between 1.34 and 4.51 only, which is indeed lower than the upper bound given by the theory, N /K i = 20, N being equal to three in this simple decomposition. Further the following observations are noteworthy:

• When the number of domains increases, the classical FETI-1 method sees its number of iteration increase significantly, whereas equipped with the GenEO coarse space, the number of iteration remains small. The dimension of the GenEO coarse spaces is roughly proportional to the number of domains in this case. • The classical FETI method convergences very slowly when the height of the domain is large compared to its width (b = 5). For that case the GenEO strategy generates only a small number of modes (43 in total) and converges very fast. • For this layered structure, the preconditioned interface problem of FETI-1 has a condition number that barely depends on the number of elements per domain, and the number of iterations is nearly invariant with respect to the discretization step. When equipped with the GenEO coarse space, a small number of modes is included in the coarse space (38 GenEO modes, independent of the discretization step), and the number of iteration is very small It is thus remarkable that the GenEO coarse space can handle automatically (once a proper threshold K has been chosen) the difficult cases of bad aspect ratios and heterogeneities along the interface.

Discontinuities along and across interfaces

In this subsection we consider the case of Figure 6 where the only difference with the previous subsection is that we have added jumps across the interfaces in subdomains 3 and 6 by inverting the soft and hard layers. The parameters are as follows: n el = 21 elements in each direction and each subdomain, N = 8 subdomains, ν = 0.3 for Poisson's ratio, E 1 /E 2 = 10 -5 for the magnitude of the jump in the coefficient, b = 1 for the aspect ratio of the subdomains and K i = 0.15 for the threshold on the GenEO eigenvalues. This is a known hard problem for FETI even with the Dirichlet preconditioner (which we use here again). In this case we show in Table V that with the K-scaling matrices (4.1) the number of bad eigenmodes is largely reduced compared to the case where multiplicity scaling is used (here multiplicity scaling reduces to setting all entries of each D i to 1/2). Indeed with K-scaling we have selected 46 modes which is only 8 more than for the same case but without the extra jumps across the interfaces (see Table IV -top -N = 8 subdomains). With the multiplicity scaling the GenEO strategy selects 173 modes. In fact, with K-scaling fewer modes are necessary because jumps across the interfaces are already accounted for in the preconditioner.

The additional modes are needed to take into account the jumps across the interfaces. This confirms that GenEO compensates for the discrepancy between the preconditioner and the actual inverse of F : when inadequate weighting is used the preconditioner is less effective and hence a larger coarse space is needed. The condition numbers for both types of scaling are almost equal when the GenEO coarse space is introduced, which confirms the theory.

CONCLUSION

We have constructed a two-level BDD method and two two-level FETI methods for which the convergence rates depend only on a chosen parameter and the maximal number of neighbours of a subdomain. The choice of this parameter is key in dimensioning the coarse space. Optimizing the choice of the parameter with respect to efficiency and the size of the coarse space is crucial. Here it has been set heuristically. For FETI the result holds for the full preconditioner based on solving Dirichlet problems in the subdomains and also on the lumped version which is a lot less expensive to implement. Compared to the Schwarz-GenEO algorithm these methods have the advantage of being non overlapping methods which means that they do not carry the extra cost of computations in the overlap.

In this paper the fundamental ideas and proofs underlying the GenEO coarse space have been explained and the numerical efficiency has been illustrated on problems hard to solve with classical FETI approaches. Future research will investigate the computational cost incurred by the GenEO coarse space (computation of the GenEO modes per domain, building and solving the coarse grid) in order to assess the overall computational efficiency of the FETI-GenEO when applied to realistic engineering problems.

Definition 1 . 3 (

 13 Local finite element spaces) For any subset D ⊂ Ω let the set of degrees of freedom in D be the set dof (D) := {k = 1, . . . , n; φ k|D = 0 |D },

Definition 1 . 8 (

 18 Partition of unity) Let µ = (µ 1 , . . . , µ N ) ∈ W be a discrete partition of unity: i=1,...,N R ⊤ i µ i = 1 | Ŵ , where 1 | Ŵ ∈ Ŵ and all vector entries are 1.

  let D : W → W be the diagonal scaling matrix D :=

  ) and(3.21)) this rewrites Ker(P * ,N P * ,0 ) ⊃ (Ker(P * ,N ) ∪ Ker(P * ,0 )) ⊃ range(G * ,0 ) ∪ range(M -1 * G N ) .

Figure 1 .

 1 Figure 1. Decomposition of the unit square into 64 regular subdomains (left) -Decomposition of the unit square into 64 subdomains using Metis (middle) -Checkerboard coefficient distribution (right)

Figure 2 .

 2 Figure 2. Checkerboard coefficient distribution -Convergence curve: primal residual versus iteration count -Left: with GenEO, Right : without GenEO -Lumped preconditioner for the Metis decomposition into 64 subdomains -(E 1 ; ν 1 ) = (10 7 ; 0.4) and (E 2 ; ν 2 ) = (10 12 ; 0.3).

Figure 3 .

 3 Figure 3. Checkerboard coefficient distribution -Spectrum of the preconditioned operator -Left: with GenEO, Right : without GenEO -Lumped preconditioner for the Metis decomposition into 64 subdomains -(E 1 ; ν 1 ) = (10 7 ; 0.4) and (E 2 ; ν 2 ) = (10 12 ; 0.3).

Figure 4 .

 4 Figure 4. Discontinuities along the interfaces

3 λ 7 = 9 . 6 • 10 -3 λ 8 = 4 . 1 • 10 -2 λ 9 = 0. 11 .Figure 5 .

 379610841109115 Figure 5. Eigenvalues and eigenmodes of the GenEO generalized eigenproblem for the geometry given in Figure 4 -dark or pink: hard material, light or yellow: soft material -The first eigenmodes (rigid body modes) are part of the natural coarse grid, and the next 6 are selected for the GenEO coarse space.

Figure 6 . 5 •

 65 Figure 6. Discontinuities across and along interfaces (subdomains 3 and 6) Table V. Geometry given in Figure 6 (discontinuities across and along the interfaces), n el = 21, N = 8, E 1 /E 2 = 10 -5 -κ : condition number; #U 0 : size of the GenEO coarse space; it: number of iterations

Table I .

 I Summary of Notations

	Function space	Description	Definition
	W h (Ω)	Global	solution space for (1.1)
	W h (Ωi) Wi W	Local Local trace Product trace	{u |Ω i ; u ∈ W h (Ω)} ((1.6); D = Ωi) {u |Γ∩∂Ω i ; u ∈ W h (Ω)} ((1.6); D = Γ ∩ ∂Ωi) W1 × . . . WN
	Ŵ	Global trace	{u |Γ ; u ∈ W h (Ω)} ((1.6) ; D = Γ)
	Stiffness matrices (defined on)	Matrix	Bilinear form
	Global (W h (Ω))	K (1.3)	â (1.1)
	Local (W h (Ωi)) Product space ( N i=1 W h (Ωi)) Lumped global (W h (Ω))	Ki K (1.11) Kbb (1.18)	aΩ i (1.7) for D = Ωi none âbb
	Lumped product space ( N i=1 W h (Ωi))	K bb (1.19)	a bb
	Schur complement (defined on)	Matrix	Bilinear form
	Global ( Ŵ )	Ŝ (1.16)	ŝ
	Local (Wi)	Si (1.13)	si (1.22)
	On the product space (W ) Weighted local (Wi)	S (1.14) Si	s si (1.23)
	Right hand sides	Notation
	Condensed onto Γ	fΓ (1.20)
	Condensed onto Γ ∩ ∂Ωi	fΓ,i (1.21)
	Condensed on product space N i=1 Γ ∩ ∂Ωi	fΓ (1.21)

Table II .

 II Checkerboard (64 regular subdomains) κ : condition number; #U 0 : size of the GenEO coarse space; it: number of iterations -For the Dirichlet preconditioner the GenEO coarse space is empty so FETI-GenEO and FETI-1 are identical

			Dirichlet			Lumped
					FETI-GenEO	FETI-1
	Coefficients	κ	#U0	it	κ	#U0	it	κ	it
	Constant	9.5	0	15	11.1	15	17	86 24
	Checkerboard	6.3	0	13	9.7	49	19	93 25

Table III .

 III Checkerboard (64 Metis subdomains) (E 1 ; ν 1 ) = (10 7 ; 0.4); κ : condition number; #U 0 : size of the GenEO coarse space; it: number of iterations. When (E 2 ; ν 2 ) = (10 7 ; 0.4) there are no jumps in the coefficients.

			Dirichlet Preconditioner			Lumped Preconditioner
		FETI-GenEO	FETI-1		FETI-GenEO	FETI-1
	(E2; ν2)	κ	#U0	it	κ	it	κ	#U0	it	κ	it
	(10 12 ; 0.3)	10.4	126	18	1.5 • 10 6 142 (1)	11.7	186	19	6.2 • 10 6 154 (2)
	(10 7 ; 0.4)	10.5	26	18	447	31	12.2	99	23	2.1 • 10 3	58
	(10 2 ; 0.49)	12.2	182	21	5.3 • 10 6 170 (3)	16.3	370	23	4.0 • 10 7 198 (4)

  1 • 10 -3 after 198 iterations.

Table IV .

 IV Three tests for the geometry in Figure4κ : condition number; #U 0 : size of the GenEO coarse space; it: number of iterations Various number of subdomains (N ), fixed aspect ratio (b = 1), fixed discretization (n el = 21),fixed jump in coefficients (E1/E2 = 10 -5 ), the problem size increases with N * ) the relative residual reaches a plateau at 2 • 10 -3 after 47 iterations.Various discretizations (n el ), fixed aspect ratios (b = 1), fixed number of subdomains (N = 8), fixed jump in coefficients (E1/E2 = 10 -5 ), the problem size increases with n el .

		FETI-GenEO	FETI-1
	N subdomains	κ	#U0 it	κ	it
	4 8 16 32 64	3 1.34 1.34 1.35 1.35	14 38 86 182 374	5 5 4 4 4	20 39 75 2.2 • 10 3 137 1.4 • 10 3 1.9 • 10 3 2.1 • 10 3 2.2 • 10 3 190
	Various aspect ratios (b), fixed number of subdomains (N = 8), fixed discretization (n el = 21),fixed jump in
	coefficients (E1/E2 = 10 -5 )
		FETI-GenEO	FETI-1
	aspect ratio b	κ	#U0	it	κ	it
	5 2 1	2.33 1.42 1.34	43 40 38	6 5 5	1.7 • 10 5 47 ( * ) 43 1.0 • 10 4 1.9 • 10 3 40
	1/2	4.51	27	9	446	33
	1/5	4.07	14	11	70	22
		FETI-GenEO	FETI-1
	n el elements	κ	#U0 it	κ	it
	21 42 70 84	1.34 1.42 1.46 1.47	38 38 38 38	5 5 5 5	1.92 • 10 3 39 1.93 • 10 3 40 1.94 • 10 3 40 1.94 • 10 3 40

(

Copyright c 0000 John Wiley & Sons, Ltd.Int. J. Numer. Meth. Engng (0000) Prepared using nmeauth.cls DOI: 10.1002/nme

† In case the global operator K is singular, a solution exists for the original problem if f is in the range of K. In that case the natural coarse grid becomes singular but the FETI approach can still be applied[START_REF] Rixen | Dual schur complement method for semi-definite problems[END_REF].

ACKNOWLEDGEMENTS

We are very grateful to Victorita Dolean, Patrice Hauret and Frédéric Nataf for many fruitful discussions and constructive comments.

Copyright c 0000 John Wiley & Sons, Ltd.

Int. J. Numer. Meth. Engng (0000) Prepared using nmeauth.cls DOI: 10.1002/nme

The FETI formulation

In the BDD section we built the coarse space for problem (2.1) which we simply recall here: find û ∈ Ŵ such that Ŝ û = fΓ , where Ŵ is the space of functions defined on the interface Γ. Instead the FETI formulation of the problem is on the product space W with an additional matching constraint at the interfaces. This constraint is ensured using matrix

which is constructed from entries 0, 1, -1 such that the components u i of a vector u in the product space W coincide on Γ when Bu = 0. More precisely each line in B corresponds to one continuity constraint for one degree of freedom and two of the subdomains to which it belongs: each line in B contains one 1 and one -1 while all other entries are zero. Denoting by λ the vector of Lagrange multipliers which is used to enforce the constraint Bu = 0 we obtain a saddle point formulation of the problem: find

We note that the solution λ of (3.2) is unique only up to an additive element of Ker(B ⊤ ) however the solution u to our problem does not depend on the choice of λ so this is not an issue in practice.

For the theoretical study we introduce the space

and will search for λ ∈ U . Given a basis for Ker(S) which consists of n K vectors, an important role is played by the prolongation operator R ⊤ N : R nK → W which columns are these basis functions. The transpose R N is a restriction operator which maps an element in W to the coordinates of its l 2orthogonal projection onto Ker(S) in the same basis. We have used the subscript N because Ker(S) is often referred to as the Natural coarse space for FETI. Going back to the system, the solution of the first equation in (3.2) can be written as

if the right-hand side associated to the operator S is such that

AUTOMATIC SPECTRAL COARSE SPACES FOR ROBUST FETI AND BDD ALGORITHMS 21

Two level FETI preconditioner with the GenEO coarse space

The proof of an upper bound for the spectrum of the preconditioned FETI system usually relies on strong assumptions on the set of equations at hand and the coefficient distribution. Once again we build a coarse space which allows us to waive all of these assumptions. The coarse space is defined next along with the two-level FETI preconditioners (projected and deflated). We use again the subscript 0 to refer to the coarse space. In order to avoid confusion with the BDD case we use calligraphic notation for the projection operator P * ,0 .

Definition 3.7 (GenEO coarse spaces for FETI) Let * denote either D (for Dirichlet) or L (for Lumped). For each subdomain i = 1, . . . , N , find the eigenpairs (q k i , Λ k i ) ∈ W i × R + of the generalized eigenvalue problem:

where M -1 * is the preconditioner defined either by (3.12) or (3.13). Next, given a threshold K i > 0 for each subdomain, define the coarse space as

Let the interpolation operator G * ,0 be the matrix whose columns are the coarse basis functions

Let the coarse solver be the exact solver on U * ,0 :

We are now ready to prove the main theorem for the GenEO FETI algorithm which is similar to Theorem 2.11. Theorem 3.14 (Main theorem for FETI with the GenEO coarse space) Let * denote either L for Lumped or D for Dirichlet. The condition number for FETI solved in range(P * ,N P * ,0 ) with the projected additive operator satisfies

.30)

As for the two-level preconditioner based on deflating the GenEO coarse space and solving in range(P * ,N ), it satisfies

(3.31) These bounds depend only on the chosen thresholds K i we use to select eigenvectors for the coarse space in Definition 3.7 and on the maximal number N of neighbours of a subdomain (including itself):

Proof From Lemma C.1 in the appendix of [START_REF] Toselli | Domain decomposition methods-algorithms and theory[END_REF], in order to prove (3.30), it is sufficient to show that, for any λ ∈ range(P * ,N P * ,0 ), the following holds: