Citation |

- Permanent Link:
- http://digital.auraria.edu/AA00003410/00001
## Material Information- Title:
- Analysis of lagrange multiplier based domain decomposition
- Creator:
- Tezaur, Radek
- Place of Publication:
- Denver, Colo.
- Publisher:
- University of Colorado Denver
- Publication Date:
- 1998
- Language:
- English
- Physical Description:
- 113 leaves : ; 28 cm
## Thesis/Dissertation Information- Degree:
- Doctorate ( Doctor of Philosophy)
- Degree Grantor:
- University of Colorado Denver
- Degree Divisions:
- Department of Mathematical and Statistical Sciences, CU Denver
- Degree Disciplines:
- Applied Mathematics
- Committee Chair:
- Mandel, Jan
- Committee Members:
- Farhat, Charbel
Franca, Leo Knyazev, Andrew Russell, Tom
## Subjects- Subjects / Keywords:
- Finite element method ( lcsh )
Boundary value problems ( lcsh ) Neumann problem ( lcsh ) Lagrangian functions ( lcsh ) Boundary value problems ( fast ) Finite element method ( fast ) Lagrangian functions ( fast ) Neumann problem ( fast ) - Genre:
- bibliography ( marcgt )
theses ( marcgt ) non-fiction ( marcgt )
## Notes- Bibliography:
- Includes bibliographical references (leaves 106-113).
- General Note:
- Department of Mathematical and Statistical Sciences
- Statement of Responsibility:
- by Radek Tezaur.
## Record Information- Source Institution:
- |University of Colorado Denver
- Holding Location:
- |Auraria Library
- Rights Management:
- All applicable rights reserved by the source institution and holding location.
- Resource Identifier:
- 40337054 ( OCLC )
ocm40337054 - Classification:
- LD1190.L622 1998d .T49 ( lcc )
## Auraria Membership |

Downloads |

## This item has the following downloads: |

Full Text |

ANALYSIS OF LAGRANGE MULTIPLIER BASED
DOMAIN DECOMPOSITION by Radek Tezaur M.S., Univerzita Karlova, Praha, Czech Republic, 1993 A thesis submitted to the University of Colorado at Denver in partial fulfillment of the requirements for the degree of Doctor of Philosophy Applied Mathematics 1998 This thesis for the Doctor of Philosophy degree by Eadek Tezaur has been approved by Jan Mandel Charbel Farhat Leo Franca Andrew Knyazev Tom Russell Date Tezaur, Eadek (Ph. D., Applied Mathematics) Analysis of Lagrange multiplier based domain decomposition Thesis directed by Professor Jan Mandel ABSTRACT The convergence of a substructuring iterative method with Lagrange multipliers known as Finite Element Tearing and Interconnecting (FETI) meth- od is analyzed in this thesis. This method, originally proposed by Farhat and Roux, decomposes finite element discretization of an elliptic boundary value problem into Neumann problems on the subdomains, plus a coarse problem for the subdomain null space components. For linear conforming elements and preconditioning by Dirichlet problems on the subdomains, the asymptotic bound on the condition number (7(1 +log(H/h))'7, where 7 = 2 or 3, is proved for a second order problem, h denoting the characteristic element size and H the size of subdomains. A similar method proposed by Park is shown to be equivalent to FETI with a special choice of some components and the bound (7(1 + log (H/h))2 on the condition number is established. Next, the origi- nal FETI method is generalized to fourth order plate bending problems. The main idea there is to enforce continuity of the transversal displacement held at the subdomain crosspoints throughout the preconditioned conjugate gradient iterations. The resulting method is shown to have a condition number that does not increase with the number of subdomains, and again grows at most 111 poly-logarithmically with the number of elements per subdomain; the condi- tion number is bounded by (7(1 + log(H/h))3. These optimal properties hold for numerous plate bending elements that are used in practice including the HCT, DKT, and a class of non-locking elements for the Reissner-Mindlin plate models. The theoretical results are confirmed by numerical experiments. This abstract accurately represents the content of the candidates thesis. I recommend its publication. Signed __________________________________________ Jan Mandel DEDICATION Dedicated to my parents. ACKNOWLEDGEMENT Foremost I would like to thank my advisor, Jan Mandel, for his encourage- ment in preparation of this thesis, his support, guidance, and time spent in discussions throughout my studies at the University of Colorado at Denver. I would also like to extend my thanks to Patrick Le Tallec for useful suggestions in early stages of this work and Charbel Farhat for computational results. CONTENTS Chapter 1. Introduction................................................. 1 2. Overview of Related Domain Decomposition Methods........................................ 7 2.1 Finite Element Approximation, Model Problem and PCG 7 2.2 Abstract Schwarz Methods............................... 12 2.3 Overlapping Domain Decomposition....................... 15 2.4 Substructuring Methods................................. 18 2.4.1 Poincare-Steklov Operators and Schur Complement . 18 2.4.2 Various Preconditioners for the Reduced Problem ... 25 2.4.3 The Neumann-Neumann and Balancing Domain Decomposition Methods................................ 27 2.4.4 The Neumann-Neumann Domain Decomposition Algorithm for Plates and Shells ..................... 29 2.4.5 Lagrange Multipliers and Poincare-Steklov Operators ........................................... 31 2.5 A Two Level Method Based on Smoothed Prolongation........................................... 33 3. Finite Element Tearing and Interconnecting (Derivation and Extensions)................. 37 3.1 Original FETI.......................................... 37 3.1.1 Problem Setting and Assumptions...................... 38 3.1.2 Derivation of the Dual Problem ...................... 42 vii 3.1.3 Saddle Point Formulation 45 3.1.4 Algebraic Properties of the Original FETI Method ... 46 3.2 Generalized FETI........................................ 48 3.2.1 Duality Derivation of the Generalized FETI.............. 50 3.2.2 Derivation of Generalized FETI Algorithm without Duality........................................ 52 3.2.3 Method Selection for Plate Problems and Other Generalizations.................................. 57 3.3 The Dirichlet Preconditioner............................ 58 3.4 Interface Formulation of Park et al..................... 59 4. Convergence Analysis....................................... 63 4.1 Preliminaries........................................... 64 4.1.1 Estimates for Pl/Ql Elements............................ 67 4.1.2 Estimates for HCT Elements.............................. 71 4.2 Abstract Analysis Framework............................. 77 4.3 Analysis of the Original FETI .......................... 81 4.3.1 Assumptions............................................. 81 4.3.2 Discrete Norm Bounds ................................... 83 4.3.3 Condition Number Estimate............................... 89 4.4 Analysis of the Method by Park.......................... 90 4.5 Convergence Estimates for Plate Bending................. 92 4.5.1 Assumptions............................................. 93 4.5.2 Discrete Norm Bounds ................................... 95 4.5.3 Condition Number Estimate............................... 99 5. Computational Results ..................................... 101 viii References 106 IX 1. Introduction In this thesis, the Finite Element Tearing and Interconnecting (FETI) method is analyzed. It is one of domain decomposition methods for solving large systems of linear equations arising from finite element discretizations of elliptic differential equations. Elliptic equations often arise in modeling physical phenomena. The Laplace equation, A0 = div V0 = / in U. (1.1) models electrostatic interactions and many other potential problems. Phys- ical quantities are often governed by systems of equations, as in the case of linear elasticity which is modeled by a system of equations for the unknown displacement vector u div A tr (e(u))S, where A and /i are called Lame constants, 8 is the Kronecker ten- sor, and e(u) is the strain, e(u) = l/2(Vtt + (Vit)T). An equation describing a plate bending problem can be obtained from the linear elasticity equation by a limit process that considers thickness of the plate to be infinitesimally small. In Kirchhoff-Love model, this yields a fourth order elliptic equation. In order to set up a well-posed problem, the equations need to be complemented by boundary conditions. A boundary value problem may be 1 reformulated in weak (variational) form [56]. Two basic kinds of boundary conditions are recognized. Essential boundary conditions, such as a Dirichlet boundary condition (f> = 0 on dQ for (1.1), need to be explicitly imposed in the weak form of the problem. Natural boundary conditions, such as a Neumann boundary condition dip/dn = 0 dQ for (1.1), are naturally incorporated into the weak form. The process of discretization may be based on a Galerkin approxima- tion. In that case, an approximate solution is sought in a finite-dimensional subspace of the space in which the weak form is posed. Finding the Galerkin approximation then requires solving a system of linear equations for the coeffi- cients of the solution relative to a basis of the subspace. Finite element spaces are widely used in this context. They are generated by basis functions that are usually polynomial on each element of a triangulation of the domain and the supports of basis functions have only small overlaps. The matrix of the discrete system has several special properties. It is symmetric and positive definite. Also, due to limited overlaps of finite element basis functions, it is sparse; the number of nonzero entries in a row of the ma- trix is a lot smaller than total number of entries in the row. The size and the condition number of this matrix are affected by fineness of the triangulation. For quasi-uniform meshes for which the finite element discretization is charac- terized by the mesh size h of the triangulation, the condition number of the matrix is of the order p for second order elliptic problems and p for fourth order elliptic problems. Two classes of methods for solving the linear system exist, direct and iterative. Direct methods, usually based on a variant of Gaussian elimination, 2 are common in everyday engineering practice. They are fairly robust; their dis- advantage lies in their memory requirements and their speed. Gaussian elimi- nation leads to fill-in and the destruction of the sparse structure of the matrix. The fill-in can be reduced by re-numbering the variables, but for problems in 3D this helps only to a certain extent. In general, the number of operations required to solve the system arising from the finite element discretization is proportional to the square of the number of unknowns or worse. Since systems arising from discretizations of problems in engineering practice constantly push available computer resources to the limits, iterative methods are often preferable for large-scale problems. Their main disadvan- tage is probably their lack of robustness when too wide a class of problems is considered. On the other hand, the memory requirements of iterative methods are typically much smaller than those of the direct methods. Usually only a small multiple of the number of nonzero entries of the original matrix needs to be stored. The number of operations can be as low as a multiple of the number of unknowns, but often the computational cost is not known in advance. The speed of convergence depends on the condition number of the matrix of the system. Since the condition number deteriorates with decreasing mesh size, so does the speed of convergence of simple iterative methods. It is therefore desirable to construct methods that overcome bad conditioning of the matrix and that are, if possible, independent of other singular perturbations such as inhomogeneities or bad Poisson ratios in the linear elasticity problem or, in the case of plate bending model, the problems arising from the thickness of a plate approaching zero and the model becoming a fourth order problem. Due to the elliptic nature of the underlying problem, this sort of stability is often possible 3 to achieve and various algorithms have been proposed. In this thesis, we will study one of domain decomposition methods. This class of methods has gained enormous popularity in the last decade. Fol- lowing a divide-and-conquer idea, domain decomposition methods divide the original problem into a number of smaller subproblems. Sometimes such a divi- sion arises from breaking up a complicated geometry. In many cases, though, it is entirely artificial. The subproblems are easier to solve because of their smaller size and often parallelism can be exploited. This is especially impor- tant with the onset of parallel computing. Domain decomposition methods can be seen from two different points of view. They may arise from separation of a physical domain into regions, where a problem is modeled by separate partial differential equations, with the interfaces between the subdomains being handled by various conditions, such as continuity. The opposite approach is to see domain decomposition methods purely as methods for solving large algebraic linear systems arising from the discretization of PDEs. In this context, the large system is subdivided into smaller problems, whose solutions can be used to produce a preconditioner for the large system. Every domain decomposition algorithm involves two principal issues. It breaks up the original problem into subproblems, that are solved by some known method, and it resolves interactions of local solutions, usually by the means of an iterative method. Many domain decomposition methods have been developed, originally, for the case of two subdomains. As the complexity of problems of interest demands splitting the problem into small enough sub- problems, a multi-domain case proves to be of more importance. Early works 4 demonstrate independence of proposed algorithms on the characteristic mesh size. However, it has been shown that to achieve independence of the number of subdomains, a coarse space needs to be introduced. Solving a small coarse space problem distributes the information about the solution globally, thus resolving global characteristics of the solution not visible by the subdomains. In the next chapter, we give an overview of several well-know domain decomposition algorithms. In particular, we will describe abstract Schwarz methods, overlapping and substructuring methods and their relationship to the finite element tearing and interconnecting method (FETI) [27, 32], the main subject of this thesis. We also briefly discuss a domain decomposition method based on smoothed prolongation [68, 64], The remaining chapters, except for Sections 3.4 and 4.4, are based on papers [50, 51, 53] by Mandel and Tezaur, and by Mandel, Tezaur and Farhat. In Chapter 3, the derivation of the original FETI method is shown. FETI tears the computational domain into non-overlapping subdomains and enforces intersubdomain continuity via Lagrange multipliers applied at the subdomain interfaces. The Lagrange multipliers are used as the unknowns and FETI formulation is obtained by solving the saddle point problem for the Lagrangian. Consistent treatment of subdomain singularities leads to a small coarse problem which is solved in each iteration of the preconditioned conjugate gradient method iteration. In Chapter 3, we also present a generalization of the original FETI algorithm. We describe its particular application to a plate bending problem. Chapter 4 is concerned with analyzing the original and generalized FETI methods. We show that the condition number of the preconditioned 5 FETI method is bounded independently of the number of subdomains and poly-logarithmically in terms of subdomain size. That is, the condition num- ber is bounded by (7(1 + log^)3, where H is the characteristic subdomain size and h the mesh size. We present a complete analysis for decompositions with crosspoints in 2D and edges and crosspoints in 3D for second order elliptic problems. Furthermore, we demonstrate that the choice of the Lagrange mul- tipliers given by Park [59, 58] assures that the condition number estimate can be improved to (7(1 + log y)2. Finally, we also show that the generalized FETI converges poly-logarithmically for a biharmonic problem. In order to illustrate the potential of the generalized FETI method, Chapter 5 summarizes some computational results. We consider the plate bending problem on a unit square and demonstrate that the condition number is almost independent of the number of subdomains and the size of the problem, as predicted by the theoretical analysis. 6 2. Overview of Related Domain Decomposition Methods 2.1 Finite Element Approximation, Model Problem and PCG In this section, we are going to summarize some concepts and algo- rithms used throughout this chapter and the rest of the thesis. We will consider the following model problem. Let fi be a bounded Lipschitz domain in IE2 or IE3. We will study the second-order elliptic equation with Dirichlet boundary condition with the matrix [a^-] symmetric, uniformly positive definite and bounded on fi. The corresponding bilinear form is then given by and is defined for all v and w in the Sobolev space if1(Q) (the space of gen- eralized functions with square integrable first derivatives). Let us denote the L2(0) inner product Au = f in fi (2,1) u = 0 on dVt where 7 By the divergence theorem (integration by parts), the problem (2.1) can be written in weak (variational) form: Find u E Hq(Q) such that a(u,v) = (f,v) (2.2) for all v in V = Hq(Q), where Hq is the completion of smooth functions with support in Q with respect to the norm in if1(Q). Since the bilinear form a(u, v) is symmetric, continuous and elliptic (coercive) [56], the problem has a unique solution by the Lax-Milgram theorem. Let us now use the standard Galerkin approximation. We look for an approximate solution of (2.2) in a finite-dimensional subspace 14 (Q) of the space V. The Galerkin approximation is the solution of the following problem: Find uh E Vh (fi) such that a(uh,vh) = ( f, vh), (2.3) for all vh E 14 (Q). This problem leads to a system of linear equation when a basis for Vh is chosen. Let i be a basis for 14(0). Assuming N uh = i=1 (2.3) becomes N YaiVuVj)ui = i=1 j = 1,...,N. Matrix K = [a((pi,(pj)\ is called the stiffness or Gramm matrix and the right hand side vector [(/, (pj)\ is called the load vector. The vector of unknowns, [iq], is also called the vector of degrees of freedom. For a wide class of approximation spaces 14, uh is a good approxima- tion of u. In the case of finite element methods, the construction of suitable 8 spaces Vh relies on triangulation, which splits the domain fi into small disjoint regions of simple geometric shape, such as triangles or quadrangles in IE2 or tetrahedrons or hexahedrons in IR3; h refers to the characteristic size of these regions. Under certain assumptions which prevent the triangulation from de- generating (e.g, the size of the angles in triangles is bounded from below so as they do not become too sharp), the finer the triangulation, the closer a finite element Galerkin solution to the exact solution, [14], Functions in a finite element space Vh usually arise from a polynomial interpolation on the elements of the triangulation. Every polynomial defined on a given region is uniquely determined by its values and perhaps also the values of its derivatives at some nodal points, usually the vertices of the region. Each function in Vh is thus determined by a set of values at nodal points, so called degrees of freedom. Simple examples of finite element spaces include spaces formed by continuous functions linear on triangle regions in IR2 or tetrahedrons in IR3 or bilinear functions on quadrangles. These spaces are being referred to as Pi and Qi, respectively, and their entries are uniquely determined by their function values at the vertices of the triangulation. Another type of finite elements we will encounter later in this thesis is HCT elements, which are C1 continuous and determined by the function values and the values of the first derivatives at the vertices of triangular regions in IR2. The standard basis for a finite element space is the one in which each of the basis functions has exactly one degree of freedom equal to 1 and the rest 0. The unknowns in the linear system arising from the discretization are di- rectly the degrees of freedom of the Galerkin approximation (hence the name). The overlaps of the supports of the basis functions are small, which causes 9 the stiffness matrix to be sparse. Due to the properties of the bilinear form, the stiffness matrix is also symmetric and positive definite. Preconditioned conjugate gradient (PCG) method is therefore often considered for solving the system. A detailed description and analysis of PCG (and other methods) can be found, e.g., in Golub and Van Loan, [38]. Here we give only a brief summary. The conjugate gradient method is based on the observation that the solution of Ax = b is the only minimum of Q(x) xTb. A sequence of vectors which converges to the solution of the system can be obtained by setting xk+l = xk + akPk> where Pk is a suitably chosen direction and vector + a^p^ minimizes Q on the line Xk + apk- In PCG, vectors pu are obtained by A-orthogonalization of residua, which leads to the following algorithm: Algorithm 1 (Conjugate Gradients) Given x0, set r0 = b Axq. For k = 1,... do: Pk Pk k Xk rl-in-i n-1 + -p-Pk-1 (Pi = r0) Pk-1 0k plApk Xk-1 + O-kPk Xk Xk-i otkApk 10 Among iteration methods the conjugate gradient method is rather exceptional. Assuming roundoff errors are not present, it gives the exact solution of the system after a finite number of iteration, which equals at most the number of unknowns. For large systems, however, performing so many iterations would be impractical. Fortunately, much smaller number of iterations is often needed to get a good enough approximate solution. This is not true if the condition number of A is large. The performance then can be improved by introducing a preconditioner M = ETE, where E is some nonsingular matrix. We point out that M defined in this way is symmetric and positive definite. Preconditioner M must have two essential properties: matrix E^TAE^1 is better conditioned than A, system Md = e can be easily solved. E^T denotes the transpose of the inverse of E. The preconditioned conjugate gradient (PCG) method is the conju- gate gradient method applied to the system E-TAE-1{Ex) = E~Tb. This leads to the following algorithm. Algorithm 2 (Preconditioned Conjugate Gradients) Given xq, set r0 = (b- Axo). For i = 1,... do: zk-i = M lrk-1 A = aT-iAc-i 11 Pk = zk^i + -p-Pk-1 (pi = Zq) Pk-l _ (h ak jr- PkAPk xk = xk-1 + akPk xk xk 1 Otk-Apk We point out that matrix M 1 is not formed explicitly and, instead, a linear system with matrix M is solved in each iteration of PCG method. Stan- dard preconditioners used in practice include diagonal scaling and incomplete factorizations (e.g., Cholesky) [38]. 2.2 Abstract Schwarz Methods Abstract Schwarz methods take their name after Herman Schwarz, a German mathematician of the past century, who used an alternating method to construct harmonic functions on regions with non-smooth boundaries. They provide an abstract framework which allows analysis of many different domain decomposition methods. Here we briefly describe the general abstract algo- rithm and refer to literature for details on analysis and practical/non-practical choices of particular methods. This presentation is based on [6] and [23]. Let V be a Hilbert space with an inner product (,). We also con- sider the inner product a(-, ) implied by the bilinear form of the variational formulation of our problem: Find u V so that a(u,v) = (f,v) Vt> e V. Let Vi, % = 0,..., n be closed subspaces of V that form a decomposi- tion of V, i.e. V W + W + W + ... + fy 12 In addition to the inner product implied by the the bilinear form a, we consider the inner products a,i(-,-) implied by some bilinear forms a* defined on V* x Vi, i = 0,... ,n. The Schwarz additive method can be then written as follows. Algorithm 3 (Additive Schwarz) Given an initial approximation u, start with k = 0 and do: 1. Find iVi E V such that adwi, Vi) = (/, Vi) a(uk, i>*), e V. 2. Deline the next iterate as n uk+1 = uk + Yswi- i=0 3. Set k = k + 1 and go to 1. The Schwarz multiplicative method is then defined by the following algorithm. Algorithm 4 (Multiplicative Schwarz) Given an initial approxi- mation u, starting with k = 0 do: 1. For / = 0.....n. (i) find iVi E Vt so that Oi(m, Vi) = (/, vt) a{uk+l/{n+l),Vi), Mvi E Vi} (ii) set yk+(i+l)/(n+l) yk+i/(n+1) + 2. Set k = k + 1 and go to step 1. The algorithms above can be combined to form hybrid methods. They are usually used as preconditioners in the conjugate gradient method. Various 13 choices of the subspaces \\ yield a wide variety of different algorithms. Since the multiplicative method above is not symmetric, it is usually symmetrized by considering the sequence of subspaces in forward and reverse order; that is F0, F,.....F, Vn, V. V ,. .... F- F0. Each sub-step of an additive and a multiplicative Schwarz method can be interpreted as an approximate projection of the error. Denoting P, the operator corresponding to the i-th sub-step of an additive method, that is the operator that satisfies a,i(PiW,v) = a(w,v) Vu e F, one step of the additive method satisfies n i=0 where ek = u uk is the error of the k-th approximation. Similarly, the reduction of the error of one step of the multiplicative method is governed by the equation ek+1 = (/ Pn)...(/ P\){I Po)ek. We refer to [25] and references included there for abstract analysis. Here we only summarize the main results. They are due to a number of authors: Dryja and Widlund [25, 24]; Nepomnyaschikh [55]; Bramble, Pasciak, Wang, Xu [9, 70]; Lions [44]; Bjprstad and Mandel [6], and others. Let Cq be the constant such that for all u E V, there exists a repre- sentation a = Y'i u 11 > ui ^ Vi, and n ^2ai(ui,Ui) < C0 a(u, u). i=0 14 Let B = [by] be the matrix of strengthened Cauchy-Schwarz coefficients, Ia(vi,Vj)\ < bija(vi,vi)1/2a(vj,vj)1/2 e V^vj G Vj,i,j = l,...,n. Furthermore, let u be the constant such that a(ui, Ui) ^ ui a,i{u,h Ui) 'tfu-i G Vj, i 0,..., n. Then the smallest eigenvalue of the operator Y,Pi of the abstract additive method used a preconditioner in the conjugate gradients method is bounded by 1/Co and the largest by u(p(B) +1). The bound on the condition number of the symmetric multiplicative method used as a preconditioner in the conjugate gradients method is C0(l + 2ujp(B))/(2 cD), where Cj = max(l,o;). The space V0 is usually regarded as a coarse space. It has been shown that without the presence of the coarse space, the condition number of an abstract Schwarz method grows with the number of subspaces. The coarse space ensures a mechanism of global exchange of information otherwise lacking in decompositions involving a large number of subdomains. 2.3 Overlapping Domain Decomposition The first known domain decomposition method is due to Schwarz [63] and is known as the Schwarz alternating method. It divides the domain fi, on which the problem is defined, into two overlapping regions fii and fi2 and performs a multiplicative Schwarz type algorithm. Let us consider the model problem (2.1) Au = f in fi, u = 0 on 15 Let fii and fi2 be overlapping subdomains of fi, fii U fi2 = fh Furthermore, let I'i = (){11 fl ilj and T2 = dfl2 H (-h- Algorithm 5 (Alternating Schwarz) Given an initial guess m2 for the values on fi2, for k = 1, 2, ... do: 1. Solve Auk = / in f^!, uk = 0 on \ r uk = |n on Tj 2. Solve Auk = / in fl2, uk = 0 on cKl2 \ T: uk _ r,.k\ i|r2 on r2. To be precise, f|ri and if2|r2 stand in the algorithm above for the trace of u\ on Id and uk on r2, respectively, rather than the simple restrictions. Approximate solutions of the model system can be assembled as uk = nt in tt2, uk = uk in fii \ fl2. The variational formulation of the algorithm above is due to Lions [44], It is the multiplicative abstract Schwarz algorithm using the subspaces Vi and V2 of V = Hi (fi) formed by functions that vanish outside of fii and fl2, respectively. The bilinear forms ai and a2 are simple restrictions of a. 16 Choosing some other bilinear forms ai and a2 can model the situation when the problems in step 1. and 2. of the algorithm above were to be solved only approximately. An additive version of the algorithm has been described by Dryja [20] and Matsokin and Nepomnyaschikh [54], In the discrete, multidomain case, overlapping methods often start with a non-overlapping partition of Q into subdomains (substructures) fij, i = 1,,NS, the diameter of the subdomains being of order H. Then, each subdomain fij is extended to a larger region fi'. It is assumed that neither Vti nor fi' cuts across any of the elements. The overlaps of subdomains fi', i = 1,..., Ns are said to be generous, if the distance between boundaries of fij and fi' is greater than some fixed fraction of H. The subspaces Vt, i = 1,..., n of the abstract Schwarz method are then chosen as restrictions of the finite element space 14(fi) to the subregions fi'. That is, \; = C,n //>( Approximate or exact solvers then can be chosen on the subspaces and addi- tive, multiplicative or hybrid two-level Schwarz methods can be used. As far as parallelization is concerned, additive versions clearly gain an edge, because they can be easily parallelized. Multiplicative versions always involve sequen- tial steps. Some parallelization is possible when mutually non-overlapping subdomains fi' are grouped together. This process is often called coloring of subdomains. Without a coarse space, the condition number of such methods is of the order 1 /H2 [24], where H is the diameter of the subdomains. Dryja and 17 Widlund in [24] propose to use the space of continuous piecewise linear func- tions on the coarse mesh defined by the subdomain Q* as the coarse space V0. They prove that, in the case of generous overlap, the condition number relevant for the conjugate gradient iteration is uniformly bounded. Mandel shows in [47] that a similar method by Cowsar [16], which uses discrete harmonic ex- tensions determined by piecewise constant values on interfaces of subdomains, leads to a poly-logarithmic bound on the condition number. In their later paper [21], Dryja and Widlund mention, that their nu- merical experiments indicate that the convergence rate is often satisfactory even for small overlaps. Running time is often smallest when the overlaps are at minimum. The number of conjugate gradient iterations is higher in such a case, because the condition number deteriorates, but this is being compen- sated for by the fact that the local problems are smaller and therefore cheaper to solve. Also, they are better conditioned and, if they are solved by iteration solvers, the rate of convergence is faster. 2.4 Substructuring Methods 2.4.1 Poincare-Steklov Operators and Schur Complement Substructuring methods borrow their name form structural engineer- ing. It was in the context of structural engineering that several substructuring algorithms have been pioneered. We shall concentrate on substructuring meth- ods on interfaces of subdomains. Let us consider the model problem (2.1) and let fi be divided into two non-overlapping subdomains fii and fi2 (Figure 2.4.1). Our goal is to solve the model problem only on the subdomains. Let T = dQ\ fi <9fi2 be the subdomain 18 interface. It is well known that the solution u E Hq(Q) exists for all f E L2(0). The solution satisfies continuity conditions on the subdomain interface: conti- nuity of fluxes of the solution and continuity of the solution. Let us introduce a few notations. For k = 1,2, we define the space Hjj(Qk) as the subspace of such that the functions in H}j(Vtk) vanish on dVtkndVt. Let Hpl(Q,k) be the dual of H}j(Vtk). Let jk be the trace operator mapping functions in H}j(Vtk) to their traces on T. Let hH2(T) be the fractional order Sobolev space on T consisting of traces of functions in HlD(Vtk) and let (hIq2{T))' denote its dual. Using the continuity of fluxes, we will split the problem into two subproblems, k = 1,2: Auk = fk in ttk n ([aijjVuk) = (-1 )k+l(f on T uk = 0 on dkl fl dQk. The vector n is the normal vector to T oriented, for example, from fii to fi2- In variational terms, these are the problems of finding uk E = 1,2 19 such that a.k(uk,Vk)= g*Vk+ fkVk, Vufe 6 HlD(Qk), (2.4) JV J where ak{v,w) = Y, / in i,j=1 dv dw hik dxi dxj a,i At this point it would be possible to design an iterative method working on the subproblems using the continuity conditions. Usually, however, we first reduce the problem to a problem on the subdomain interface using Poincare-Steklov operators. First, we will look for the unknown Neumann data g* on T. We define the Poincare-Steklov operators Qk : (hIq2{T))' ifdt)2(r),Â£; = 1,2 by Qk9* = IkUki (2.5) where, for g* e is the solution of (2.4) with fk = 0. Such uk is the harmonic function satisfying the Neumann condition given by g*. In other words, the Poincare-Steklov operator maps the Neumann boundary condition into the corresponding Dirichlet boundary condition duk Qk ' ----^ TWfc- on Furthermore, we define Rk : H]J1(Qk) Hq12(T), k = 1,2 by the equation Rkfk where for fk e L2(f2fe), uk is the solution of (2.4) with g* = 0. In terms of the Poincare-Steklov operators, the problem is to find the solution g* such that (Qi + Q2)g* = R2f2-Rifi (2.6) 20 That is to find the Neumann data g* on T such that the traces of the solutions Uk, k = 1, 2 of (2.4) coincide on T. The second possibility is to enforce continuity of the solution on the boundary a priori. This yields a dual formulation to (2.6). Consider the Dirichlet problems, k = 1,2, Auk = fk in klk uk = 9 on r uk = 0 on dn n dttk We will look for the unknown Dirichlet data g on T such that the fluxes are continuous for the solutions Uk,k = 1,2 of the problems above. Using the Poincare-Steklov operators, this problem can be written as (Or1 + orOs = or'-Ri/i + Q^Rih- (2.7) The equations above and iterative methods for solving them have been studied, for example, by Agoshkov [2] and in a mixed-method setting by Glowinski and Wheller [36]. The paper by Bakhvalov and Knyazev [4] is concerned with highly discontinuous coefficients between the subdomains. As the case of two subdomains is hardly of practical importance, we conclude this section by introducing a multi-domain discrete analog of (2.7). An analog of (2.6) is described in the next chapter. Let us recall now the concept of Schur complement. Definition 6 (Schur complement) Let A be an n x n matrix and a, fi C {l,...,n} be index sets. We denote A(a, (3) the submatrix that lies in the rows of A indexed by a and the columns indexed by /3, and A{a!, (Â¥) the sub- matrix given by deleting the rows indexed by a and the columns given by (3. 21 Let A(a, a) be nonsingular. Then the matrix A(ex, ex ) A(ex, ex^jA^ex^ o;) A^cXj ex ) is called the Schur complement of A(a, a) in A. We note that the Schur complement is the matrix that arises when eliminating the unknowns x{a) from the equation Ax = b. Then, b(ex') A(ex', ex)A(ex, a:)_1&(a:) = (A(ex ex ) A(ex cx^)A(cX] ex'} A^ex^ ex }}x(ex ). We also observe that the Schur complement S of A(ex, a) in A has the property (Sy,y) inf (Ax,x). xmn,x(r)=y (2.8) Let Q be a domain in IR2 or IR3 decomposed into Ns non-overlapping subdomains fii, , &ns- Let iq be the vector of degrees of freedom for subdomain Qj corresponding to a conforming hnite element discretization of the second order elliptic problem (2.1) defined on fi, such that each subdomain is a union of some of the elements. Let u,. !\ ,. and /*, be the vector of degrees of freedom, the local stiffness matrix, and the load vectors, respectively, associated with the subdomain Q*. Let Li denote the zero-one assembly matrix mapping the subdomain degrees of freedom iq into global degrees of freedom u, that is iq = Lju. The stiffness matrix then is Ns K = jTLlKlLj i=1 and the load vector Ns / = Â£Â£./. i=1 22 The problem to be solved, Ku = /, can be reduced to an interface problem by splitting the degrees of freedom into interface and interior degrees of freedom. Let us assume that the interior degrees of freedom are listed first. The subdomain stiffness matrices Kt and the restriction matrices Lt can be split accordingly: Ui Ui Ki kf Ki = Ui ki k Li \Li, Li By eliminating the interior degrees of freedom, we obtain the Schur comple- ments of Ki in Ki ~ ' (2.9) Si = Ki^ KjKS K' The problem then reduces to the interface problem Su = f for the interface unknowns u with the global Schur complement Ns s = y,l1s1zj i=1 and the interface right hand side Ns f = Y^Li{fi- KiKf'fi). i=1 In most cases, explicit computation of Si would be too expensive, but for gradient-descent methods, only evaluation of the action of Si is necessary. This evaluation can be performed efficiently by factorizing and it corre- sponds to solving a Dirichlet problem on every subdomain. The subdomain problems can be solved in parallel. 23 Let us demonstrate that the action of S) can be evaluated by solving a Dirichlet problem on the subdomain Q*. Consider the discretized problem on the subdomain Ki kT Ui 0 Ki ki Ui fi where the Dirichlet data tq is known. From the first equation, we find -tq = l\; ll\jd,. This the solution of the Dirichlet problem given by the data tq. Then, we substitute into the second equation and obtain {Ki ~ kikr'kDiii = SiUi = Many preconditioners for the reduced problem are based on compu- tations of the action of Sf1. This can be evaluated by solving a Neumann problem on the subdomain, for unknown tq and -tq and given rq, k kT i Ui 0 k k Ui 9i where rq corresponds to the Neumann data. Then, tq = Sil9i- Indeed, from the first equation, we find -tq = lI\ jH,. Substituting into the second equation, we find {Ki ~ kik^kT)ui = SiUi = g,. We note that the Schur complement Si is the discrete analog of the inverse of the Poincare-Steklov operator Q^1 : u|r |^, and the inverse of the Schur complement is the analog of Qit : u|r (2.5). One of the advantages of the reduced problem is that it improves the condition number of the original problem. For our model problem, assuming 24 the Poincare-Friedrichs inequality \\u/||^ 2 n. < cp\v/l\\2Cl. holds, the subdomain stiffness matrix Kt satisfies r\U,h\\2 < < r\Uih\\2 HI ui 111,2,Hi 22: ai 22: Ly\\ai 111,2,0,) where u/ is the finite element function corresponding to the vector of degrees of freedom iq. Then, the subdomain Schur complement S) satisfies (cf. Lemma 26 and Lemma 33) c\\u h ||2 i 11-5,2,(90, < U TSm, < C\\uh t 11 ^ ,2,(902; This implies that, for a second order problem, the condition number of the reduced problem is of order 1/h for triangulations of characteristic mesh size h, while the condition number of the original problem is of the order 1/h2. 2.4.2 Various Preconditioners for the Reduced Problem Instead of solving the reduced problem directly, iterative methods such as the preconditioned conjugate gradient method can be applied in which only matrix vector products are required. Since each iteration is quite ex- pensive (it requires solving a Dirichlet problem on each subdomain), efficient preconditioning is important to keep the number of iterations small. Diagonal or block diagonal preconditioning of the reduced problem requires knowing the diagonal of S. Instead of computing it directly, Chan [13] proposes a boundary probing technique to construct an approximate diagonal by evaluating actions of S on carefully selected vectors. For geometrically simple subdomain in 2D, Bramble, Pasciak and Schatz in [10] propose a method for preconditioning the original problem by 25 first splitting the functions in the finite element space into discrete harmon- ics on each subdomain and functions which vanish on the boundary of each subdomain. This is equivalent to eliminating the interior degrees of freedom as in the reduced system. Then the space of discrete harmonic functions is divided into functions linear on interfaces and functions whose values are zero at subdomain crosspoints. The preconditioner is then based on a bilinear form defined on the space splitting. The condition number of the method is shown to be bounded by (7(1 + log H/h)2. A similar idea is exploited in [3] for the p-version of finite elements and a bound (7(1 + logp)2 is obtained. A different algorithm [8] by Bramble, Pasciak and Schatz is described for the case of two subdomains fii and fi2 separated by the interface T. They propose to precondition the original problem in the following way: First split the local solution
homogeneous equation with the zero boundary condition on <9fi2 and a discrete |

Full Text |

PAGE 1 ANAL YSIS OF LA GRANGE MUL TIPLIER BASED DOMAIN DECOMPOSITION b y Radek T ezaur M.S., Univ erzita Karlo v a, Praha, Czec h Republic, 1993 A thesis submitted to the Univ ersit y of Colorado at Den v er in partial fulllmen t of the requiremen ts for the degree of Doctor of Philosoph y Applied Mathematics 1998 PAGE 2 This thesis for the Doctor of Philosoph y degree b y Radek T ezaur has been appro v ed b y Jan Mandel Charbel F arhat Leo F ranca Andrew Kn y azev T om Russell Date PAGE 3 T ezaur, Radek (Ph. D., Applied Mathematics) Analysis of Lagrange m ultiplier based domain decomposition Thesis directed b y Professor Jan Mandel ABSTRA CT The con v ergence of a substructuring iterativ e method with Lagrange m ultipliers kno wn as Finite Elemen t T earing and In terconnecting (FETI) method is analyzed in this thesis. This method, originally proposed b y F arhat and Roux, decomposes nite elemen t discretization of an elliptic boundary v alue problem in to Neumann problems on the subdomains, plus a coarse problem for the subdomain n ull space componen ts. F or linear conforming elemen ts and preconditioning b y Diric hlet problems on the subdomains, the asymptotic bound on the condition n um ber C (1+log( H=h )) r where r = 2 or 3, is pro v ed for a second order problem, h denoting the c haracteristic elemen t size and H the size of subdomains. A similar method proposed b y P ark is sho wn to be equiv alen t to FETI with a special c hoice of some componen ts and the bound C (1 + log( H=h )) 2 on the condition n um ber is established. Next, the original FETI method is generalized to fourth order plate bending problems. The main idea there is to enforce con tin uit y of the transv ersal displacemen t eld at the subdomain crosspoin ts throughout the preconditioned conjugate gradien t iterations. The resulting method is sho wn to ha v e a condition n um ber that does not increase with the n um ber of subdomains, and again gro ws at most iii PAGE 4 poly-logarithmicall y with the n um ber of elemen ts per subdomain; the condition n um ber is bounded b y C (1 + log( H=h )) 3 These optimal properties hold for n umerous plate bending elemen ts that are used in practice including the HCT, DKT, and a class of non-loc king elemen ts for the Reissner-Mindlin plate models. The theoretical results are conrmed b y n umerical experimen ts. This abstract accurately represen ts the con ten t of the candidate's thesis. I recommend its publication. Signed Jan Mandel iv PAGE 5 DEDICA TION Dedicated to m y paren ts. PAGE 6 A CKNO WLEDGEMENT F oremost I w ould lik e to thank m y advisor, Jan Mandel, for his encouragemen t in preparation of this thesis, his support, guidance, and time spen t in discussions throughout m y studies at the Univ ersit y of Colorado at Den v er. I w ould also lik e to extend m y thanks to P atric k Le T allec for useful suggestions in early stages of this w ork and Charbel F arhat for computational results. PAGE 7 CONTENTS Chapter 1. In troduction . . . . . . . . . . . . . . 1 2. Ov erview of Related Domain Decomposition Methods . . . . . . . . . . . 7 2.1 Finite Elemen t Appro ximation, Model Problem and PCG 7 2.2 Abstract Sc h w arz Methods . . . . . . . . 12 2.3 Ov erlapping Domain Decomposition . . . . . . 15 2.4 Substructuring Methods . . . . . . . . . 18 2.4.1 P oincar e-Steklo v Operators and Sc h ur Complemen t . 18 2.4.2 V arious Preconditioners for the Reduced Problem . 25 2.4.3 The Neumann-Neumann and Balancing Domain Decomposition Methods . . . . . . . . . 27 2.4.4 The Neumann-Neumann Domain Decomposition Algorithm for Plates and Shells . . . . . . 29 2.4.5 Lagrange Multipliers and P oincar e-Steklo v Operators . . . . . . . . . . . . 31 2.5 A Tw o Lev el Method Based on Smoothed Prolongation . . . . . . . . . . . . 33 3. Finite Elemen t T earing and In terconnecting (Deriv ation and Extensions) . . . . . 37 3.1 Original FETI . . . . . . . . . . . . 37 3.1.1 Problem Setting and Assumptions . . . . . . 38 3.1.2 Deriv ation of the Dual Problem . . . . . . 42 vii PAGE 8 3.1.3 Saddle P oin t F orm ulation . . . . . . . . 45 3.1.4 Algebraic Properties of the Original FETI Method . 46 3.2 Generalized FETI . . . . . . . . . . . 48 3.2.1 Dualit y Deriv ation of the Generalized FETI . . . 50 3.2.2 Deriv ation of Generalized FETI Algorithm without Dualit y . . . . . . . . . . . 52 3.2.3 Method Selection for Plate Problems and Other Generalizations . . . . . . . . . 57 3.3 The Diric hlet Preconditioner . . . . . . . . 58 3.4 In terface F orm ulation of P ark et al. . . . . . . 59 4. Con v ergence Analysis . . . . . . . . . . . 63 4.1 Preliminaries . . . . . . . . . . . . 64 4.1.1 Estimates for P1/Q1 Elemen ts . . . . . . . 67 4.1.2 Estimates for HCT Elemen ts . . . . . . . 71 4.2 Abstract Analysis F ramew ork . . . . . . . . 77 4.3 Analysis of the Original FETI . . . . . . . 81 4.3.1 Assumptions . . . . . . . . . . . . 81 4.3.2 Discrete Norm Bounds . . . . . . . . . 83 4.3.3 Condition Num ber Estimate . . . . . . . 89 4.4 Analysis of the Method b y P ark . . . . . . . 90 4.5 Con v ergence Estimates for Plate Bending . . . . 92 4.5.1 Assumptions . . . . . . . . . . . . 93 4.5.2 Discrete Norm Bounds . . . . . . . . . 95 4.5.3 Condition Num ber Estimate . . . . . . . 99 5. Computational Results . . . . . . . . . . . 101 viii PAGE 9 References . . . . . . . . . . . . . . . . . 106 ix PAGE 10 1. In troduction In this thesis, the Finite Elemen t T earing and In terconnecting (FETI) method is analyzed. It is one of domain decomposition methods for solving large systems of linear equations arising from nite elemen t discretizations of elliptic dieren tial equations. Elliptic equations often arise in modeling ph ysical phenomena. The Laplace equation, = div r = f in n ; (1.1) models electrostatic in teractions and man y other poten tial problems. Ph ysical quan tities are often go v erned b y systems of equations, as in the case of linear elasticit y whic h is modeled b y a system of equations for the unkno wn displacemen t v ector u div ( u ) = f : In the isotropic case, the stress tensor can be written as ( u ) = 2 ( u ) + tr ( ( u )) where and are called Lam e constan ts, is the Kronec k er tensor, and ( u ) is the strain, ( u ) = 1 = 2( r u + ( r u ) T ). An equation describing a plate bending problem can be obtained from the linear elasticit y equation b y a limit process that considers thic kness of the plate to be innitesimally small. In Kirc hho-Lo v e model, this yields a fourth order elliptic equation. In order to set up a w ell-posed problem, the equations need to be complemen ted b y boundary conditions. A boundary v alue problem ma y be 1 PAGE 11 reform ulated in w eak (v ariational) form [56]. Tw o basic kinds of boundary conditions are recognized. Essen tial boundary conditions, suc h as a Diric hlet boundary condition = 0 on @ n for (1.1), need to be explicitly imposed in the w eak form of the problem. Natural boundary conditions, suc h as a Neumann boundary condition @=@ n = 0 @ n for (1.1), are \naturally" incorporated in to the w eak form. The process of discretization ma y be based on a Galerkin appro ximation. In that case, an appro ximate solution is sough t in a nite-dimensional subspace of the space in whic h the w eak form is posed. Finding the Galerkin appro ximation then requires solving a system of linear equations for the coecien ts of the solution relativ e to a basis of the subspace. Finite elemen t spaces are widely used in this con text. They are generated b y basis functions that are usually polynomial on eac h elemen t of a triangulation of the domain and the supports of basis functions ha v e only small o v erlaps. The matrix of the discrete system has sev eral special properties. It is symmetric and positiv e denite. Also, due to limited o v erlaps of nite elemen t basis functions, it is sparse; the n um ber of nonzero en tries in a ro w of the matrix is a lot smaller than total n um ber of en tries in the ro w. The size and the condition n um ber of this matrix are aected b y neness of the triangulation. F or quasi-uniform meshes for whic h the nite elemen t discretization is c haracterized b y the mesh size h of the triangulation, the condition n um ber of the matrix is of the order 1 h 2 for second order elliptic problems and 1 h 4 for fourth order elliptic problems. Tw o classes of methods for solving the linear system exist, direct and iterativ e. Direct methods, usually based on a v arian t of Gaussian elimination, 2 PAGE 12 are common in ev eryda y engineering practice. They are fairly robust; their disadv an tage lies in their memory requiremen ts and their speed. Gaussian elimination leads to ll-in and the destruction of the sparse structure of the matrix. The ll-in can be reduced b y re-n um bering the v ariables, but for problems in 3D this helps only to a certain exten t. In general, the n um ber of operations required to solv e the system arising from the nite elemen t discretization is proportional to the square of the n um ber of unkno wns or w orse. Since systems arising from discretizations of problems in engineering practice constan tly push a v ailable computer resources to the limits, iterativ e methods are often preferable for large-scale problems. Their main disadv antage is probably their lac k of robustness when too wide a class of problems is considered. On the other hand, the memory requiremen ts of iterativ e methods are t ypically m uc h smaller than those of the direct methods. Usually only a small m ultiple of the n um ber of nonzero en tries of the original matrix needs to be stored. The n um ber of operations can be as lo w as a m ultiple of the n um ber of unkno wns, but often the computational cost is not kno wn in adv ance. The speed of con v ergence depends on the condition n um ber of the matrix of the system. Since the condition n um ber deteriorates with decreasing mesh size, so does the speed of con v ergence of simple iterativ e methods. It is therefore desirable to construct methods that o v ercome bad conditioning of the matrix and that are, if possible, independen t of other singular perturbations suc h as inhomogeneities or bad P oisson ratios in the linear elasticit y problem or, in the case of plate bending model, the problems arising from the thic kness of a plate approac hing zero and the model becoming a fourth order problem. Due to the elliptic nature of the underlying problem, this sort of stabilit y is often possible 3 PAGE 13 to ac hiev e and v arious algorithms ha v e been proposed. In this thesis, w e will study one of domain decomposition methods. This class of methods has gained enormous popularit y in the last decade. F ollo wing a divide-and-conquer idea, domain decomposition methods divide the original problem in to a n um ber of smaller subproblems. Sometimes suc h a division arises from breaking up a complicated geometry In man y cases, though, it is en tirely articial. The subproblems are easier to solv e because of their smaller size and often parallelism can be exploited. This is especially importan t with the onset of parallel computing. Domain decomposition methods can be seen from t w o dieren t poin ts of view. They ma y arise from separation of a ph ysical domain in to regions, where a problem is modeled b y separate partial dieren tial equations, with the in terfaces bet w een the subdomains being handled b y v arious conditions, suc h as con tin uit y The opposite approac h is to see domain decomposition methods purely as methods for solving large algebraic linear systems arising from the discretization of PDE's. In this con text, the large system is subdivided in to smaller problems, whose solutions can be used to produce a preconditioner for the large system. Ev ery domain decomposition algorithm in v olv es t w o principal issues. It breaks up the original problem in to subproblems, that are solv ed b y some kno wn method, and it resolv es in teractions of local solutions, usually b y the means of an iterativ e method. Man y domain decomposition methods ha v e been dev eloped, originally for the case of t w o subdomains. As the complexit y of problems of in terest demands splitting the problem in to small enough subproblems, a m ulti-domain case pro v es to be of more importance. Early w orks 4 PAGE 14 demonstrate independence of proposed algorithms on the c haracteristic mesh size. Ho w ev er, it has been sho wn that to ac hiev e independence of the n um ber of subdomains, a coarse space needs to be in troduced. Solving a small coarse space problem distributes the information about the solution globally th us resolving global c haracteristics of the solution not visible b y the subdomains. In the next c hapter, w e giv e an o v erview of sev eral w ell-kno w domain decomposition algorithms. In particular, w e will describe abstract Sc h w arz methods, o v erlapping and substructuring methods and their relationship to the nite elemen t tearing and in terconnecting method (FETI) [27, 32 ], the main subject of this thesis. W e also briery discuss a domain decomposition method based on smoothed prolongation [68, 64]. The remaining c hapters, except for Sections 3.4 and 4.4, are based on papers [50 51, 53] b y Mandel and T ezaur, and b y Mandel, T ezaur and F arhat. In Chapter 3, the deriv ation of the original FETI method is sho wn. FETI tears the computational domain in to non-o v erlapping subdomains and enforces in tersubdomain con tin uit y via Lagrange m ultipliers applied at the subdomain in terfaces. The Lagrange m ultipliers are used as the unkno wns and FETI form ulation is obtained b y solving the saddle poin t problem for the Lagrangian. Consisten t treatmen t of subdomain singularities leads to a small coarse problem whic h is solv ed in eac h iteration of the preconditioned conjugate gradien t method iteration. In Chapter 3, w e also presen t a generalization of the original FETI algorithm. W e describe its particular application to a plate bending problem. Chapter 4 is concerned with analyzing the original and generalized FETI methods. W e sho w that the condition n um ber of the preconditioned 5 PAGE 15 FETI method is bounded independen tly of the n um ber of subdomains and poly-logarithmicall y in terms of subdomain size. That is, the condition n umber is bounded b y C (1+log H h ) 3 where H is the c haracteristic subdomain size and h the mesh size. W e presen t a complete analysis for decompositions with crosspoin ts in 2D and edges and crosspoin ts in 3D for second order elliptic problems. F urthermore, w e demonstrate that the c hoice of the Lagrange m ultipliers giv en b y P ark [59, 58] assures that the condition n um ber estimate can be impro v ed to C (1+log H h ) 2 Finally w e also sho w that the generalized FETI con v erges poly-logarithmicall y for a biharmonic problem. In order to illustrate the poten tial of the generalized FETI method, Chapter 5 summarizes some computational results. W e consider the plate bending problem on a unit square and demonstrate that the condition n um ber is almost independen t of the n um ber of subdomains and the size of the problem, as predicted b y the theoretical analysis. 6 PAGE 16 2. Ov erview of Related Domain Decomposition Methods 2.1 Finite Elemen t Appro ximation, Model Problem and PCG In this section, w e are going to summarize some concepts and algorithms used throughout this c hapter and the rest of the thesis. W e will consider the follo wing model problem. Let n be a bounded Lipsc hitz domain in I R 2 or I R 3 W e will study the second-order elliptic equation with Diric hlet boundary condition A u = f in n u = 0 on @ n ; (2.1) where A u = d X i;j =1 @ @x i a ij @u @x j : with the matrix [ a ij ] symmetric, uniformly positiv e denite and bounded on n. The corresponding bilinear form is then giv en b y a ( v; w ) = d X i;j =1 Z n a ij @v @x i @w @x j and is dened for all v and w in the Sobolev space H 1 (n) (the space of generalized functions with square in tegrable rst deriv ativ es). Let us denote the L 2 (n) inner product ( v; w ) = Z n vw: 7 PAGE 17 By the div ergence theorem (in tegration b y parts), the problem (2.1) can be written in w eak (v ariational) form: Find u 2 H 1 0 (n) suc h that a ( u; v ) = ( f; v ) (2.2) for all v in V = H 1 0 (n), where H 1 0 is the completion of smooth functions with support in n with respect to the norm in H 1 (n). Since the bilinear form a ( u; v ) is symmetric, con tin uous and elliptic (coerciv e) [56], the problem has a unique solution b y the Lax-Milgram theorem. Let us no w use the standard Galerkin appro ximation. W e look for an appro ximate solution of (2.2) in a nite-dimensional subspace V h (n) of the space V The Galerkin appro ximation is the solution of the follo wing problem: Find u h 2 V h (n) suc h that a ( u h ; v h ) = ( f; v h ) ; (2.3) for all v h 2 V h (n). This problem leads to a system of linear equation when a basis for V h is c hosen. Let f g N i =1 be a basis for V h (n). Assuming u h = N X i =1 u i i ; (2.3) becomes N X i =1 a ( i ; j ) u i = ( f; j ) ; j = 1 ; : : :; N Matrix K = [ a ( i ; j )] is called the stiness or Gramm matrix and the righ t hand side v ector [( f; j )] is called the load v ector. The v ector of unkno wns, [ u i ], is also called the v ector of degrees of freedom. F or a wide class of appro ximation spaces V h u h is a good appro ximation of u In the case of nite elemen t methods, the construction of suitable 8 PAGE 18 spaces V h relies on triangulation, whic h splits the domain n in to small disjoin t regions of simple geometric shape, suc h as triangles or quadrangles in I R 2 or tetrahedrons or hexahedrons in I R 3 ; h refers to the c haracteristic size of these regions. Under certain assumptions whic h prev en t the triangulation from degenerating (e.g, the size of the angles in triangles is bounded from belo w so as they do not become too sharp), the ner the triangulation, the closer a nite elemen t Galerkin solution to the exact solution, [14]. F unctions in a nite elemen t space V h usually arise from a polynomial in terpolation on the elemen ts of the triangulation. Ev ery polynomial dened on a giv en region is uniquely determined b y its v alues and perhaps also the v alues of its deriv ativ es at some nodal poin ts, usually the v ertices of the region. Eac h function in V h is th us determined b y a set of v alues at nodal poin ts, so called degrees of freedom. Simple examples of nite elemen t spaces include spaces formed b y con tin uous functions linear on triangle regions in I R 2 or tetrahedrons in I R 3 or bilinear functions on quadrangles. These spaces are being referred to as P 1 and Q 1 respectiv ely and their en tries are uniquely determined b y their function v alues at the v ertices of the triangulation. Another t ype of nite elemen ts w e will encoun ter later in this thesis is HCT elemen ts, whic h are C 1 con tin uous and determined b y the function v alues and the v alues of the rst deriv ativ es at the v ertices of triangular regions in I R 2 The standard basis for a nite elemen t space is the one in whic h eac h of the basis functions has exactly one degree of freedom equal to 1 and the rest 0. The unkno wns in the linear system arising from the discretization are directly the degrees of freedom of the Galerkin appro ximation (hence the name). The o v erlaps of the supports of the basis functions are small, whic h causes 9 PAGE 19 the stiness matrix to be sparse. Due to the properties of the bilinear form, the stiness matrix is also symmetric and positiv e denite. Preconditioned conjugate gradien t (PCG) method is therefore often considered for solving the system. A detailed description and analysis of PCG (and other methods) can be found, e.g., in Golub and V an Loan, [38 ]. Here w e giv e only a brief summary The conjugate gradien t method is based on the observ ation that the solution of Ax = b is the only minim um of Q ( x ) = 1 2 x T Ax x T b: A sequence of v ectors whic h con v erges to the solution of the system can be obtained b y setting x k +1 = x k + k p k ; where p k is a suitably c hosen direction and v ector x k + k p k minimizes Q on the line x k + p k In PCG, v ectors p k are obtained b y A -orthogonalization of residua, whic h leads to the follo wing algorithm: Algorithm 1 (Conjugate Gradien ts) Giv en x 0 set r 0 = b Ax 0 : F or k = 1 ; : : : do: k = r T k 1 r k 1 p k = r k 1 + k k 1 p k 1 ( p 1 = r 0 ) k = k p T k Ap k x k = x k 1 + k p k r k = r k 1 k Ap k 10 PAGE 20 Among iteration methods the conjugate gradien t method is rather exceptional. Assuming roundo errors are not presen t, it giv es the exact solution of the system after a nite n um ber of iteration, whic h equals at most the n um ber of unkno wns. F or large systems, ho w ev er, performing so man y iterations w ould be impractical. F ortunately m uc h smaller n um ber of iterations is often needed to get a good enough appro ximate solution. This is not true if the condition n um ber of A is large. The performance then can be impro v ed b y in troducing a preconditioner M = E T E where E is some nonsingular matrix. W e poin t out that M dened in this w a y is symmetric and positiv e denite. Preconditioner M m ust ha v e t w o essen tial properties: matrix E T AE 1 is better conditioned than A system Md = e can be easily solv ed. E T denotes the transpose of the in v erse of E The preconditioned conjugate gradien t (PCG) method is the conjugate gradien t method applied to the system E T AE 1 ( Ex ) = E T b: This leads to the follo wing algorithm. Algorithm 2 (Preconditioned Conjugate Gradien ts) Giv en x 0 set r 0 = ( b Ax 0 ) : F or i = 1 ; : : : do: z k 1 = M 1 r k 1 k = r T k 1 z k 1 11 PAGE 21 p k = z k 1 + k k 1 p k 1 ( p 1 = z 0 ) k = k p T k Ap k x k = x k 1 + k p k r k = r k 1 k Ap k W e poin t out that matrix M 1 is not formed explicitly and, instead, a linear system with matrix M is solv ed in eac h iteration of PCG method. Standard preconditioners used in practice include diagonal scaling and incomplete factorizations (e.g., Cholesky) [38]. 2.2 Abstract Sc h w arz Methods Abstract Sc h w arz methods tak e their name after Herman Sc h w arz, a German mathematician of the past cen tury who used an alternating method to construct harmonic functions on regions with non-smooth boundaries. They pro vide an abstract framew ork whic h allo ws analysis of man y dieren t domain decomposition methods. Here w e briery describe the general abstract algorithm and refer to literature for details on analysis and practical/non-practical c hoices of particular methods. This presen tation is based on [6 ] and [23 ]. Let V be a Hilbert space with an inner product h ; i W e also consider the inner product a ( ; ) implied b y the bilinear form of the v ariational form ulation of our problem: Find u 2 V so that a ( u; v ) = h f; v i 8 v 2 V: Let V i ; i = 0 ; : : :; n be closed subspaces of V that form a decomposition of V i.e. V = V 0 + V 1 + V 2 + : : : + V n : 12 PAGE 22 In addition to the inner product implied b y the the bilinear form a w e consider the inner products a i ( ; ) implied b y some bilinear forms a i dened on V i V i ; i = 0 ; : : :; n The Sc h w arz additiv e method can be then written as follo ws. Algorithm 3 (Additiv e Sc h w arz) Giv en an initial appro ximation u 0 start with k = 0 and do: 1. Find w i 2 V i suc h that a i ( w i ; v i ) = h f; v i i a ( u k ; v i ) ; 8 v i 2 V i : 2. Dene the next iterate as u k +1 = u k + n X i =0 w i : 3. Set k = k + 1 and go to 1. The Sc h w arz m ultiplicativ e method is then dened b y the follo wing algorithm. Algorithm 4 (Multiplicativ e Sc h w arz) Giv en an initial appro ximation u 0 starting with k = 0 do: 1. F or i = 0 ; : : :; n (i) nd w i 2 V i so that a i ( w i ; v i ) = h f; v i i a ( u k + i= ( n +1) ; v i ) ; 8 v i 2 V i ; (ii) set u k +( i +1) = ( n +1) = u k + i= ( n +1) + w i : 2. Set k = k + 1 and go to step 1. The algorithms abo v e can be com bined to form h ybrid methods. They are usually used as preconditioners in the conjugate gradien t method. V arious 13 PAGE 23 c hoices of the subspaces V i yield a wide v ariet y of dieren t algorithms. Since the m ultiplicativ e method abo v e is not symmetric, it is usually symmetrized b y considering the sequence of subspaces in forw ard and rev erse order; that is V 0 ; V 1 ; : : :; V n 1 ; V n ; V n ; V n 1 ; : : :; V 1 ; V 0 Eac h sub-step of an additiv e and a m ultiplicativ e Sc h w arz method can be in terpreted as an appro ximate projection of the error. Denoting P i the operator corresponding to the i -th sub-step of an additiv e method, that is the operator that satises a i ( P i w; v ) = a ( w; v ) 8 v 2 V i ; one step of the additiv e method satises e k +1 = ( I n X i =0 P i ) e k ; where e k = u u k is the error of the k -th appro ximation. Similarly the reduction of the error of one step of the m ultiplicativ e method is go v erned b y the equation e k +1 = ( I P n ) : : : ( I P 1 )( I P 0 ) e k : W e refer to [25] and references included there for abstract analysis. Here w e only summarize the main results. They are due to a n um ber of authors: Dryja and Widlund [25 24 ]; Nepomn y asc hikh [55]; Bram ble, P asciak, W ang, Xu [9, 70]; Lions [44 ]; Bjrstad and Mandel [6 ], and others. Let C 0 be the constan t suc h that for all u 2 V there exists a represen tation u = P n i =0 u i ; u i 2 V i and n X i =0 a i ( u i ; u i ) C 0 a ( u; u ) : 14 PAGE 24 Let B = [ b ij ] be the matrix of strengthened Cauc h y-Sc h w arz coecien ts, j a ( v i ; v j ) j b ij a ( v i ; v i ) 1 = 2 a ( v j ; v j ) 1 = 2 8 v i 2 V i ; 8 v j 2 V j ; i; j = 1 ; : : :; n: F urthermore, let be the constan t suc h that a ( u i ; u i ) a i ( u i ; u i ) 8 u i 2 V i ; i = 0 ; : : :; n: Then the smallest eigen v alue of the operator P P i of the abstract additiv e method used a preconditioner in the conjugate gradien ts method is bounded b y 1 =C 0 and the largest b y ( ( B )+1). The bound on the condition n um ber of the symmetric m ultiplicativ e method used as a preconditioner in the conjugate gradien ts method is C 0 (1 + 2^ ( B )) = (2 ^ ), where ^ = max(1 ; ). The space V 0 is usually regarded as a coarse space. It has been sho wn that without the presence of the coarse space, the condition n um ber of an abstract Sc h w arz method gro ws with the n um ber of subspaces. The coarse space ensures a mec hanism of global exc hange of information otherwise lac king in decompositions in v olving a large n um ber of subdomains. 2.3 Ov erlapping Domain Decomposition The rst kno wn domain decomposition method is due to Sc h w arz [63] and is kno wn as the Sc h w arz alternating method. It divides the domain n, on whic h the problem is dened, in to t w o o v erlapping regions n 1 and n 2 and performs a m ultiplicativ e Sc h w arz t ype algorithm. Let us consider the model problem (2.1) A u = f in n ; u = 0 on @ n : 15 PAGE 25 Let n 1 and n 2 be o v erlapping subdomains of n, n 1 [ n 2 = n. F urthermore, let 1 = @ n 1 \ n 2 and 2 = @ n 2 \ n 1 Algorithm 5 (Alternating Sc h w arz) Giv en an initial guess u 0 2 for the v alues on n 2 for k = 1 ; 2 ; : : : do: 1. Solv e A u k 1 = f in n 1 ; u k 1 = 0 on @ n 1 n 1 ; u k 1 = u k 1 2 j 1 on 1 ; 2. Solv e A u k 2 = f in n 2 ; u k 2 = 0 on @ n 2 n 2 ; u k 2 = u k 1 j 2 on 2 : T o be precise, u k 1 j 1 and u k 2 j 2 stand in the algorithm abo v e for the trace of u k 1 on 1 and u k 2 on 2 respectiv ely rather than the simple restrictions. Appro ximate solutions of the model system can be assem bled as u k = u k 2 in n 2 ; u k = u k 1 in n 1 n n 2 : The v ariational form ulation of the algorithm abo v e is due to Lions [44]. It is the m ultiplicativ e abstract Sc h w arz algorithm using the subspaces V 1 and V 2 of V = H 1 0 (n) formed b y functions that v anish outside of n 1 and n 2 respectiv ely The bilinear forms a 1 and a 2 are simple restrictions of a 16 PAGE 26 Choosing some other bilinear forms a 1 and a 2 can model the situation when the problems in step 1. and 2. of the algorithm abo v e w ere to be solv ed only appro ximately An additiv e v ersion of the algorithm has been described b y Dryja [20] and Matsokin and Nepomn y asc hikh [54]. In the discrete, m ultidomain case, o v erlapping methods often start with a non-o v erlapping partition of n in to subdomains (substructures) n i i = 1 ; : : :; N s the diameter of the subdomains being of order H Then, eac h subdomain n i is extended to a larger region n 0 i It is assumed that neither n i nor n 0 i cuts across an y of the elemen ts. The o v erlaps of subdomains n 0 i i = 1 ; : : :; N s are said to be generous, if the distance bet w een boundaries of n i and n 0 i is greater than some xed fraction of H The subspaces V i ; i = 1 ; : : :; n of the abstract Sc h w arz method are then c hosen as restrictions of the nite elemen t space V h (n) to the subregions n 0 i That is, V i = V h \ H 1 0 (n 0 i ) ; i = 1 ; : : :; N s : Appro ximate or exact solv ers then can be c hosen on the subspaces and additiv e, m ultiplicativ e or h ybrid t w o-lev el Sc h w arz methods can be used. As far as parallelization is concerned, additiv e v ersions clearly gain an edge, because they can be easily parallelized. Multiplicativ e v ersions alw a ys in v olv e sequential steps. Some parallelization is possible when m utually non-o v erlapping subdomains n 0 i are grouped together. This process is often called coloring of subdomains. Without a coarse space, the condition n um ber of suc h methods is of the order 1 =H 2 [24], where H is the diameter of the subdomains. Dryja and 17 PAGE 27 Widlund in [24] propose to use the space of con tin uous piecewise linear functions on the coarse mesh dened b y the subdomain n i as the coarse space V 0 They pro v e that, in the case of generous o v erlap, the condition n um ber relev an t for the conjugate gradien t iteration is uniformly bounded. Mandel sho ws in [47] that a similar method b y Co wsar [16], whic h uses discrete harmonic extensions determined b y piecewise constan t v alues on in terfaces of subdomains, leads to a poly-logarithmic bound on the condition n um ber. In their later paper [21], Dryja and Widlund men tion, that their n umerical experimen ts indicate that the con v ergence rate is often satisfactory ev en for small o v erlaps. Running time is often smallest when the o v erlaps are at minim um. The n um ber of conjugate gradien t iterations is higher in suc h a case, because the condition n um ber deteriorates, but this is being compensated for b y the fact that the local problems are smaller and therefore c heaper to solv e. Also, they are better conditioned and, if they are solv ed b y iteration solv ers, the rate of con v ergence is faster. 2.4 Substructuring Methods 2.4.1 P oincar e-Steklo v Operators and Sc h ur Complemen t Substructuring methods borro w their name form structural engineering. It w as in the con text of structural engineering that sev eral substructuring algorithms ha v e been pioneered. W e shall concen trate on substructuring methods on in terfaces of subdomains. Let us consider the model problem (2.1) and let n be divided in to t w o non-o v erlapping subdomains n 1 and n 2 (Figure 2.4.1). Our goal is to solv e the model problem only on the subdomains. Let = @ n 1 \ @ n 2 be the subdomain 18 PAGE 28 W 2 W G 1 nFigure 2.1. Model problem in terface. It is w ell kno wn that the solution u 2 H 1 0 (n) exists for all f 2 L 2 (n). The solution satises con tin uit y conditions on the subdomain in terface: con tin uit y of ruxes of the solution and con tin uit y of the solution. Let us in troduce a few notations. F or k = 1 ; 2, w e dene the space H 1 D (n k ) as the subspace of H 1 (n k ) suc h that the functions in H 1 D (n k ) v anish on @ n k \ @ n. Let H 1 D (n k ) be the dual of H 1 D (n k ). Let r k be the trace operator mapping functions in H 1 D (n k ) to their traces on Let H 1 = 2 00 () be the fractional order Sobolev space on consisting of traces of functions in H 1 D (n k ) and let ( H 1 = 2 00 ()) 0 denote its dual. Using the con tin uit y of ruxes, w e will split the problem in to t w o subproblems, k = 1 ; 2: A u k = f k in n k n ([ a ij ] r u k ) = ( 1) k +1 g ? on u k = 0 on @ n \ @ n k : The v ector n is the normal v ector to orien ted, for example, from n 1 to n 2 In v ariational terms, these are the problems of nding u k 2 H 1 D (n k ) ; k = 1 ; 2 19 PAGE 29 suc h that a k ( u k ; v k ) = Z g ? v k + Z n k f k v k ; 8 v k 2 H 1 D (n k ) ; (2.4) where a k ( v; w ) = d X i;j =1 Z n k a ij @v @x i @w @x j : A t this poin t it w ould be possible to design an iterativ e method w orking on the subproblems using the con tin uit y conditions. Usually ho w ev er, w e rst reduce the problem to a problem on the subdomain in terface using P oincar e-Steklo v operators. First, w e will look for the unkno wn Neumann data g ? on W e dene the P oincar e-Steklo v operators Q k : ( H 1 = 2 00 ()) 0 H 1 = 2 00 () ; k = 1 ; 2 b y Q k g ? = r k u k ; (2.5) where, for g ? 2 ( H 1 = 2 00 ()) 0 u k is the solution of (2.4) with f k = 0. Suc h u k is the harmonic function satisfying the Neumann condition giv en b y g ? In other w ords, the P oincar e-Steklo v operator maps the Neumann boundary condition in to the corresponding Diric hlet boundary condition Q k : @u k @ n r k u k : F urthermore, w e dene R k : H 1 D (n k ) H 1 = 2 00 () ; k = 1 ; 2 b y the equation R k f k = r k u k ; where for f k 2 L 2 (n k ), u k is the solution of (2.4) with g ? = 0. In terms of the P oincar e-Steklo v operators, the problem is to nd the solution g ? suc h that ( Q 1 + Q 2 ) g ? = R 2 f 2 R 1 f 1 : (2.6) 20 PAGE 30 That is to nd the Neumann data g ? on suc h that the traces of the solutions u k ; k = 1 ; 2 of (2.4) coincide on The second possibilit y is to enforce con tin uit y of the solution on the boundary a priori. This yields a dual form ulation to (2.6). Consider the Diric hlet problems, k = 1 ; 2, A u k = f k in n k u k = g on u k = 0 on @ n \ @ n k : W e will look for the unkno wn Diric hlet data g on suc h that the ruxes are con tin uous for the solutions u k ; k = 1 ; 2 of the problems abo v e. Using the P oincar e-Steklo v operators, this problem can be written as ( Q 1 1 + Q 1 2 ) g = Q 1 1 R 1 f 1 + Q 1 2 R 2 f 2 : (2.7) The equations abo v e and iterativ e methods for solving them ha v e been studied, for example, b y Agoshk o v [2] and in a mixed-method setting b y Glo winski and Wheller [36]. The paper b y Bakh v alo v and Kn y azev [4] is concerned with highly discon tin uous coecien ts bet w een the subdomains. As the case of t w o subdomains is hardly of practical importance, w e conclude this section b y in troducing a m ulti-domain discrete analog of (2.7). An analog of (2.6) is described in the next c hapter. Let us recall no w the concept of Sc h ur complemen t. Denition 6 (Sc h ur complemen t) Let A be an n n matrix and ; f 1 ; : : :; n g be index sets. W e denote A ( ; ) the submatrix that lies in the ro ws of A indexed b y and the columns indexed b y and A ( 0 ; 0 ) the submatrix giv en b y deleting the ro ws indexed b y and the columns giv en b y 21 PAGE 31 Let A ( ; ) be nonsingular. Then the matrix A ( 0 ; 0 ) A ( 0 ; ) A ( ; ) 1 A ( ; 0 ) is called the Sc h ur complemen t of A ( ; ) in A W e note that the Sc h ur complemen t is the matrix that arises when eliminating the unkno wns x ( ) from the equation Ax = b Then, b ( 0 ) A ( 0 ; ) A ( ; ) 1 b ( ) = ( A ( 0 ; 0 ) A ( 0 ; ) A ( ; ) 1 A ( ; 0 )) x ( 0 ) : W e also observ e that the Sc h ur complemen t S of A ( ; ) in A has the propert y h Sy; y i = inf x 2 I R n ;x ( 0 )= y h Ax; x i : (2.8) Let n be a domain in I R 2 or I R 3 decomposed in to N s non-o v erlapping subdomains n 1 ; n 2 ; :::; n N s Let u i be the v ector of degrees of freedom for subdomain n i corresponding to a conforming nite elemen t discretization of the second order elliptic problem (2.1) dened on n, suc h that eac h subdomain is a union of some of the elemen ts. Let u i K i and f i be the v ector of degrees of freedom, the local stiness matrix, and the load v ectors, respectiv ely associated with the subdomain n i Let L i denote the zero-one assem bly matrix mapping the subdomain degrees of freedom u i in to global degrees of freedom u that is u i = L T i u The stiness matrix then is K = N s X i =1 L i K i L T i and the load v ector f = N s X i =1 L i f i : 22 PAGE 32 The problem to be solv ed, Ku = f; can be reduced to an in terface problem b y splitting the degrees of freedom in to in terface and in terior degrees of freedom. Let us assume that the in terior degrees of freedom are listed rst. The subdomain stiness matrices K i and the restriction matrices L i can be split accordingly: u i = 2 6 6 4 u i u i 3 7 7 5 ; K i = 2 6 6 4 K i ~ K T i ~ K i K i 3 7 7 5 ; L i = [ L i ; L i ] : By eliminating the in terior degrees of freedom, w e obtain the Sc h ur complemen ts of K i in K i S i = K i ~ K i K 1 i ~ K T i : (2.9) The problem then reduces to the in terface problem S u = f for the in terface unkno wns u with the global Sc h ur complemen t S = N s X i =1 L i S i L T i and the in terface righ t hand side f = N s X i =1 L i ( f i ~ K i K 1 i f i ) : In most cases, explicit computation of S i w ould be too expensiv e, but for gradien t-descen t methods, only ev aluation of the action of S i is necessary This ev aluation can be performed ecien tly b y factorizing K i and it corresponds to solving a Diric hlet problem on ev ery subdomain. The subdomain problems can be solv ed in parallel. 23 PAGE 33 Let us demonstrate that the action of S i can be ev aluated b y solving a Diric hlet problem on the subdomain n i Consider the discretized problem on the subdomain 2 6 6 4 K i ~ K T i ~ K i K i 3 7 7 5 2 6 6 4 u i u i 3 7 7 5 = 2 6 6 4 0 f i 3 7 7 5 ; where the Diric hlet data u i is kno wn. F rom the rst equation, w e nd u i = K 1 i ~ K T i u i This the solution of the Diric hlet problem giv en b y the data u i Then, w e substitute in to the second equation and obtain ( K i ~ K i K 1 i ~ K T i ) u i = S i u i = f i : Man y preconditioners for the reduced problem are based on computations of the action of S 1 i This can be ev aluated b y solving a Neumann problem on the subdomain, for unkno wn u i and u i and giv en g i 2 6 6 4 K i ~ K T i ~ K i K i 3 7 7 5 2 6 6 4 u i u i 3 7 7 5 = 2 6 6 4 0 g i 3 7 7 5 ; where g i corresponds to the Neumann data. Then, u i = S 1 i g i Indeed, from the rst equation, w e nd u i = K 1 i ~ K T i u i Substituting in to the second equation, w e nd ( K i ~ K i K 1 i ~ K T i ) u i = S i u i = g i : W e note that the Sc h ur complemen t S i is the discrete analog of the in v erse of the P oincar e-Steklo v operator Q 1 i : u j @u @ and the in v erse of the Sc h ur complemen t is the analog of Q i : @u @ u j (2.5). One of the adv an tages of the reduced problem is that it impro v es the condition n um ber of the original problem. F or our model problem, assuming 24 PAGE 34 the P oincar e-F riedric hs inequalit y k u h i k 2 1 ; 2 ; n i c p j u h i j 2 1 ; 2 ; n i holds, the subdomain stiness matrix K i satises c k u h i k 2 1 ; 2 ; n i u T i K i u i C k u h i k 2 1 ; 2 ; n i ; where u h i is the nite elemen t function corresponding to the v ector of degrees of freedom u i Then, the subdomain Sc h ur complemen t S i satises (cf. Lemma 26 and Lemma 33) c k u h i k 2 1 2 ; 2 ;@ n i u T i S i u i C k u h i k 2 1 2 ; 2 ;@ n i : This implies that, for a second order problem, the condition n um ber of the reduced problem is of order 1 =h for triangulations of c haracteristic mesh size h while the condition n um ber of the original problem is of the order 1 =h 2 2.4.2 V arious Preconditioners for the Reduced Problem Instead of solving the reduced problem directly iterativ e methods suc h as the preconditioned conjugate gradien t method can be applied in whic h only matrix v ector products are required. Since eac h iteration is quite expensiv e (it requires solving a Diric hlet problem on eac h subdomain), ecien t preconditioning is importan t to k eep the n um ber of iterations small. Diagonal or bloc k diagonal preconditioning of the reduced problem requires kno wing the diagonal of S Instead of computing it directly Chan [13] proposes a \boundary probing" tec hnique to construct an appro ximate diagonal b y ev aluating actions of S on carefully selected v ectors. F or geometrically simple subdomain in 2D, Bram ble, P asciak and Sc hatz in [10 ] propose a method for preconditioning the original problem b y 25 PAGE 35 rst splitting the functions in the nite elemen t space in to discrete harmonics on eac h subdomain and functions whic h v anish on the boundary of eac h subdomain. This is equiv alen t to eliminating the in terior degrees of freedom as in the reduced system. Then the space of discrete harmonic functions is divided in to functions linear on in terfaces and functions whose v alues are zero at subdomain crosspoin ts. The preconditioner is then based on a bilinear form dened on the space splitting. The condition n um ber of the method is sho wn to be bounded b y C (1 + log H=h ) 2 A similar idea is exploited in [3] for the p-v ersion of nite elemen ts and a bound C (1 + log p ) 2 is obtained. A dieren t algorithm [8] b y Bram ble, P asciak and Sc hatz is described for the case of t w o subdomains n 1 and n 2 separated b y the in terface They propose to precondition the original problem in the follo wing w a y: First split the local solution on n 2 in to the componen t P that satises the nonhomogeneous equation with the zero boundary condition on @ n 2 and a discrete harmonic function H that satises the corresponding homogeneous equation. A Diric hlet problem is solv ed on n 2 to obtain P This is follo w ed b y a solution of a mixed Neumann-Diric hlet problem on n 1 using the Neumann data obtained b y solving the problem on n 2 Finally another Diric hlet problem is solv ed to obtain H a discrete harmonic function, that satises the corresponding homogeneous equation. This is sometimes called Neumann-Diric hlet decomposition. This algorithm has also been described in [24 ] and it is in fact preconditioning the reduced system b y either S 1 1 or S 1 2 The in v erses of the Sc h ur complemen ts need not be explicitly computed; the action corresponds to solving a Neumann (or mixed Neumann-Diric hlet) problem on one subdomain and then extending the solution to the other subdomain b y solving a Diric hlet 26 PAGE 36 problem with the Diric hlet data found from the Neumann problem. The idea of the Neumann-Diric hlet preconditioner can be extended to m ultiple subdomain case when there is a red-blac k ordering of the subdomains. Sev eral methods solving an equation with the global Sc h ur complemen t are proposed b y Mandel [47]. The methods are set in the space of all discrete harmonic functions on the union of subdomain in terfaces = [ i @ n i @ n. The methods are h ybrid Sc h w arz methods using the coarse space V 0 in a m ultiplicativ e fashion and other spaces are treated additiv ely A method of this kind can also be in terpreted as a t w o lev el v ariational m ultigrid. The methods use V i that are associated with subdomains or globs. A glob is a v ertex, an edge that does not con tain its endpoin ts, or a face of a subdomain in terface. Piecewise linear functions dened on the subdomain triangulation or a piecewise constan t space based on glob-wise a v eraging of subdomain v alues is used as the coarse space. P oly-logarithmic bounds are obtained. Man y other c hoices of the coarse space are discussed b y Dryja, Smith and Widlund [23 ]. Their exhaustiv e in v estigation comprises v ertex based coarse spaces, based on piecewise linear functions on substructures used as elemen ts, wire bask et algorithms, that use a v erages on substructures, and face based algorithms. 2.4.3 The Neumann-Neumann and Balancing Domain Decomposition Methods The balancing domain decomposition (BDD) is based on the so called Neumann-Neumann preconditioner that preconditions the reduced problem b y a w eigh ted sum of in v erses of the local Sc h ur complemen ts. It is called 27 PAGE 37 Neumann-Neumann because it corresponds to solving Neumann problems only on subdomains (as opposed to the Diric hlet-Neumann preconditioner, that uses both Diric hlet and Neumann problems). Since this is a preconditioner in a w a y dual to the method that is the subject of this thesis, w e describe it in a little more detail. T o describe the preconditioner, w eigh t matrices D i ; i = 1 ; : : :; N s satisfying the decomposition of unit y I = N s X i =1 L i D i L T i are used. A simple c hoice for D i is a diagonal matrix with the diagonal elemen ts being the reciprocals of the n um ber of subdomains the degree of freedom is associated with. The Neumann-Neumann method used as a preconditioner of the problem is as follo ws: Giv en the residual r distribute it to subdomains r i = D T i L T i r solv e the local problems S i u i = r i and a v erage the results M 1 r = P N s i =1 L i D i u i [7, 43 ]. The dra wbac k of the Neumann-Neumann preconditioner is that it lac ks a mec hanism of distributing the error globally This has been resolv ed b y adding a coarse space to the problem. The resulting method is called the balancing domain decomposition. Let Z i be the matrix with linearly independen t columns that generate the k ernel of K i Im Z i = Ker K i If K i is regular, Z i is a v oid matrix. The balancing preconditioner is as follo ws: 1. Balance the original residual r b y solving the auxiliary problems Z T i D T i L T i ( r S N s X j =1 L j D j Z j v j ) = 0 ; i = 1 ; : : :; N s 2. Distribute the balanced residual to the subdomains and nd a solution of 28 PAGE 38 the local problems S i u i = D T i L T i ( r S N s X j =1 L j D j Z j v j ) ; i = 1 ; : : :; N s 3. Balance b y solving the auxiliary problems Z T i D T i L T i ( r S N s X j =1 L j D j ( u j + Z j w j ) = 0 ; i = 1 ; : : :; N s 4. Av erage the result on the in terface M 1 r = N s X i =1 L i D i ( u i + Z i w i ) The BDD is pro v en to be independen t of n um ber of subdomains, and its condition n um ber is independen t of jumps of coecien ts bet w een subdomains with appropriate w eigh t matrices D i It can be written also as an abstract additiv e Sc h w arz algorithm [24 41]. 2.4.4 The Neumann-Neumann Domain Decomposition Algorithm for Plates and Shells Let us in troduce the Kirc hho-Lo v e model of plate bending follo wing [41]. W e consider a plate occup ying a domain in I R 2 whic h is clamped on the part of boundary @ n 00 and simply supported on @ n 0 n @ n 00 The Kirc hhoLo v e model plate model c haracterizes the v ertical displacemen t u 2 V of the plate as the solution of the v ariational problem a ( u; v ) = F ( v ) ; 8 v 2 V; (2.10) where the bilinear form on the righ t hand side is symmetric, con tin uous and coerciv e a ( u; v ) = Z n ( ( u )) : K : ( ( v )) ; 29 PAGE 39 the functional on the righ t hand side is giv en b y F ( v ) = Z n fv + Z @ n @ n 00 m g @v @ n + Z @ n @ n 0 gv: The space of kinematically admissible elds is V = f v 2 H 2 (n) ; v = 0 on @ n 0 ; v = @v @ n = 0 on @ n 00 g : In the denitions abo v e f is the densit y of v ertical forces, m g the densit y of rexion momen ts applied on the part of the boundary where the plate is free to rotate, and g is the densit y of v ertical boundary loading. The sym bol ":" denotes a tensor product, ( ) = 1 2 ( r + ( r ) T ) is the curv ature tensor, ( u ) = r u represen ts the in-plane notation of the plate and K is the plate rexural stiness. K is symmetric, elliptic and con tin uous in the sense that ( ( u )) : K : ( ( u )) ct 3 j ( ( u )) j 2 ( ( u )) : K : ( ( v )) Ct 3 j ( ( u )) jj ( ( v )) j ; where t is the plate thic kness. F or a simple case of an isotropic plate made of an homogeneous material with Y oung modulus E and P oisson coecien t it is giv en b y ( ( u )) : K : ( ( u )) = Et 3 12(1 2 ) ((1 ) r 2 u : r 2 v + u v ) : Since V is a subspace of H 2 (n), the appropriate nite elemen t spaces are C 1 con tin uous. Examples of suc h elemen ts include the discrete Kirc hho triangle (DKT) and HCT elemen ts. The Neumann-Neumann and BDD preconditioner do not perform w ell for the plate bending problem [41]. The condition n um ber estimate of BDD is 30 PAGE 40 based on the estimate [45 ] sup 8 < : P N s j =1 k L T j P N s i =1 L i D i u i k 2 S j P N s i =1 k u i k 2 S i : u i ? Ker S i ; S i u i ? Im Z i 9 = ; (2.11) A t crosspoin ts of subdomains, the v ectors in the estimate are constructed from con tributions sev eral subdomains dieren t from those that share an edge leading to the crosspoin t. This leads to a discon tin uit y that is inappropriate for a fourth order problem. The BDD method for plates [42] a v oids this problem arising at subdomain crosspoin ts b y enhancing the coarse space of the balancing domain decomposition algorithm. The coarse space is again Im Z i where Z i = [ x i 1 ; : : :; x in i ; y i 1 ; : : :y im i ] : f x i 1 ; : : :; x in i g is a basis for Ker S i and for eac h crosspoin t j = 1 ; : : :; m i of the subdomain n i y ij is the solution of the problem S i y ij = e ij with e ij the v ector corresponding to the unit normal load at the crosspoin t j With this c hoice, the normal displacemen t componen t of the v ectors u i coming out of the coarse space problem is zero since S i u i ? y ij implies that u i ? S i y ij = e ij Then, the suprem um (2.11) is tak en o v er functions with zero at endpoin ts of the edges, whic h mak es it possible to pro v e a poly-logarithmic bound. F or another approac h, imposing zeros at crosspoin t directly see [41]. 2.4.5 Lagrange Multipliers and P oincar e-Steklo v Operators F ollo wing [19], w e sho w ho w Lagrange m ultiplier approac h to enforcing solution con tin uit y is related to in terface form ulations using P oincar eSteklo v operators dened in Section 2.4.1. The discrete m ulti-domain case will be treated in the next c hapter. 31 PAGE 41 W e consider the problem from Section 2.4.1 and the notation in troduced there. W e reform ulate our problem as a constrained minimization problem: Find the solution ( u 1 ; u 2 ) 2 H 1 D (n 1 ) H 1 D (n 2 ) that minimizes 1 2 2 X k =1 a k ( u k ; u k ) Z n k fv k subject to the condition r 1 u 1 = r 2 u 2 Then, for eac h ( u 1 ; u 2 ; g ? ) 2 H 1 D (n 1 ) H 1 D (n 2 ) ( H 1 = 2 00 ()) 0 w e dene the Lagrangian ( u 1 ; u 2 ; ? ) = 1 2 2 X k =1 a k ( u k ; u k ) Z n k fu k Z ? ( u 1 u 2 ) : The Lagrange m ultiplier ? is used to to enforce the con tin uit y of the solution on the boundary The rest of the Lagrangian is the usual quadratic functional implied b y the w eak form of the problem. The critical poin ts ( u 1 ; u 2 ; ? ) of the saddle poin t problem of no w m ust satisfy the v ariational equalit y 2 X i =1 a k ( u k ; v k ) Z ( ? ( v 1 v 2 ) + ? ( u 1 u 2 ) d = 2 X i =1 Z n i fv k ; for all ( v 1 ; v 2 ; ? ) 2 H 1 D (n 1 ) H 1 D (n 2 ) ( H 1 = 2 00 ()) 0 This problem has a unique solution [19] whic h is the solution of the problem ( Q 1 + Q 2 ) ? = R 2 f 2 R 1 f 1 This is the equation (2.6). It sho ws that solving the Lagrange m ultiplier form ulation is equiv alen t to nding Neumann in terface data on (cf. [18]). The paper [19] uses the Lagrange form ulation to in troduce nite elemen t spaces of Lagrangians of small dimension per in terface for regular meshes. This can reduce the size of the problem substan tially but it is restricted to regular meshes. The space of Lagrangians can be c hosen as the restriction of the nite elemen ts space on the subdomains to the in terface. This is similar to the approac h tak en b y FETI as explained in the next c hapter. 32 PAGE 42 2.5 A Tw o Lev el Method Based on Smoothed Prolongation The t w o-lev el method described in [68 66 ] dev elops a simple abstract framew ork based on the concept of smoothed ten tativ e prolongator in troduced in [65]. The ten tativ e prolongator is deriv ed from a system of nono v erlapping subdomains. As opposed to the previously described methods, the union of all subdomains n i i = 1 ; : : :; N s does not co v er whole n. Instead, there is a la y er one elemen t wide bet w een eac h t w o subdomains. Our algorithm can be written as a v ariational t w o-lev el m ultigrid with a special c hoice of componen ts. W e will rst describe componen ts of the method and abstract assumptions that ensure coarse space size independen t con v ergence. Let us consider the system of linear algebraic equations Ku = f; where K is an n n symmetric positiv e denite stiness matrix arising from a discretization of a second order elliptic PDE, for example (2.1). In terpolation from the coarse space to the ne space is represen ted b y an operator MP composed from a ten tativ e prolongator P and prolongator smoother M Let us denote K M = M 2 K; M 0 = I K M ; where 2 (0 ; 2) is a giv en constan t and ( K M ) is an estimated upper bound for ( K M ), the spectral radius of K M The follo wing assumption species requiremen ts on M P and ( K M ). are needed for pro ving coarse space size independen t con v ergence. 33 PAGE 43 Assumption 7 There exist positiv e constan ts C 1 C 2 independen t of n and N s and a positiv e constan t C D ( N s ; n ) suc h that (1) F or ev ery u 2 I R n there exists v 2 I R N s suc h that k u Pv k C 1 C D ( N s ; n ) 1 = 2 ( K ) k u k K : (2.12) (2) The prolongator smoother M is symmetric, comm utes with K ( M ) < 1, and ( K M ) ( K M ) C 2 2 C 2 D ( N s ; n ) ( K ) : (2.13) W e also need to dene a pre-smoother S M and post-smoother S M 0 ; S M ( u; f ) = Mu + Nf; and analogously for M 0 The matrix N is c hosen so the smoother is consisten t with the system Ku = f The smoothers S M and S M 0 are relaxation operators used to smooth the appro ximate solution corrected b y an in terpolated coarse lev el error. The algorithm can no w be written as follo ws. Algorithm 8 Giv en an initial guess u repeat un til con v ergence: (1) u S M ( u; f ), (2) solv e ( P T MKMP ) w = P T M ( Ku f ) ; (3) u u MPw (4) u S M 0 ( u; f ). P ost process u S M ( u; f ). The ten tativ e prolongator P ma y be dened as P = D 1 = 2 ~ P; where D = diag( K ) and ~ P ij = 1 ; if the node corresponding to u i belongs to subdomain n j ; = 0 ; otherwise. 34 PAGE 44 The most practical c hoice of M is a polynomial in K The c hoice described belo w is a polynomial in K that attempts to minimize ( K M ) giv en an upper bound on ( K ) and a c hosen degree of the polynomial. F or certain degrees, this polynomial can be found explicitly b y the follo wing recursiv e algorithm: Algorithm 9 Let ^ be the estimate of ( K ) satisfying ( K ) ^ C ( K ), with a giv en constan t C Set ^ i = ^ 9 i A 0 = K and for i = 1 ; 2 ; : : : do: (1) A i = ( I 4 3 ^ 1 i 1 A i 1 ) 2 A i 1 (2) M i = Q i 1 j =0 ( I 4 3 ^ 1 j A j ) Notice that deg( M i ) 3 i 1 2 F or a 2 D problem, w e c hoose the prolongator smoother M = M k where deg( M k +1 ) qN 1 = 2 es deg( M k ) ; (2.14) where q 2 (0 ; 1] is a giv en parameter and N es the a v erage a v erage n um ber of degrees of freedom per subdomain, and M k are the polynomials constructed b y the algorithm abo v e. With the c hoices of the ten tativ e prolongator and prolongator smoother described abo v e, the con v ergence of our method can be sho wn to be independen t of the meshsize, N es inhomogeneities bet w een subdomains, and boundary conditions [68 ]. It th us o v ercomes one of the main disadv an tages of standard t w o-lev el m ultigrid methods, suering from the dependence of the rate of conv ergence on the size of the coarse space. On the other hand, its computational complexit y is lo w er in comparison with domain decomposition methods using direct local solv ers. This is true ev en for simple direct solv ers as in ternal solv ers in our method (whic h are preferred, as iterativ e solv ers tend to ha v e negativ e impact on the robustness of the method). 35 PAGE 45 Application of the method to solid problems is treated in [64] and similar results are obtained. The method is pro v ably robust with respect to jumps in coecien ts [67]. 36 PAGE 46 3. Finite Elemen t T earing and In terconnecting (Deriv ation and Extensions) 3.1 Original FETI The FETI (Finite Elemen t T earing and In terconnecting) method is a domain decomposition algorithm deriv ed from a h ybrid v ariational principle and designed for the iterativ e solution of systems of equations arising from the nite elemen t discretization of self-adjoin t elliptic partial dieren tial equations. It w as dev eloped in [28 27 33], and also discussed in detail in monograph [34]. In this method, a giv en spatial domain is \torn" in to non-overlapping subdomains where an incomplete solution of the primary eld is rst ev aluated using a direct solv er. Next, in tersubdomain eld con tin uit y is enforced via Lagrange m ultipliers applied at the subdomain in terfaces. This \gluing" phase generates a smaller size symmetric dual problem where the unkno wns are the Lagrange m ultipliers, and whic h is best solv ed b y a preconditioned conjugate gradien t (PCG) algorithm. This idea is related to the ctitious domain method where the Lagrange m ultipliers enforce boundary conditions as in Dinh et al. [17]. In con trast with other related domain decomposition methods using Lagrange m ultipliers as unkno wns [36 62], the FETI method distinguishes itself with the treatmen t of the n ull spaces of the subdomain stiness matrices (rigid 37 PAGE 47 body modes) associated with the so-called roating subdomains, i.e., subdomains without a sucien t n um ber of essen tial boundary conditions to prev en t the local stiness matrix from being singular. Resolving the rigid body modes leads to a small \coarse" problem that is solv ed in eac h PCG iteration. This is an added complication, but also a blessing. F arhat, Mandel, and Roux [32] ha v e sho wn n umerically and pro v ed for the FETI method without preconditioning, that the auxiliary problem pla ys the role of a coarse problem, namely it causes the condition n um ber to be bounded independen tly of the n um ber of subdomains. In [26], F arhat, Chen, and Mandel extended to time-dependen t problems, whic h lac k the naturally occurring coarse problem. The FETI method is in a sense dual to the Neumann-Neumann method with a coarse problem, dev eloped b y Mandel under the name Balancing Domain Decomposition [45 ] (BDD), whic h w e described in Chapter 2 and whic h is based on an earlier method of de Roec k and Le T allec [61 ]. A modied method w as analyzed b y Dryja and Widlund [25]. It should be noted that while the underlying ideas of FETI and BDD are in a w a y dual, FETI is not the BDD method applied to the dual problem. 3.1.1 Problem Setting and Assumptions Let us rst in troduce some notation used throughout this c hapter and the rest of this thesis. F or u; v 2 I R n the inner product h u; v i = u T v serv es also as dualit y pairing. The ` 2 norm is denoted k u k = h u; u i 1 = 2 F or a symmetric positiv e semidenite matrix A k u k A = h Au; u i 1 = 2 the induced seminorm. This is a norm if A is positiv e denite. The superscript + denotes pseudoin v erse, dened as follo ws. 38 PAGE 48 Denition 10 (Pseudoin v erse) Let A be a linear operator. The pseudoin v erse A + is an y linear operator suc h that if a 2 Im A then AA + a = a In general, a pseudoin v erse is not unique. Our algorithms will be in v arian t to a specic c hoice of the pseudoin v erse. If A is a symmetric operator on a nite dimensional space, A + can be c hosen to be also symmetric. Considering the spectral decomposition of A A = X v v T ; Av = v ; v T v = 1 ; (3.1) a pseudoin v erse A + can be c hosen as A + = X 6 =0 1 v v T : F or A positiv e semidenite, w e denote A = X > 0 v v T : In particular, with this notation w e ha v e, A = A 1 = 2 A 1 = 2 and A + = A 1 = 2 A 1 = 2 W e poin t out that Ker A = Ker A for an y real W e refer to Section 2.2 and the model problem (2.1) for explanation of some of the terminology used next. Let n be a domain in I R 2 or I R 3 decomposed in to N s non-o v erlapping subdomains n 1 ; n 2 ; :::; n N s Let u i be the v ector of degrees of freedom for subdomain n i corresponding to a conforming nite elemen t discretization of an elliptic problem (e.g. linear elasticit y plate bending) dened on n, suc h that eac h subdomain is a union of some of the elemen ts. Let u i K i and f i be the v ector of degrees of freedom, the local stiness matrix, and the load v ectors, 39 PAGE 49 respectiv ely associated with the subdomain n i W e will use the bloc k notation u = 2 6 6 6 6 6 6 6 6 6 6 4 u 1 u 2 . u N s 3 7 7 7 7 7 7 7 7 7 7 5 ; f = 2 6 6 6 6 6 6 6 6 6 6 4 f 1 f 2 . f N s 3 7 7 7 7 7 7 7 7 7 7 5 ; K = 2 6 6 6 6 6 6 6 6 6 6 4 K 1 0 : : : 0 0 K 2 : : : 0 : : : : : : : : : : : : 0 0 : : : K N s 3 7 7 7 7 7 7 7 7 7 7 5 : Depending on the imposed boundary conditions and the location of the subdomain, the local stiness matrix K i is positiv e denite or positiv e semidenite. A subdomain without sucien t essen tial boundary conditions to prev en t the subdomain stiness matrix K i from being singular is called a ro ating sub domain Let Z i be the matrix with linearly independen t columns that generate the k ernel of K i Im Z i = Ker K i If K i is regular, Z i is a v oid matrix. Denote Z = 2 6 6 6 6 6 6 6 6 6 6 4 Z 1 0 : : : 0 0 Z 2 : : : 0 : : : : : : : : : : : : 0 0 : : : Z N s 3 7 7 7 7 7 7 7 7 7 7 5 ; th us, Im Z = Ker K; Ker Z = f 0 g : A single mesh poin t x 2 n has sev eral degrees of freedom associated with it if it lies on the in terface (in tersection of the boundaries) of t w o or more subdomains, see Figure 3.1. Let B be a giv en matrix suc h that Bu = 0 expresses the condition that for eac h mesh node shared b y t w o or more subdomains the v alues of the degrees of freedom associated with that node coincide. W e denote b y W the space of all v ectors of degrees of freedom, and 40 PAGE 50 W W j i l l u u i m j n xFigure 3.1. Model domain: the mesh node x shared b y subdomains n i and n j and the degrees of freedom u i m and u j n associated with the mesh node x 41 PAGE 51 b y the space of the v ectors of v alues of the con tin uit y constrain t; th us, K : W W ; B : W : The problem to be solv ed is the minimization of the energy of the system subject to in tersubdomain con tin uit y conditions E ( u ) = 1 2 u T Ku f T u min subject to Bu = 0 ; u 2 W : (3.2) W e assume that the global structure is not roating, that is, the solution of (3.2) is unique. F rom (3.2), this is equiv alen t to the assumption that Ker K \ Ker B = f 0 g : (3.3) The follo wing notation will be needed to write the FETI algorithm concisely: G = BZ; F = BK + B T ; d = BK + f; (3.4) e = Z T f; P = I G ( G T G ) 1 G T : Note that P is an ` 2 orthogonal projection on (Im G ) ? = Ker G T Lemma 12 justies the expression for P b y sho wing that G T G is in v ertible. 3.1.2 Deriv ation of the Dual Problem The method, as originally deriv ed b y F arhat and Roux [27 ], in troduces Lagrange m ultipliers to enforce the con tin uit y of the solution. Solving the 42 PAGE 52 problem (3.2) leads to the system of equations Ku + B T = f Bu = 0 (3.5) A solution u of the rst equation in (3.5) exists if and only if f B T 2 Im K: (3.6) It m ust then ha v e the form u = K + ( f B T ) + Z: (3.7) where is to be determined. Substituting u from (3.7) in to the second equation of (3.5) yields BK + ( f B T ) + BZ = 0 : (3.8) No w (3.8) m ultiplied b y P together with (3.6) sho w that satises the system of equations P ( F d ) = 0 G T = e; (3.9) where w e ha v e used notation (3.4). The rst of the equations (3.9) is solv ed b y a projected preconditioned conjugate gradien t method using an initial appro ximation 0 that satises the second equation. The conjugate gradien t algorithm requires ev aluating only the actions of PF Most of the computational w ork in F = BK + B T is concentrated in the action of K + Since K + is subdomain bloc k diagonal, its action can be computed in parallel and in v olv es solving subdomain problems only Application of P requires solving a small coarse problem. F or a scalar problem, the size of this problem is less than the n um ber of subdomains N s F or a 43 PAGE 53 linear elasticit y problem, the size of the problem is the n um ber of rigid body modes, whic h again does not exceed a small m ultiple of N s A preconditioned projected CG for the equation PF = d using a symmetric preconditioner D can be written as follo ws. Algorithm 11 (FETI) Giv en an initial 0 compute the initial estimate 0 = G ( G T G ) 1 e + P 0 and the initial residual r 0 = P ( F 0 d ) : Repeat for k = 1 ; 2 ; : : : un til con v ergence: z k 1 = Dr k 1 y k 1 = Pz k 1 k = r T k 1 y k 1 p k = y k 1 + k k 1 p k 1 ( p 1 = y 0 ) k = k p T k PFp k k = k 1 + k p k r k = r k 1 k PFp k The c hoice of the preconditioner is discussed in Section 3.3. Lagrange m ultipliers in FETI can be seen as interfac e for c es and moments in the ph ysical system. F rom (3.7) and the denition of F in (3.4), the residual P ( F d ) = Bu has the in terpretation of jumps of the values of de gr e es of fr e e dom bet w een subdomains. The condition f B T ? Ker K 44 PAGE 54 means that the action of the loads and in tersubdomain forces and momen ts does not excite rigid body motions. 3.1.3 Saddle P oin t F orm ulation Let us explain in more detail ho w enforcing the con tin uit y using Lagrange m ultipliers yields the system (3.9). This approac h is presen ted b y F arhat, Mandel and Roux [32], and it will be later used to deriv e a generalized FETI method. F or the Lagrangian of the minimization problem (3.2), L ( u; ) = 1 2 u T Ku f T u + T Bu; u 2 W ; 2 ; w e solv e the dual problem: nd suc h that C ( ) = max 2 C ( ) max 2 inf u 2 W L ( u; ) : (3.10) By a direct computation, C ( ) = 8 > > < > > : 1 if f B T 6? Ker K; 1 2 h K + ( f B T ) ; f B T i otherwise : (3.11) The dual problem (3.10) is th us equiv alen t to maximizing C ( ) on the admissible set A = f 2 j C ( ) > 1g : The space of admissible incremen ts is f 1 2 j 1 2 A ; 2 2 Ag = f 2 j B T ? Ker K g = Ker G T : (3.12) A t the maxim um of C ( ), 2 A the deriv ativ e of C D C ( ; ), is zero in all directions in 2 Ker G T : D C ( ; ) = 0 8 2 Ker G T : 45 PAGE 55 By a straigh tforw ard computation, this becomes h BK + B T + BK + f; i = 0 ; 8 2 Ker G T : (3.13) In order to express (3.13) as a linear equation in the space Ker G T w e use the l 2 orthogonal projection dened b y (3.4). Then for 2 Ker G T h BK + B T + BK + f; i = h BK + B T + BK + f; P i = h P ( BK + B T + BK + f ) ; i Therefore, the dual problem (3.10) is equiv alen t to the linear equation in Ker G T for the unkno wn 2 Ker G T ; P ( BK + B T ( + 0 ) + BK + f ) = 0 ; (3.14) where 0 is an arbitrary starting feasible solution, that is, 0 2 A Denoting = 0 + w e obtain the system (3.9). 3.1.4 Algebraic Properties of the Original FETI Method Lemma 12 ( G T G ) 1 exists. Pr o of. Let G = BZ = 0. Then Z 2 Ker B and Z 2 Ker K b y the denition of Z It follo ws from (3.3) that Z = 0, hence = 0, since Z w as assumed to be of full rank. See also [34, Theorem 5.4]. Theorem 13 The solution of (3.9) is unique up to addition of a v ector from Ker B T An y solution of (3.9) yields the same solution u of the minimization problem (3.2), using (3.7) with = ( G T G ) 1 G T ( d F ). Pr o of. The relation bet w een and follo ws b y a direct computation. T o pro v e uniqueness of it is sucien t to sho w that Ker PF \ Ker G T = Ker F \ Ker G T = Ker B T : (3.15) 46 PAGE 56 First, B T = 0 implies F = 0 and G T = 0, so Ker B T Ker F \ Ker G T Ker PF \ Ker G T : Con v ersely assume G T = 0 and PF = 0. Then, since G T = Z T B T the rst equation implies B T ? Im Z = Ker K: Th us B T 2 Im K F rom the denition of P G T = 0 is equiv alen t to P = F rom PF = 0, w e obtain 0 = T PF = T F = ( B T ) T K + B T : Since K + is positiv e semidenite, this implies that B T 2 Ker K + = Ker K T ogether, B T 2 Ker K \ Im K = f 0 g : Non-uniqueness of the m ultipliers corresponds to redundan t in tersubdomain con tin uit y constrain ts, whic h occur naturally at crosspoin ts of more than t w o subdomains. Dene the space of the Lagrange m ultipliers as the factorspace ~ = = Ker B T : Since for an y 2 Ker B T P F and G T are the zero v ectors, the operators P F and G T induce operators on ~ whic h will be denoted b y the same sym bols. T o a v oid confusion, all n ull spaces will refer to the space not the factorspace ~ F or example, the n ull space of the induced operator G T on ~ will be denoted b y Ker G T = Ker B T It is easy to see that F is symmetric and positiv e semi-denite and hence so is PFP The next lemma sho ws that the operator PF restricted to Ker G T = Ker B T whic h coincides with the restriction of PFP is positiv e denite. 47 PAGE 57 Lemma 14 The operator PF is symmetric and positiv e denite on the factorspace Ker G T = Ker B T Pr o of. F or u; v 2 Ker G T w e ha v e from the denition of P h PFPu; v i = h PFu; v i = h Fu; v i ; whic h pro v es that PF is symmetric positiv e semidenite on Ker G T T o pro v e that PF is nonsingular on Ker G T = Ker B T let PFu = 0 and G T u = 0. Since F = BKB T and G T = Z T B T b y denition, it follo ws that B T u 2 Ker K and B T u ? Ker K hence B T u = 0. The original FETI method is therefore the method of preconditioned conjugate gradien ts in the factorspace ~ for the operator equation, equiv alen t to (3.9), PF = Pd; 2 0 + Ker G T ; (3.16) where 0 is an initial appro ximation to the conjugate gradien t method c hosen so that G T 0 = e and all searc h directions are in the space Ker G T 3.2 Generalized FETI There ha v e been sev eral extensions to the original method. Extension to time-dependen t problems w as done in [26]. Generalization to plate-bending problems is discussed in [51, 53] and, in detail, in this thesis. F or the original method for plate bending problems, the condition n um ber w as observ ed to gro w fast with the n um ber of elemen ts per subdomain [32]. This is caused b y the fact that plate bending is a fourth order problem, while the FETI domain decomposition method \tears" the appro ximate solution at subdomain crosspoin ts, whic h is suitable only for second order problems. 48 PAGE 58 The limitation of method can be cured to extend the FETI methodology to obtain a non-o v erlapping domain decomposition method for plate bending problems. This new method has the properties one usually looks for in iterativ e substructuring methods: the condition n um ber can be bounded independen tly of the n um ber of subdomains, and it gro ws only poly-logarithmicall y with the n um ber of elemen ts per subdomain. The computational cost per iteration is proportional to the solution of a boundary v alue problem in eac h subdomain, plus the solution of a sparse c o arse pr oblem with only few v ariables per subdomain. Suc h methods are commonly referred to as sc alable and quasi-optimal though, of course, for v ery large n um ber of subdomains, the solution of the coarse problem w ould dominate. The k ey idea of our method is to enforce the con tin uit y of the appro ximate solution at the subdomain crosspoin ts throughout the iterations b y adding the corresponding Lagrange m ultipliers to the coarse problem. A similar idea w as emplo y ed in the Balancing Domain Decomposition (BDD) method for plates [42], where appro ximate con tin uit y at crosspoin ts is enforced b y adding new basis functions to the original coarse space [46 49 ] in order to k eep the energy of the appro ximate solution minimal with respect to displacemen ts that are solutions for poin t loads at the subdomain crosspoin ts (cf. Section 2.4.4). The distinguishing features of both the presen t method and the method from [42] is that they are non-overlapping and w ork for standard nite elemen ts used in ev eryda y engineering practice. F or other domain decomposition methods for the biharmonic equation and plate bending see, for example, [12, 71]. Extensions to shells, implementation issues, and further computational can be found in [29, 31 ]. 49 PAGE 59 3.2.1 Dualit y Deriv ation of the Generalized FETI In this section, w e presen t a deriv ation and form ulation [52 ] of the generalized FETI, based on the concept of coarse optimalit y of a dual objectiv e function. Preserving the condition G T k = e throughout the iterations of the Algorithm 11 can be in terpreted as enforcing that ev ery k be optimal with respect to all possible incremen ts of the form G : G T k = e ( ) C ( k ) C ( k G ) ; 8 ; (3.17) where C is dened b y (3.11). The k ey to the generalization of the FETI method to plate bending problems is to mak e all k optimal in more directions. Let C be some giv en matrix with the same n um ber of ro ws as G Eac h column of C will giv e rise to an additional v ariable in the coarse problem. W e shall satisfy in eac h iteration the c o arse optimality pr op erty C ( ) C ( G C ) ; 8 ; : (3.18) with = k T o satisfy this propert y consider an auxiliary problem: F or a giv en ~ nd and so that C ( ) max ; = ~ G C (3.19) Since w e only need solutions satisfying C ( ) > 1 w e consider the maximization problem (3.19) along with the constrain t G T = e (3.20) 50 PAGE 60 In troducing new Lagrange m ultipliers for (3.20), w e get that , and solv e the saddle poin t problem inf ; sup ~ L ( ~ G C; ) = sup inf ; ~ L ( ~ G C; ) where ~ L ( ; ) = 1 2 T F + T d + T ( G T e ) F rom the optimalit y conditions @ ~ L ( ; ) @ = 0 ; @ ~ L ( ; ) @ = 0 ; @ ~ L ( ; ) @ = 0 ; = ~ G C; w e obtain that (3.19) is equiv alen t to the bloc k linear system M 2 6 6 6 6 6 6 4 3 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 4 G T F C T F G T 3 7 7 7 7 7 7 5 ~ 2 6 6 6 6 6 6 4 G T d C T d e 3 7 7 7 7 7 7 5 (3.21) where M = 2 6 6 6 6 6 6 4 G T FG G T FC G T G C T FG C T FC C T G G T G G T C 0 3 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 4 G T 0 C T 0 0 G T 3 7 7 7 7 7 7 5 2 6 6 4 F I I 0 3 7 7 5 2 6 6 4 G C 0 0 0 G 3 7 7 5 : (3.22) The solution of (3.19) is unique up to the addition of a v ector in Ker F \ Ker G T and w e write it as = ~ G C 0 M + 0 B B B B B B @ 2 6 6 6 6 6 6 4 G T F C T F G T 3 7 7 7 7 7 7 5 ~ 2 6 6 6 6 6 6 4 G T d C T d e 3 7 7 7 7 7 7 5 1 C C C C C C A : (3.23) 51 PAGE 61 Note that since (3.19) has a solution for an y ~ so does (3.21), hence Im 2 6 6 6 6 6 6 4 G T F C T F G T 3 7 7 7 7 7 7 5 Im M: (3.24) F urther, if is coarse optimal, then + is also coarse optimal if and only if G C 0 M + 2 6 6 6 6 6 6 4 G T F C T F G T 3 7 7 7 7 7 7 5 2 Ker F \ Ker G T : (3.25) The generalized FETI algorithm is th us obtained from Algorithm 11 b y prescribing the initial appro ximation 0 b y 0 = 0 G C 0 M + 0 B B B B B B @ 2 6 6 6 6 6 6 4 G T F C T F G T 3 7 7 7 7 7 7 5 0 2 6 6 6 6 6 6 4 G T d C T d e 3 7 7 7 7 7 7 5 1 C C C C C C A and b y replacing the step y k 1 = Pz k 1 b y y k 1 = z k 1 [ G; C; 0] M + 2 6 6 6 6 6 6 4 G T F C T F G T 3 7 7 7 7 7 7 5 z k 1 ; cf. Algorithm 20 in the next section. 3.2.2 Deriv ation of Generalized FETI Algorithm without Dualit y The generalized FETI method can be also obtained b y forcing the iterates to satisfy also a w eigh ted residual condition [53]. That is, w e require throughout the iterations that C T P ( F d ) = 0 ; G T = e; (3.26) 52 PAGE 62 where C is some giv en matrix. Searc h directions that preserv e (3.26) form the space V 0 = f 2 j G T = 0 ; C T PF = 0 g : (3.27) The corresponding space of residuals is V = f u 2 Im B j G T u = 0 ; C T u = 0 g : (3.28) Denote the associated subspaces of the factorspace ~ as ~ V 0 = f u + Ker B T j u 2 V 0 g = V 0 = Ker B T ; ~ V = f u + Ker B T j u 2 V g : Note that V and ~ V are isomorphic: since V Im B and Im B \ Ker B T = f 0 g eac h class of ~ V con tains exactly one elemen t of V W e will need sev eral properties of the spaces ~ V and ~ V 0 Lemma 15 ~ V = PF ~ V 0 : Pr o of. F rom the denition, clearly PFV 0 V hence PF ~ V 0 ~ V T o sho w that ~ V PF ~ V 0 let u + Ker B T 2 ~ V Since PF is a bijection on Ker G T = Ker B T b y Lemma 14, there is ~ 2 Ker G T = Ker B T suc h that PF ~ = u + Ker B T It follo ws that ~ 2 V 0 Lemma 16 The space ~ V 0 is the dual of ~ V with the dualit y pairing h ; i Pr o of. An y 2 ~ V 0 denes a linear functional 0 on ~ V b y 0 ( v ) = h ; v i Let 0 be an arbitrary linear functional on ~ V F rom Lemma 15, it follo ws that 0 PF is a linear functional on ~ V 0 and, from Riesz represen tation theorem, there is a unique 2 ~ V 0 suc h that 0 ( PFv ) = h ; v i for all v 2 ~ V 0 F rom Lemma 14, PF is a bijection on ~ V and it follo ws that the mapping bet w een a linear functional 0 and its represen tation 2 ~ V 0 is an isomorphism. 53 PAGE 63 Lemma 17 ~ = ~ V 0 ~ V ? : Pr o of. F rom Lemma 14, PF is symmetric positiv e denite on Ker G T = Ker B T Hence, ( PF ) 1 is also symmetric, positiv e denite on Ker G T = Ker B T and, using Lemma 15, it follo ws that ~ = ~ V 0 ( ~ V 0 ) ? ( PF ) 1 = ~ V 0 (( PF ) ~ V ) ? ( PF ) 1 = ~ V 0 ~ V ? ; (3.29) whic h w as to be pro v ed. W e no w dene a projection operator Q b y Q : ~ ~ ; Q 2 = Q; Im Q = ~ V 0 ; Ker Q = ~ V ? ; (3.30) and compute the matrix represen tation of Q This represen tation denes also an operator on whic h will be denoted b y the same sym bol Q Lemma 18 The projector Q is giv en b y the form ula Q = I G C 0 M + 2 6 6 6 6 6 6 4 G T F C T F G T 3 7 7 7 7 7 7 5 ; where M is dened b y (3.22). Pr o of. Let ; 2 Q ( + Ker B T ) = + Ker B T Since 2 V 0 from the denition of P there exists suc h that PF = F + G F rom the denition of V 0 G T = 0 and C T PF = 0, hence G T F + G T G = 0 ; C T F + C T G = 0 ; G T = 0 : (3.31) 54 PAGE 64 F rom the denition of Q w e ha v e ? V so for an y u 2 Im B G T u = 0 and C T u = 0 implies h ; u i = 0. Consequen tly there exist and suc h that for all u 2 Im B h ; u i = h G; u i + h C; u i ; whic h implies that = + G + C + r; r 2 (Im B ) ? = Ker B T (3.32) Since 2 V 0 substituting (3.32) in to (3.31) giv es that , satisfy the linear system G T F ( + G + C ) + G T G = 0 C T F ( + G + C ) + C T G = 0 G T ( + G + C ) = 0 : (3.33) On the other hand, it is easy to see from (3.31) and the denitions of Q and V that if , satisfy (3.33), then = + G + C 2 Q + Ker B T Similarly as in the case of the original FETI, w e need to nd an initial appro ximation 0 satisfying certain conditions, in this case (3.26). Lemma 19 F or an y 0 2 the system of equations G T F ( 0 + G + C ) + G T G = G T d C T F ( 0 + G + C ) + C T G = C T d G T ( 0 + G + C ) = e (3.34) has a solution ; ; and 0 = 0 + G + C (3.35) satises C T P ( F 0 d ) = 0 ; G T 0 = e: (3.36) 55 PAGE 65 Pr o of. As in the preceding proof, 0 satises (3.36) if and only if there is a suc h that G T F 0 + G T G = G T d C T F 0 + C T G = C T d G T 0 = e (3.37) Substituting 0 from (3.35) in to (3.37) yields the system (3.34), whic h can be written as M 2 6 6 6 6 6 6 4 3 7 7 7 7 7 7 5 = X T 2 6 6 4 d G ( G T G ) 1 e 3 7 7 5 ; X = 2 6 6 4 G C 0 0 0 G 3 7 7 5 : F rom the factorization (3.22), it follo ws that Ker M = Ker X Using symmetry of M w e ha v e Im M = (Ker M ) ? = (Ker X ) ? = Im X T ; hence, (3.34) has a solution. The generalized FETI method is the conjugate gradien ts method for the operator PF : V 0 V preconditioned b y QDQ T : V V 0 where the D : is a giv en operator symmetric on V Since Q T is also a projection and Im Q T = (Ker Q ) ? = V the application of Q T on V can be omitted, and one obtains the follo wing algorithm similar to the algorithm of the original FETI (Algorithm 11). Algorithm 20 (Generalized FETI) Giv en an initial 0 compute the initial 0 from (3.34) and (3.35), and compute the initial residual r 0 = P ( F 0 d ) : Repeat for k = 1 ; 2 ; : : : un til con v ergence: z k 1 = Dr k 1 56 PAGE 66 y k 1 = Qz k 1 k = r T k 1 y k 1 p k = y k 1 + k k 1 p k 1 ( p 1 = y 0 ) k = k p T k PFp k k = k 1 + k p k r k = r k 1 k PFp k The c hoice of the preconditioner is discussed in Section 3.3. 3.2.3 Method Selection for Plate Problems and Other Generalizations W e c hoose the columns of matrix C whic h appears in the description of the generalized FETI method abo v e, as v ectors with a one at the position of a Lagrange m ultiplier that enforces the con tin uit y of the transv ersal displacemen t at a crosspoin t, and zeros ev erywhere else. By a crosspoin t w e understand an in terface node adjacen t to at least three subdomains or to t w o subdomains and the complemen t of n. F or a more precise form ulation, see Section 4.5.1. A similar idea can be exploited for shell problems. The con tin uit y of the normal displacemen t needs to be enforced. T o a v oid nding normals and the added complexit y of enforcing con tin uit y of the normal displacemen t, one ma y enforce con tin uit y of all displacemen t degrees of freedom. This, ho w ev er, increases the size of the coarse space. Other possibilities and computational issues are discussed in [30]. F arhat et al. in [26 ] studies the case of time-dependen t problems where the subdomain stiness matrices K i are perturbed b y the addition of 57 PAGE 67 a m ultiple of the subdomain mass matrix, th us making the new local matrix positiv e denite. Consequen tly all matrices Z i are v oid and the natural coarse problem is lost in time-dependen t applications. The methodology dev eloped in [26] for rein troducing a coarsening operator in the FETI algorithm for dynamics problems is a special case of the presen t generalization where C is tak en to be the matrix G b efor e the p erturb ation that is, C = [ B i ~ Z i ] where the columns of ~ Z i are the basis of the k ernel of the local stiness matrix of the subdomain n i The reason wh y the preconditioner w orks for the dynamics problems is quite dieren t from that for the plate bending. 3.3 The Diric hlet Preconditioner Let us decompose the space of all degrees of freedom W in to the space of in ternal degrees of freedom and the degrees of freedom on subdomain in terfaces, W = W W: In the corresponding bloc k notation, B = [ 0 B ] ; B : W ; since B has nonzero en tries for in terface degrees of freedom only Also, Z = 2 6 6 4 Z Z 3 7 7 5 ; and w e ha v e G = BZ = B Z; Ker B T = Ker B T : 58 PAGE 68 Let S be the Sc h ur complemen t of K obtained b y elimination of degrees of freedom in ternal to all subdomains. Then F = BK + B T = BS + B T (3.38) and Ker S = Im Z Ev aluation of the matrix-v ector product S + u reduces to the solution of independen t Neumann pr oblems on all subdomains, cf. Section 2.4.1. Inspired b y (3.38), w e c hoose D = BS B T giving the preconditioner QD = Q BS B T : (3.39) This preconditioner is called the Dirichlet pr e c onditioner since ev aluating the matrix-v ector product Sr is equiv alen t to solving independen t Dirichlet pr oblems on all subdomains, cf. Section 2.4.1. 3.4 In terface F orm ulation of P ark et al. In [59, 58 ], K.C. P ark has dev eloped a similar substructuring method. W e will sho w that this method can essen tially be written as FETI with a special c hoice of the in terface con tin uit y operator B The method augmen ts the Lagrangian b y in troducing another v ariable u g the global v ector of degrees of freedom, on the in terfaces bet w een subdomains. This v ariable is redundan t and is later eliminated, yielding an analog of FETI. W e will sho w that the method can be written as FETI and discuss the resulting c hoice of B and its implications. Let L be the subdomain assem bly matrix, that is a zero-one matrix mapping the local subdomain degrees of freedom to global ones. Using this matrix, the con tin uit y of the solution is enforced through the constrain t u = Lu g : 59 PAGE 69 The Lagrangian can be written as L ( u; u g ; ) = 1 2 u T Ku f T u + T ( u Lu g ) : Decomposing u = Z + in to a rigid body mode componen t and its complemen t (pure deformation componen t), w e obtain L ( ; ; u g ; ) = 1 2 T K f T f T Z + T ( Z + Lu g ) : The stationarit y condition leads to the system of equations 2 6 6 6 6 6 6 6 6 6 6 4 K I 0 0 I 0 Z L 0 Z T 0 0 0 L 0 0 3 7 7 7 7 7 7 7 7 7 7 5 2 6 6 6 6 6 6 6 6 6 6 4 u g 3 7 7 7 7 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 6 6 6 6 4 f 0 Zf 0 3 7 7 7 7 7 7 7 7 7 7 5 Eliminating the deformation componen t, = K + ( f ), the system becomes 2 6 6 6 6 6 6 4 K + Z L Z T 0 0 L T 0 0 2 6 6 6 6 6 6 4 2 6 6 6 6 6 6 4 u g 3 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 4 K + f Zf 0 3 7 7 7 7 7 7 5 The rst equation of the system is K + + Z Lu g = K + f: (3.40) W e rst eliminate the term Z b y m ultiplying (3.40) b y the orthogonal projection P Z on to (Im Z ) ? = Ker Z T P Z = I Z ( Z T Z ) 1 Z T and obtain the equation P Z K + P Z Lu g = P Z K + f; (3.41) 60 PAGE 70 equiv alen t to (3.40) holding for some Next, using the orthogonal projection P l on to (Im P Z L ) ? = Ker( P Z L ) T P l = I P Z L ( L T P Z L ) 1 L T P Z ; in a similar w a y P l P Z K + = P l P Z K + f: This is the equation (15) in [58 ]. W e ma y rev erse the order of elimination of the terms, and w e obtain from (3.40) P z P L K + = P z P L K + f; (3.42) where the orthogonal projections P L and P z are giv en b y P L = I L ( L T L ) 1 L T P z = I P L Z ( Z T P L Z ) 1 Z T P L : Since Im P z Im P L = Ker L T w e ha v e = P L for 2 Im( P z ). Th us, P L K + = P L K + P L : No w c hoosing B = P L ; (3.43) w e ma y write P z = I P L Z ( Z T P L Z ) 1 Z T P L = I G ( G T G ) 1 G T = P; where w e ha v e used the notation (3.4). This sho ws that the method with the equation (3.42) can be written as FETI with the special c hoice of B (3.43). 61 PAGE 71 Let us illustrate this c hoice. If a mesh node is on the in terface of t w o subdomains, then, omitting all other degrees of freedom and Lagrange m ultipliers, the corresponding bloc k of the assem bly matrix L is 2 6 6 4 1 1 3 7 7 5 and the corresponding bloc k of the projection P L is 1 2 2 6 6 4 1 1 1 1 3 7 7 5 : Similarly for a three node and four node in terfaces, w e get 1 3 2 6 6 6 6 6 6 4 2 1 1 1 2 1 1 1 2 3 7 7 7 7 7 7 5 ; and 1 4 2 6 6 6 6 6 6 6 6 6 6 4 3 1 1 1 1 3 1 1 1 1 3 1 1 1 1 3 3 7 7 7 7 7 7 7 7 7 7 5 ; respectiv ely The c hoice of B = P L has some in teresting implications on the analysis, whic h will be discussed in Section 4.4. 62 PAGE 72 4. Con v ergence Analysis In this c hapter, w e pro v e that for a second order problem the condition n um ber of the preconditioned FETI method is bounded independen tly of the n um ber of subdomains and poly-logarithmical ly in terms of subdomain size, similarly as it is in the case for other optimal non-o v erlapping domain decomposition methods [11 23 25 49 48 ]. W e pro v e that the c hoice of the in terface con tin uit y operator as giv en b y P ark et al. assures that the condition n um ber is bounded b y C (1+log H h ) 2 Finally w e also sho w that the generalized FETI con v erges poly-logarithmicall y for the plate bending problem problem. Analysis of domain decomposition methods t ypically demonstrates spectral equiv alence of the quadratic form that denes the problem in a v ariational setting and the quadratic form that denes the preconditioner, often b y w a y of P .L. Lions lemma [6, 23, 24 44 ]. Since the preconditioner in the FETI method is quite complicated and is not dened in terms of a quadratic form, w e proceed dieren tly W e nd a bound on the norm of the product of the system operator and the preconditioner, so as to bound the maximal eigenv alue, and a bound on the in v erse, to bound the minimal eigen v alue. Related analyses w ere previously done for methods without crosspoin ts bet w een the subdomains, or done formally in functional spaces, cf., for example, Glo winski and Wheeler [37 ] and Bakh v alo v and Kn y azev [4]. W e presen t a complete analysis in terms of upper and lo w er bound on the preconditioned operator for decompositions with crosspoin ts in 2D and edges and crosspoin ts in 3D for 63 PAGE 73 second order elliptic problems. W e sho w that the condition n um ber is bounded b y C (1 + log H h ) r where r = 2 ; 3, h is the c haracteristic mesh size, and H the diameter of the subdomains. 4.1 Preliminaries In this section, w e presen t some results that will be used in analysis of FETI methods. The results are mostly concerned with estimates for nite elemen t functions used to discretize our model problems. W e will assume that a domain n is divided in to a set of nono v erlapping subdomains n i ; i = 1 ; : : :; N s n = n 1 [ : : : [ n N s The subdomains are assumed to be shape regular of diameter H according to the denition belo w. W e will form ulate all results for one of the subdomains n i and assume that constan ts in the estimates do not depend on the index of the subdomain. Denition 21 A subdomain n i I R d is said to be shape-regular of diameter O ( H ) if it can be generated from a reference domain (square or cube) ^ n of unit diameter b y a mapping F i suc h that n i = F i ( ^ n). The mapping is assumed to satisfy k @ F i k c s H; k @ F 1 i k c s H 1 (4.1) where @ F i is the Jacobian of the mapping, k : k is the Euclidean I R d matrix norm, and c s is a positiv e constan t. W e will be using the Sobolev spaces W k;p F or p = 2, the spaces are Hilbert spaces and w e denote them H k W k; 2 Denitions of Sobolev spaces can be found, for example, in [1]. F or k 6 = 1 = 2 ; 1, w e use k k k;p; n and j j k;p; n to denote the standard Sobolev norms and seminorms of functions in W k;p F ollo wing [22 ], w e dene the scaled Sobolev norm for a scalar function 64 PAGE 74 u 2 H 1 (n i ) k u k 2 1 ; 2 ; n i = j u j 2 1 ; 2 ; n i + 1 H 2 k u k 2 0 ; 2 ; n i ; where the Sobolev seminorm is dened b y j u j 2 1 ; 2 ; n i = Z n i kr u ( x ) k 2 dx: W e dene the scaled Sobolev norm for a scalar function u 2 H 1 = 2 ( @ n i ) k u k 2 1 = 2 ; 2 ;@ n i = j u j 2 1 = 2 ; 2 ;@ n i + 1 H k u k 2 0 ; 2 ;@ n i ; where the Sobolev seminorm is dened b y j u j 2 1 = 2 ; 2 ;@ n i = Z @ n i Z @ n i j u ( x ) u ( y ) j 2 k x y k d dxdy: Here k k is the Euclidean norm in I R d W e note that the space H 1 = 2 ( @ n i ) is the space of traces of functions in the space H 1 (n i ). If u = ( u 1 ; u 2 ) is a v ector function, then w e dene, for example, k u k 2 1 ; 2 ; n i = k u 1 k 2 1 ; 2 ; n i + k u 2 k 2 1 ; 2 ; n i : Other norms and seminorms of v ector functions are dened analogously W e presen t t w o v arian ts of the P oincar e-F riedric hs inequalit y W e pro v e the rst for the sak e of completeness and to demonstrate the tec hnique that is used to pro v e its v arian ts for nite elemen t functions. W e note that P k ; k = 0 ; 1 ; : : : denotes the space of polynomials of degree at most k Lemma 22 Let G be a con tin uous linear functional on H 1 = 2 ( @ n i ) suc h that for all u 2 P 0 G ( u ) = 0 implies u = 0. Then there exists a constan t c independen t of H suc h that for all u 2 H 1 = 2 ( @ n i ) ; u 2 Ker G k u k 1 2 ; 2 ;@ n i c j u j 1 2 ; 2 ;@ n i : 65 PAGE 75 Pr o of. W e pro v e the theorem for a reference domain of diameter 1 b y con tradiction. The result in the scaled norms then follo ws. If the inequalit y is not true, there exists a sequence f u n g Ker G suc h that k u n k 1 2 ; 2 ;@ ^ n = 1 and j u n j 1 2 ; 2 ;@ ^ n 0 : Due to the compact im bedding of H 1 = 2 ( @ ^ n) in L 2 ( @ ^ n), there exists a subsequence of f u n g that con v erges in L 2 ( @ ^ n). Since j u n j 1 2 ; 2 ;@ ^ n 0, the subsequence f u n k g is Cauc h y in H 1 = 2 ( @ ^ n). Therefore it con v erges in H 1 = 2 ( @ ^ n) to some u 2 H 1 = 2 ( @ ^ n). The con tin uit y of the norm and the seminorm on H 1 = 2 ( @ ^ n) implies k u k 1 2 ; 2 ;@ ^ n = 1 and j u j 1 2 ; 2 ;@ ^ n = 0 : (4.2) Therefore, b y the second equation in (4.2), u = k almost ev erywhere, where k is some real n um ber. Let us assume without loss of generalit y that u = k ev erywhere; that is u 2 P 0 Since G is con tin uous G ( u n k ) G ( u ) = 0. This yields, b y the assumption, u = 0 whic h con tradicts the rst equation in (4.2). The proof of this P oincar e-F riedric hs inequalit y in scaled norms is similar to the proof of the previous lemma. Lemma 23 Let G be a con tin uous linear functional on H 1 (n i ) suc h that for all u 2 P 0 G ( u ) = 0 implies u = 0. Then there exists a constan t c independen t of H suc h that for all u 2 H 1 (n i ) ; u 2 Ker G k u k 1 ; 2 ; n i c j u j 1 ; 2 ; n i : 66 PAGE 76 4.1.1 Estimates for P1/Q1 Elemen ts Let V P 1 h (n i ) be a conforming nite elemen t space of P1 or Q1 elemen ts [14] satisfying the usual regularit y and in v erse properties and possibly some essen tial boundary conditions. That is, for example in I R 2 w e assume that n i = [ K2T h;i K ; where eac h elemen t K of the triangulation T h;i is a triangle or a rectangle. F urthermore, for all K 2 T h;i c 1 d ( K ) h c 2 ( K ) ; (4.3) where d ( K ) is the diameter of K and ( K ) the diameter of the circle inscribed in K Then h is called the c haracteristic mesh size. A v ertex of an elemen t K 2 T h;i will be referred to as a mesh poin t, a nodal poin t, or just a v ertex. If essen tial boundary conditions are prescribed on @ n i w e assume that () c c ( @ n i ), where ( ) denotes the measure. W e note that functions in V P 1 h (n i ) are con tin uous. As in the previous c hapter, the corresponding space of v ectors of degrees of freedom is denoted W i Let I P 1 : W i V P 1 h (n i ) denote the linear one-to-one transformation that maps a v ector of degrees of freedom to the corresponding nite elemen t function. This transformation is often called nite elemen t in terpolation. W e decompose a v ector of degrees of freedom U 2 W i in to in ternal and boundary degrees of freedom assuming the boundary degrees of freedom are listed last. That is, in bloc k notation, U = [ U T ; U T ] T where U is the v ector of boundary degrees of freedom. W e denote W i the space of boundary degrees 67 PAGE 77 of freedom and dene the matrix T i so that U = T i U Then, I P 1 U denes a function on @ n i whic h depends on U only By abuse of notation, w e will write I P 1 U = I P 1 [0 T ; U T ] T on @ n i W e summarize some w ell kno wn results and inequalities in a form suitable for our purposes. The next lemma summarizes the fact that the H 1 = 2 norm of a zero extension of a "piece" of a function can be bounded b y the norm of the function. Lemma 24 Let @ n i be a v ertex, edge, or face (if d = 3) of subdomain n i A face is understood not to con tain adjacen t edges, and an edge does not con tains its endpoin ts. F or for all z 2 V P 1 h ( @ n i ), dene w 2 V P 1 h ( @ n i ) b y w ( x ) = z ( x ) on all nodes of triangulation x 2 w ( x ) = 0 on all other nodes of @ n i Then k w k 2 1 2 ; 2 ;@ n i C (1 + log H h ) ( j z j 2 1 2 ; 2 ;@ n i + 1 H k z k 2 0 ; 2 ;@ n i ) ; where = 1 if d = 2 and is a v ertex, or d = 3 and is an edge or a v ertex; = 2 if d = 2 and is an edge, or d = 3 and is a face. Pr o of. The inequalit y for d = 2 w as pro v ed in [49, 48]. The case when d = 3 follo ws from Lemmas 4.1 and 4.2 in [11] if is an edge or a v ertex, and Lemma 4.3 in [11] if is a face. Cf. also [23]. The follo wing lemma can be pro v ed b y using Lemmas 4.1 and 4.2 in [11 ] and estimates the H 1 = 2 seminorm of a spik e function. Lemma 25 There exists a constan t c independen t of h and H suc h that for all u 2 V P 1 h (n i ), u ( x ) = 0 for all mesh poin ts x 2 @ n i x 6 = x 0 j u j 2 1 2 ; 2 ;@ n i c k u k 2 0 ; 1 ;@ n i 68 PAGE 78 The follo wing lemma sho ws the equiv alence of the discrete seminorm dened b y the Sc h ur complemen t and the H 1 = 2 seminorm. In the next section w e will apply this result to the stiness matrix K i whic h is wh y w e the same sym bol here. This result is standard [10 69 ] and it is pro v ed here for the sak e of completeness. Lemma 26 Let K i be a symmetric matrix that satises c j I P 1 U j 2 1 ; 2 ; n i h K i U; U i C j I P 1 U j 2 1 ; 2 ; n i 8 U 2 W i : Let U = T i U and the matrix K i be decomposed so that K i U = 2 6 6 4 K ~ K T ~ K K 3 7 7 5 2 6 6 4 U U 3 7 7 5 : Let S i be the Sc h ur complemen t of K in K i Then, there exist constan ts c 0 and C 0 independen t of h and H suc h that c 0 j I P 1 U j 2 1 2 ; 2 ;@ n i h S i U; U i C 0 j I P 1 U j 2 1 2 ; 2 ;@ n i 8 U 2 W i : Pr o of. F rom the denition of the Sc h ur complemen t (Denition 6) and the propert y (2.8), it follo ws that h S i U; U i = h K U; U i h K 1 ~ K T U; ~ K T U i = inf U =[ U T ; U T ] T h K i U; U i Hence, from the assumption b y in v ariance of the j j 1 ; 2 ; n i seminorm to adding a constan t and b y the discrete extension theorem [22] for the scaled norms, w e obtain h S i U; U i C inf U =[ U T ; U T ] T j I P 1 U + k j 1 ; 2 ; n i Cc e k I P 1 U + k k 1 2 ; 2 ;@ n i ; 69 PAGE 79 where k 2 I R. Choosing k = R @ n i I P 1 Udx R @ n i dx yields the equiv alence of the scaled norm to the seminorm b y the P oincar e-F riedric hs inequalit y (Lemma 22). On the other hand, b y the trace theorem, in v ariance of the seminorms to adding a constan t, and the P oincar e-F riedric hs inequalit y (Lemma 23), w e ha v e 1 c t c p c j I P 1 U j 2 1 2 ; 2 ;@ n i c inf U =[ U T ; U T ] T j I P 1 U j 1 ; 2 ; n i h S i U; U i ; where c t and c p are constan ts arising from the trace and P oincar e-F riedric hs inequalit y respectiv ely The follo wing theorem sho ws a v arian t of the P oincar e-F riedric hs inequalit y W e note that the dimension of W i implicitly depends on h and for a uniform mesh, the lemma follo ws from Lemma 22. Lemma 27 Let S i be the Sc h ur complemen t from Lemma 26. There exists a constan t C independen t of H and h suc h that for all U 2 W i U ? Ker S i k I P 1 U k 1 2 ; 2 ;@ n i C j I P 1 U j 1 2 ; 2 ;@ n i : Pr o of. If there is an essen tial boundary condition imposed on @ n i then Ker S i = f 0 g and the statemen t follo ws from the P oincar e-F riedric hs inequalit y Lemma 22. Otherwise, Ker S i is spanned b y the v ector of ones. Then the condition U ? Ker S i implies P N j =1 U j = 0 for U = [ U 1 ; : : :; U N ] T F ollo wing along the lines of the proof of Lemma 22, w e obtain a sequence of v ectors U n 2 W i I R N ( n ) suc h that ( I P 1 U n ) F i u in L 2 ( @ ^ n) ; and k u k 1 2 ; 2 ;@ ^ n = 1 and j u j 1 2 ; 2 ;@ ^ n = 0 (4.4) 70 PAGE 80 Therefore, u = k almost ev erywhere, where k 2 I R. Th us, ( I P 1 U n ) F i k in L 2 ( @ ^ n). Then, from the quasi-uniformit y of the triangulation, w e obtain k ( I P 1 U n ) F i k k 2 0 ; 2 ;@ ^ n = c 1 N ( n ) k U n k k 2 c N ( n ) 2 0 @ N ( n ) X j =1 ( U n j k ) 1 A 2 Since U n ? Ker S i it holds that P N ( n ) j =1 U n j = 0. Th us, ck 2 k ( I P 1 U n ) F i k k 2 0 ; 2 ;@ ^ n 0 : This is possible only if k = 0 whic h is a con tradiction with (4.4). 4.1.2 Estimates for HCT Elemen ts Let V HCT h (n i ) be the nite elemen t space of HCT elemen ts satisfying the usual regularit y and in v erse properties and possibly some essen tial boundary conditions. That is, w e assume that n i = [ K 2T h;i K ; where eac h elemen t K of the triangulation T h;i is a triangle. F urthermore, condition (4.3) is satised and h denotes the c haracteristic mesh size. W e note that functions in V HCT h (n i ) are C 1 con tin uous and the space V HCT h (n i ) can be written as V HCT h (n i ) = f v 2 C 1 (n i ) : 8K 2 T h;i ; v j K j 2 P 3 ( K j ) for all subtriangles K j of K ; @v @n j a k a j 2 P 1 ( a k a j ) for all sides ( a k ; a j ) of the triangle Kg On eac h triangle, a function v in V HCT h (n i ) is determined b y the v alues v ( a i ) and the v alues of its deriv ativ es @v @x j ( a i ) at the v ertices of the triangle. The 71 PAGE 81 in ternal v ertex of subtriangles is not arbitrary but it is determined to guaran tee the unisolv ence of the elemen t. This is as m uc h detail as w e will need in our considerations. W e refer to [14] for more details about the HCT elemen t. The corresponding space of v ectors of degrees of freedom is denoted W i Let I P 1 : W i V HCT h (n i ) denote the linear one-to-one transformation that maps a v ector of degrees of freedom to the corresponding nite elemen t function. As in the previous section, w e decompose a v ector of degrees of freedom U in to in ternal and boundary degrees of freedom assuming the boundary degrees of freedom are listed last. That is, in bloc k notation, U = [ U T ; U T ] T where U is the v ector of boundary degrees of freedom. W e denote W i the space of boundary degrees of freedom and dene the matrix T i so that U = T i U Then, I HCT U denes a function on @ n i whic h depends on U only By abuse of notation, w e will write I HCT U = I HCT [0 T ; U T ] T on @ n i W e summarize here some w ell kno wn results and inequalities in a form suitable for our purposes. The follo wing lemma estimates the norm of a \spik e" function and is pro v ed in [42 ]. Lemma 28 Let x be a v ertex of a subdomain n i F or u 2 V HCT h (n i ) suc h that u ( x ) = 0, dene z 2 V HCT h (n i ) b y z ( x ) = u ( x ) = 0 ; r z ( x ) = r u ( x ), and z ( y ) = 0 ; r z ( y ) = 0 at all other nodes y of @ n i Then kr z k 2 1 2 ; 2 ;@ n i C 1 + log H h kr u k 2 1 2 ; 2 ;@ n i + 1 H kr u k 2 0 ; 2 ;@ n i : The follo wing estimate of the trace norm of the extension b y zero is pro v ed as in [10 Lemma 3.5]. Lemma 29 Let u 2 V HCT h (n i ). Then there exists a constan t C suc h that if 72 PAGE 82 supp u \ @ n i is con tained in a segmen t of @ n i of length then j u j 2 1 2 ; 2 ;@ n i j u j 2 1 2 ; 2 ; + C 1 + log h k u k 2 0 ; 1 ; : The follo wing estimate is a modication of the previous lemma. Lemma 30 Let u 2 V HCT h (n i ) and supp u \ @ n i be con tained in a segmen t of @ n i and x 1 and y 1 be its endpoin ts. Let x 2 and y 2 be the nodal poin ts next to x 1 and y 1 respectiv ely in the segmen t. Let ( x 2 ; y 2 ) denote the segmen t bet w een x 2 and y 2 Let the length of this segmen t be h Assume that the function u satises the condition k u k 0 ; 1 ; ( x 1 ;y 1 ) c k u k 0 ; 1 ; ( x 2 ;y 2 ) Then j u j 2 1 2 ; 2 ;@ n i j u j 2 1 2 ; 2 ; ( x 2 ;y 2 ) + C 1 + log h k u k 2 0 ; 1 ; ( x 2 ;y 2 ) : Pr o of. By the previous lemma, the inequalit y holds with ( x 2 ; y 2 ) replaced b y the segmen t ( x 1 ; y 1 ). By the denition of the H 1 = 2 seminorm, w e ha v e j u j 2 1 2 ; 2 ; ( x 1 ;y 1 ) = j u j 2 1 2 ; 2 ; ( x 2 ;y 2 ) + Z ( x 1 ;x 2) Z ( x 1 ;y 2) j u ( x ) u ( y ) j 2 k x y k 2 dxdy + Z ( y 1 ;y 2) Z ( x 1 ;y 2) j u ( x ) u ( y ) j 2 k x y k 2 dxdy: Using the fact that j u ( x ) u ( y ) j min f c 1 c h k u k 0 ; 1 ; ( x 2 ;y 2 ) k x y k ; 2 c k u k 0 ; 1 ; ( x 2 ;y 2 ) g ; the last t w o terms are bounded b y c 2 k u k 0 ; 1 ; ( x 2 ;y 2 ) W e will also need a straigh tforw ard extension of the discrete Sobolev inequalit y of Dryja [22 ] to piecewise polynomial functions of order p > 1 [57]. Lemma 31 F or ev ery p 1 there exists a constan t C = C ( p ) suc h that for ev ery u con tin uous on @ n i of length suc h that u 2 P p on the side of ev ery elemen t, k u k 2 0 ; 1 ; C 1 + log h j u j 2 1 2 ; 2 ; + 1 k u k 2 0 ; 2 ; : 73 PAGE 83 The follo wing lemma establishes a useful inequalit y bet w een discrete and Sobolev norms. Lemma 32 Let u 2 V HCT h (n i ) and u ( x 0 ) = 0 for some x 0 2 @ n i Let U = I 1 HCT u and U = T i U Then, for all suc h U h k U k 2 c (1 + H 2 ) kr u k 2 0 ; 2 ;@ n i ; where c is independen t of h and H Pr o of. Let E be an elemen t edge on @ n i and x 1 and x 2 be its endpoin ts. Eac h componen t of r u is a polynomial of degree at most t w o on E Since all norms in a nite dimensional space are equiv alen t, it holds that h ( kr u ( x 1 ) k 2 + kr u ( x 2 ) k 2 ) c 1 min f 2 P 2 P 2 ;f ( x 1 )= r u ( x 1 ) ;f ( x 2 )= r u ( x 2 ) k f k 2 0 ; 2 ;E c 1 kr u k 2 0 ; 2 ;E Summing o v er all edges of the boundary w e obtain h k U k 2 c 1 2 kr u k 2 0 ; 2 ;@ n i + h k U 0 k 2 ; where U 0 denotes the displacemen t degrees of freedom of U W e sho w that h k U 0 k 2 < c 2 H 2 kr u k 2 0 ; 2 ;@ n i Since u ( x 0 ) = 0, w e can write for an y mesh poin t x on the boundary u ( x ) = Z ( x 0 ;x ) r u ( y ) ( y ) dy; where ( x 0 ; x ) @ n i is the part of the boundary bet w een x 0 and x and is the v ector tangen tial to the boundary Squaring, using the Cauc h y inequalit y and considering that the length of ( x 0 ; x ) is bounded b y c 2 H w e obtain u ( x ) 2 c 2 H kr u k 2 0 ; 2 ;@ n i : 74 PAGE 84 The statemen t of the lemma then follo ws b y summing o v er all mesh poin ts on the boundary The follo wing lemma sho ws the equiv alence of the discrete seminorm dened b y the Sc h ur complemen t and the H 1 = 2 seminorm of the gradien t of the corresponding nite elemen t function and is an analog of Lemma 26 of the previous section. Lemma 33 Let K i be a matrix that satises c j I HCT U j 2 2 ; 2 ; n i h K i U; U i C j I HCT U j 2 2 ; 2 ; n i 8 U 2 W i : Let U = T i U and let the matrix K i be decomposed so that K i U = 2 6 6 4 K ~ K T ~ K K 3 7 7 5 2 6 6 4 U U 3 7 7 5 : Let S i be the Sc h ur complemen t of K in K i Then, there exist constan ts c 0 and C 0 independen t of h and H suc h that c 0 jr I HCT U j 2 1 2 ; 2 ;@ n i h S i U; U i C 0 jr I HCT U j 2 1 2 ; 2 ;@ n i 8 U 2 W i : Pr o of. F rom the denition of the Sc h ur complemen t it follo ws that h S i U; U i = h K U; U i h K 1 ~ K T U; ~ K T U i = inf U =[ U T ; U T ] T h K i U; U i Hence, b y the discrete extension theorem [41], w e obtain h S i U; U i C inf U =[ U T ; U T ] T j I P 1 U j 2 ; 2 ; n i Cc e jr I P 1 U j 1 2 ; 2 ;@ n i : On the other hand, the trace theorem and in v ariance of the j j 1 2 ; 2 ; n i seminorm to adding a linear function, and the P oincar e-F riedric hs theorem, jr I P 1 U j 2 1 2 ; 2 ;@ n i c t c p inf U =[ U T ; U T ] T j I P 1 U j 2 ; 2 ; n i : 75 PAGE 85 The follo wing theorem sho ws another v arian t of the P oincar e-F riedric hs inequalit y W e again emphasize that the dimension of W i implicitly depends on h and is an analog of Lemma 27 of the previous section. Lemma 34 Let S i be the Sc h ur complemen t from Lemma 33. There exists a constan t C independen t of H and h suc h that for all U 2 W i U ? Ker S i kr I HCT U k 1 2 ; 2 ;@ n i C jr I HCT U j 1 2 ; 2 ;@ n i : Pr o of. Let us pro v e the case when Ker S is non trivial. In the other case, the statemen t follo ws from the standard P oincar e-F riedric hs inequalit y F ollo wing along the lines of the proof of Lemma 22, w e obtain a sequence of v ectors U n 2 W i I R N ( n ) suc h that ( r I HCT U n ) F i u in L 2 ( @ ^ n) L 2 ( @ ^ n) ; and k u k 1 2 ; 2 ;@ ^ n = 1 and j u j 1 2 ; 2 ;@ ^ n = 0 (4.5) Therefore, u = ( k 1 ; k 2 ) 2 I R I R almost ev erywhere. Th us, ( r I HCT U n ) F i ( k 1 ; k 2 ) in L 2 ( @ ^ n) L 2 ( @ ^ n). W e decompose the v ector U n in to t w o v ectors corresponding to the partial deriv ativ e with respect to x U n 1 and the partial deriv ativ e with respect to y U n 2 Then, b y the proof of Lemma 32, w e obtain k ( r I HCT U n ) F i ( k 1 ; k 2 ) k 2 0 ; 2 ;@ ^ n c 1 N ( n ) ( k U n 1 k 1 k 2 + k U n 2 k 2 k 2 ) c N ( n ) 2 0 @ N ( n ) X i =1 ( U n 1 ;i k 1 ) 1 A 2 + c N ( n ) 2 0 @ N ( n ) X i =1 ( U n 2 ;i k 2 ) 1 A 2 : Since U n ? Ker S i it holds that P N ( n ) i =1 U n 1 ;i = 0 and P N ( n ) i =1 U n 2 ;i = 0. Th us, c ( k 2 1 + k 2 2 ) k ( r I HCT U n ) F i ( k 1 ; k 2 ) k 2 0 ; 2 ;@ ^ n 0 : 76 PAGE 86 This is possible only if ( k 1 ; k 2 ) = (0 ; 0) whic h is a con tradiction with (4.5). 4.2 Abstract Analysis F ramew ork This section dev elops a framew ork for estimating the condition n umber of our method. W e estimate the minim um and maxim um eigen v alues of the matrix of the FETI form ulation and the preconditioner. The follo wing simple lemma will be the basis of our estimates. It will allo w to reduce estimates of norms to estimates of boundedness and coercivit y The proof follo ws a standard argumen t and it is presen ted for completeness only Lemma 35 Let X be a Banac h space, X 0 the dual of X and A : X X 0 a linear operator that is on to and satises the conditions h y; Ax i = h x; Ay i ; 8 x; y 2 X (4.6) c A k x k 2 X h x; Ax i C A k x k 2 X ; 8 x 2 X (4.7) with constan ts C A ; c A > 0. Then k A k X X 0 C A ; k A 1 k X 0 X 1 c A : Pr o of. F rom (4.6), k A k X X 0 = sup x 2 X k Ax k X 0 k x k X = sup x; ~ x 2 X h x; A ~ x i k x k X k ~ x k X = sup x 2 X h x; Ax i k x k 2 X C A : F rom (4.7) and the fact that A is on to, w e obtain 1 k A 1 k X 0 X = inf x 2 X k Ax k X 0 k x k X = inf x 2 X sup ~ x 2 X h x; A ~ x i k x k X k ~ x k X inf x 2 X h x; Ax i k x k 2 X c A ; concluding the proof. 77 PAGE 87 Assume that the preconditioner of a FETI method can be written in the form RDR T where R is a projection, and RDR T is an isomorphism from H on to H 0 Let H and its dual H 0 are subspaces of Ker G T F urthermore, let H = Im R T The original FETI satises this with R = P and H = Ker G T = Ker B T The generalized FETI satises the assumption with R = Q and H = ~ V F or the purpose of analysis, w e equip the space H with the norm k v k H = k B T v k S = h S B T v; B T v i 1 = 2 : (4.8) Since B T v ? Ker S for v 2 H (4.8) indeed denes a norm rather than only seminorm. The dual space H 0 is equipped with the dual norm k k H 0 = sup v 2H h ; v i k v k H : (4.9) Since H ~ it immediately follo ws k k H 0 = sup v 2H h ; v i k B T v k S = sup w 2 W; Bw 2H h ; Bw i k B T Bw k S : (4.10) The norm on H w as c hosen so that the preconditioner RD is trivially coerciv e and bounded. Lemma 36 F or all v 2 H h v; RDv i = k v k 2 H : Pr o of. Let v 2 H = Im R T Since R T is a projection, w e ha v e b y denition of the preconditioner D h v; RDv i = h R T v; BS B T v i = h v; BS B T v i = h B T v; S B T v i = k v k 2 H ; whic h w as to be sho wn. 78 PAGE 88 Coercivit y and boundedness of the system operator PF on H 0 will be estimated using the follo wing lemma. Lemma 37 F or all 2 H 0 h ; F i = sup w 2 W;w ? Ker S h ; Bw i 2 k w k 2 S = sup w 2 W h ; Bw i 2 k w k 2 S : Pr o of. Let 2 H 0 Then h ; F i = h S + B T ; B T i = h S 1 = 2 B T ; S 1 = 2 B T i = k S 1 = 2 B T k 2 = sup x 2 W h S 1 = 2 B T ; x i 2 k x k 2 = sup x 2 W; x = x 1 + x 2 x 1 2 Ker S x 2 ? Ker S h B T ; S 1 = 2 x i 2 k x 1 + x 2 k 2 = sup x 2 2 W; x 2 ? Ker S h B T ; S 1 = 2 x 2 i 2 k x 2 k 2 since S 1 = 2 x 1 = 0 and k x k 2 = k x 1 k 2 + k x 2 k 2 No w write an y w 2 W as w = w 1 + w 2 ; w 1 2 Ker S; w 2 = S 1 = 2 x 2 ? Ker S: 2 H 0 implies that h B T ; w 1 i = 0, hence h B T ; w 2 i = h B T ; w i = h ; Bw i It follo ws that h ; F i = sup w 2 2 W; w 2 ? Ker S h B T ; w 2 i 2 h w 2 ; Sw 2 i = sup w 2 W h ; Bw i 2 k w k 2 S ; whic h w as to be pro v ed. It is w ell kno wn [39] that after k iterations of the preconditioned conjugate gradien t method, the energy norm or the error jjj e jjj = h PFe; e i 1 = 2 is reduced b y a factor of at least 2(( p 1) = ( p +1)) k ; where is the condition n um ber. The condition n um ber is giv en in our case b y = ( RDPF ) = max ( RDPF ) min ( RDPF ) ; (4.11) where max and min are the maxim um and minim um eigen v alue, respectiv ely W e are no w ready to pro v e an abstract bound on 79 PAGE 89 Theorem 38 Assume there exist constan ts C 1 C 2 suc h that (i) for an y 2 H 0 and w 2 W suc h that Bw 2 H there is ~ w 2 W suc h that h ; B ~ w i = h ; Bw i ; and k ~ w k 2 S C 1 k B T Bw k 2 S ; (ii) for an y 2 H 0 and w 2 W w ? Ker S there is ~ w 2 W suc h that B ~ w 2 H h ; B ~ w i = h ; Bw i ; and k B T B ~ w k 2 S C 2 k w k 2 S : Then ( RDPF ) C 1 C 2 Pr o of. Lemma 35 applied to the operator RD together with Lemma 36 giv e k QD k 2 H!H 0 1 ; k ( RD ) 1 k 2 H 0 !H 1 : (4.12) F rom assumption (i), w e ha v e sup w 2 W h ; Bw i 2 k w k 2 S 1 C 1 sup w 2 W; Bw 2H h ; Bw i 2 k B T Bw k 2 S (4.13) while assumption (ii) giv es the con v erse bound sup w 2 W;w 2 Ker S h ; Bw i 2 k w k 2 S C 2 sup w 2 W; Bw 2H h ; Bw i 2 k B T Bw k 2 S : (4.14) Using (4.10) and Lemma 37, w e see that the inequalities (4.13) and (4.14) imply the inequalit y 1 C 1 k k 2 H 0 h ; PF i C 2 k k 2 H 0 ; 8 2 H 0 : Applying Lemma 35 to the operator PF w e obtain k PF k 2 H 0 !H C 2 ; k ( PF ) 1 k 2 H!H 0 C 1 : (4.15) F rom (4.12) and (4.15), w e ha v e k QDPF k H 0 !H 0 k QD k H!H 0 k PF k H 0 !H C 2 80 PAGE 90 and k ( QDPF ) 1 k H 0 !H 0 k ( QD ) 1 k H 0 !H k ( PF ) 1 k H!H 0 C 1 : The result follo ws. 4.3 Analysis of the Original FETI In this section, w e use the abstract analysis framew ork from the previous section to pro v e a bound on the condition n um ber of FETI for a second elliptic problem, utilizing results from Section 4.1. 4.3.1 Assumptions W e need more specic assumptions in order to be able to pro v e a bound on the condition n um ber Let us recall the model problem (2.1). Assume w e are solving the boundary v alue problem A u = g in n ; u = 0 on @ n ; where A is a second order elliptic operator A v = d X i;j =1 @ @x i ( x ) @v ( x ) @x j ; with ( x ) a measurable function suc h that 0 < 0 ( x ) 1 a.e. in n. The domain n is assumed to be divided in to non-o v erlapping subdomains n i i = 1 ; :::; N s that are shape regular of diameter H (Denition 21). Assume that V h (n) is a conforming P1 or Q1 nite elemen t space on a triangulation of n, whic h satises the standard regularit y and in v erse assumptions (cf. Section 4.1). In particular, w e recall that the degrees of freedom are v alues at nodes of the triangulation h denotes the c haracteristic elemen t size. Eac h subdomain n i is assumed to be a union of some of the elemen ts, and all functions in V h (n) are zero on @ n. 81 PAGE 91 Let us dene the restriction operator R i : W W i b y the equation R i w = w i ; where w = [ w T 1 w T 2 : : :w T N s ] T w i 2 W i ; i = 1 ; : : :; N s W e will assume that the in terface con tin uit y operator B is dened as follo ws. Let w r ( x ) and w s ( x ) denote the pair of degrees of freedom corresponding to the mesh node x 2 @ n r \ @ n s and let ( B w ) rs ( x ) be the en try of v ector B w that corresponds to the mesh poin t x and subdomains n r and n s W e require eac h suc h ( B w ) rs ( x ) to ha v e form ( Bw ) rs ( x ) = rs ( w r ( x ) w s ( x )) ; (4.16) where rs is either 1 or -1. Note that rs does not depend on x that is, coecien ts rs are uniform along edges (and faces, in 3D) bet w een t w o subdomains. F or eac h node x that belongs to the in terface of t w o and more subdomains, x 2 T @ n i ; i = s 1 ; s 2 ; : : :; s n v ector Bw con tains n 1 en tries, ( Bw ) s k s k +1 ( x ) ; k = 1 ; : : :; s n 1 F or an example of the denition of the v alues of Bw for ( s 1 ; s 2 ; s 3 ) = (1 ; 3 ; 2) in 2D around a crosspoin t, see Fig. 4.1. W e poin t out that B c hosen in this w a y has full ro w rank, that is, this denition implies that there are no redundan t constrain ts in enforcing the con tin uit y of the solution at the nodes where more than t w o subdomains meet. Only the impro v ed estimate in statemen t 3 of Lemma 39 will require the specic denition of this section. If redundan t constrain ts are allo w ed, the estimates, with the exceptions of the impro v ed estimates, still hold whic h can be sho wn b y minor modications of the proofs. See Section 4.4. An additional connectivit y assumption on the decomposition is needed. Let N rs ; r ; s = 1 ; : : :; N s ; r 6 = s be the n um ber of in terface conditions 82 PAGE 92 bet w een subdomains n r and n s N rr the n um ber of degrees of freedom of subdomain n r and N = max N s i =1 N ii W e will assume that there exists constan ts c and n 0 independen t of h and H suc h that for all r ; s = 1 ; : : :; N s ; r 6 = s for whic h N rs > 0, there exists a sequence of indices f r i g r 0 = r r k = s k < n 0 suc h that N r i 1 r i cN 8 i = 1 ; : : :; k: Throughout the next section, c; C; c 1 ; c 2 ; c 3 ; c 4 ; c 5 and c 6 denote positiv e constan ts independen t of H and h 4.3.2 Discrete Norm Bounds The stiness matrices K i ; i = 1 ; : : :; N s satisfy the assumption of Lemma 26. Therefore, using Lemma 26 for eac h subdomain and summing o v er all subdomains, w e obtain c 1 N s X i =1 j I P 1 R i w j 2 1 2 ; 2 ;@ n i k w k 2 S c 2 N s X i =1 j I P 1 R i w j 2 1 2 ; 2 ;@ n i 8 w 2 W ; (4.17) where the positiv e constan ts c 1 ; c 2 are independen t of the c haracteristic mesh size h and the subdomain diameter H The constan ts c 1 ; c 2 ma y depend on the regularit y of the shape of subdomains. The follo wing lemmas v erify the assumptions of Theorem 38. Lemma 39 F or all 2 Ker G T and all w 2 W suc h that B T Bw ? Ker S there exists ~ w 2 W suc h that h ; B ~ w i = h ; Bw i and k ~ w k 2 S C (1 + log H=h ) k B T Bw k 2 S : where = 1, and = 0 in the follo wing special cases: (1) 1 2 B B T = I whic h means that there are no nodes shared b y more than t w o subdomains. 83 PAGE 93 Figure 4.1. Denition of B A A A A A A A v v v v v v v v v r r r r r n 1 n 2 n 3 @ n 3 @ n 2 @ n 1 w 3 ( x ) w 3 ( x ) w 3 ( x r ) w 2 ( x r ) w 1 ( x ) w 1 ( x ) w 2 ( x ) w 2 ( x ) w 1 ( x ) w 1 ( x ) w 3 ( x ) w 1 ( x ) w 3 ( x ) w 2 ( x ) w 3 ( x ) w 2 ( x r ) w 3 ( x r ) w 1 ( x ) w 2 ( x ) 84 PAGE 94 (2) d = 2, and the matrix B has the follo wing propert y: If w 2 Im B T x is a crosspoin t (node shared b y more than t w o subdomains), I P 1 R i w ( x ) = I P 1 R i w ( y ) for all i suc h that x 2 @ n i and all nodes y that are adjacen t to x on @ n i then I P 1 R i w ( x ) = 0 for all i suc h that x 2 @ n i (3) d = 2, B is dened as in Section 4.3.1 and all nodes in the triangulation belong to either one, t w o, or an odd n um ber of subdomains. Pr o of. Let us rst pro v e that, in the general case, w e obtain 1. Let w 2 W and dene ~ w = B T ( B B T ) 1 Bw The triangle inequalit y sho ws that k ~ w k S r r r r 1 2 B T Bw r r r r S + r r r r r 1 2 B T I 1 2 B B T 1 Bw r r r r r S : (4.18) Denote z = 1 2 B T ( I ( 1 2 B B T ) 1 ) Bw F rom the denition of B in (4.16), z is zero at all nodes that belong to at most t w o subdomains. The remaining nodes lie on crosspoin ts or edges (in the 3D case) of subdomains. F rom the denition of B at ev ery suc h node x z is a linear com bination of the en tries of B T Bw that correspond to the same node x and the coecien ts of the linear com binations are bounded only in terms of the n um ber of subdomains the node belongs to. Using Lemma 24 for the crosspoin ts of subdomains and the equiv alence (4.17), w e obtain for the 2D case that k z k 2 S c N s X i =1 j I P 1 R i z j 2 1 2 ; 2 ;@ n i C (1 + log( H=h )) N s X i =1 k I P 1 R i B T Bw k 2 1 2 ; 2 ;@ n i : (4.19) The P oincar e inequalit y Lemma 27, yields N s X i =1 k I P 1 R i B T Bw k 2 1 2 ; 2 ;@ n i k B T Bw k 2 S In the 3D case, the argumen t for subdomain crosspoin ts is the same. In addition, w e note that the coecien ts of the linear com bination do not c hange along a subdomain edge, so it remains to apply Lemma 24 on ev ery edge. 85 PAGE 95 Let us no w turn to the special cases that giv e = 0. If 1 2 B B T = I w e c hoose ~ w = 1 2 B T Bw = w whic h pro v es the special case 1. No w w e pro v e special case 2. F rom the denition of the H 1 = 2 norm, the equiv alence of norms (4.17), and the fact that I P 1 R i B T Bw is a piecewise linear function, it follo ws that k B T Bw k 2 S c N s X i =1 j I P 1 R i B T Bw j 2 1 2 ; 2 ;@ n i (4.20) c X x crosspoin t ;x 2 @ n i y adjacen t to x;y 2 @ n i I P 1 R i B T Bw ( x ) I P 1 R i B T Bw ( y ) 2 : F or an y crosspoin t x it follo ws from the assumption that for ev ery w 2 Im B T X i;@ n i 3 x y adjacen t to x; y 2 @ n i ( I P 1 R i w ( x ) I P 1 R i w ( y )) 2 = 0 ) X i;@ n i 3 x ( I P 1 R i w ( x )) 2 = 0 : Consequen tly b y compactness, and since there are only nitely man y dieren t n um bers of subdomains sharing a crosspoin t, for all w 2 Im B T X i;@ n i 3 x ( I P 1 R i w ( x )) 2 C X i;@ n i 3 x y adjacen t to x; y 2 @ n i ( I P 1 R i w ( x ) I P 1 R i w ( y )) 2 : By summation o v er all crosspoin ts x using Lemma 25 and (4.20), w e get k z k 2 S C k B T u k 2 S ; whic h concludes the proof of this case. In order to pro v e case 3, w e v erify the assumptions of case 2. W e form ulate the proof only for a crosspoin t shared b y three subdomains (Fig. 4.1). The proof is similar for a dieren t odd n um ber of subdomains. Let w 2 Im B T Since I P 1 R 1 w ( x ) I P 1 R 1 w ( x ) = 0, and I P 1 R 1 w ( x ) I P 1 R 1 w ( x ) = 0, w e ha v e I P 1 R 1 w ( x ) = I P 1 R 1 w ( x ). Similarly w e obtain I P 1 R 2 w ( x ) = 86 PAGE 96 I P 1 R 2 w ( x r ), and I P 1 R 3 w ( x ) = I P 1 R 3 w ( x r ). Moreo v er w 2 Im B T implies I P 1 R 1 w ( x ) = I P 1 R 2 w ( x ), I P 1 R 2 w ( x r ) = I P 1 R 3 w ( x r ), and I P 1 R 3 w ( x ) = I P 1 R 1 w ( x ), whic h can be satised only if I P 1 R 1 w ( x ) = I P 1 R 1 w ( x ) = ::: = 0. Remark 40 In general, the exponen t = 1 in Lemma 39 cannot be impro v ed. T o see that, let us consider the conguration with v alues of u and Bu in the neigh borhood of a crosspoin t as in Fig. 4.2, whic h violate the assumptions of special case 2. Extending the v alues of u in Fig. 4.2 to deca y as log r ( t=H ), r < 1 = 2, where t is the distance from the crosspoin t, w e obtain a v ector u 2 suc h that k B T u k S C; k I P 1 u k 1 2 ; 2 ;@ n 1 \ @ n 2 j log h=H j r : If u = Bw then on @ n 1 \ @ n 2 u = w 2 w 1 whic h giv es j I P 1 u j 1 2 ; 2 ;@ n 1 \ @ n 2 j I P 1 w 1 j 1 2 ; 2 ;@ n 1 \ @ n 2 + j I P 1 w 2 j 1 2 ; 2 ;@ n 1 \ @ n 2 j I P 1 w 1 j 1 2 ; 2 ;@ n 1 + j I P 1 w 2 j 1 2 ; 2 ;@ n 2 k w k S ; and therefore k w k S C ( r ) j log h=H j r for all r < 1 = 2. Lemma 41 Let 2 Ker G T Then for all w 2 W w ? Ker S there is a ~ w 2 W suc h that B T B ~ w ? Ker S h ; Bw i = h ; B ~ w i ; and k B T B ~ w k 2 S C (1 + log H=h ) 2 k w k 2 S Pr o of. Let w 2 W be arbitrary and put B ~ w = PBw Then h ; Bw i = h ; B ~ w i : (4.21) Denote Bz = Bw B ~ w z 2 Ker S Then, since P is an orthogonal projection, k Bz k k Bw k : 87 PAGE 97 Figure 4.2. Coun ter-example t t t t t t t t t t t t q q q q q q q n 1 n 2 n 3 n 4 1 1 1 1 1 0 1 1 = 2 1 = 2 1 = 2 1 = 2 1 = 2 1 = 2 1 = 2 1 = 2 1 = 2 1 = 2 1 = 2 1 = 2 88 PAGE 98 No w, from the denition of B and the P oincar e inequalit y Lemma 27, w e obtain h k Bw k 2 2 h k w k 2 CH k w k 2 S : Also, since z 2 Ker S it is constan t on eac h @ n i using the connectivit y assumption, w e ha v e the follo wing b y Lemma 24: k B T Bz k 2 S C h H k Bz k 2 (1 + log H=h ) 2 : T ogether this yields k B T Bz k 2 S C (1 + log H=h ) 2 k w k 2 S : By the denition of B R i B T Bw on @ n i [ @ n j is a linear com bination (with bounded coecien ts) of (a bounded n um ber of) w k from all @ n k adjacen t to @ n i [ @ n j F rom Lemma 24 and the P oincar e-F riedric hs inequalit y Lemma 27, w e obtain k B T Bw k S C (1 + log( H=h )) k w k S ; 8 w 2 W: Finally summarizing, k B T B ~ w k S k B T Bw k S + k B T Bz k S C (1 + log H=h ) k w k S : F rom this and (4.21), the statemen t of the lemma follo ws. 4.3.3 Condition Num ber Estimate W e are no w ready to pro v e the nal result. It follo ws from the abstract estimate in Lemma 35 with the assumptions v eried b y Lemma 39 and Lemma 41. 89 PAGE 99 Theorem 42 Under the assumptions of section 4.3.1, the condition n um ber of the FETI method with the Diric hlet preconditioner satises = max ( PDPF ) min ( PDPF ) C (1 + log H h ) r with r = 3, and r = 2 in the special cases listed in Lemma 39. 4.4 Analysis of the Method b y P ark W e consider the same assumptions as for the original FETI for second order problems except for the c hoice of B W e consider B giv en b y projection matrix P L (3.43). The essen tial propert y of this matrix B used here is that B = B T B = B T since B is an orthogonal projection. W e presen t t w o lemmas that v erify assumptions of Theorem 38. Lemma 43 F or all 2 Ker G T = Ker B T and all w 2 W suc h that Bw 2 Ker G T = Ker B T there exists ~ w 2 W suc h that h ; B ~ w i = h ; Bw i ; k ~ w k S = k B T Bw k S Pr o of. Since B T = B and B 2 = B w e ma y c hoose ~ w = w Lemma 44 Let 2 Ker G T = Ker B T Then for all w 2 W w ? Ker S there is a ~ w 2 W suc h that B T B ~ w ? Ker S and h ; Bw i = h ; B ~ w i ; k B T B ~ w k 2 S C (1 + log H=h ) 2 k w k 2 S Pr o of. Let w 2 W and B ~ w = PBw Then h ; Bw i = h ; B ~ w i (4.22) Let z 2 Ker S Bz = Bw B ~ w Then, since P is an orthogonal projection, k Bz k k Bw k : 90 PAGE 100 No w, from the denition of B and the P oincar e inequalit y w e obtain h k Bw k 2 2 h k w k 2 CH k w k 2 S : Also, since z 2 Ker S it is constan t on eac h @ n i and using the connectivit y assumption, w e ha v e the follo wing b y Lemma 24 k Bz k 2 S C h H k Bz k 2 (1 + log H=h ) 2 : T ogether this yields k Bz k 2 S C (1 + log H=h ) 2 k w k 2 S : By the denition of B ( Bw ) i on @ n i [ @ n j is a linear com bination (with bounded coecien ts) of (a bounded n um ber of) w k from all @ n k adjacen t to @ n i [ @ n j F rom Lemma 24 and the P oincar e-F riedric hs inequalit y (Lemma 27), k Bw k S C (1 + log( H=h )) k w k S ; 8 w 2 W: Finally summarizing, k B ~ w k S k Bw k S + k B T Bz k S C (1 + log H=h ) k w k S : F rom this, (4.22), and B T B = B the result follo ws. Theorem 45 Under the assumption of Section 4.3.1, the condition n um ber of the P ark's v arian t of FETI method with the Diric hlet preconditioner (3.39) satises = max ( QDPF ) min ( QDPF ) C 1 + log H h 2 : Pr o of. Lemmas 43 and 44 v erify the assumptions (i) and (ii) of Theorem 38 with C 1 = 1 and C 2 = C (1 + log H=h ) 2 respectiv ely 91 PAGE 101 4.5 Con v ergence Estimates for Plate Bending The follo wing appr oximate p ar ametric variational principle form ulated in [42 ] will allo w us to considerate to the plate bending problem and the biharmonic equation. Assumption 46 ([41, 42]) W e consider elemen ts with displacemen t and rotation degrees of freedom at the v ertices only and assume that there exist constan ts c 1 > 0, c 2 suc h that if the plate thic kness t satises 0 < t h then for eac h elemen t T the local stiness matrix K T satises c 1 K HCT T K T c 2 K HCT T (4.23) where K HCT T is the HCT elemen t lev el stiness matrix of the biharmonic equation [15], with the rotations in terpreted as deriv ativ es of the transv ersal displacemen t in the HCT elemen t. That is, as the thic kness of the plate goes to zero, the stiness matrix of the elemen t should be spectrally equiv alen t to that of the HCT elemen t for the biharmonic equation. W e refer to [14 ] and Section 4.1.2 for more details on HCT elemen ts. Here w e just summarize that HCT elemen ts are in C 1 and they use cubic splines for v alues on elemen t sides, linear in terpolation for normal deriv ativ es on the sides, and piecewise polynomial extension in to the elemen t in terior. In [41 ], Assumption 46 is v eried for the particular case of the DKT elemen t [5]. Assumption 46 also holds for the follo wing general class of non-loc king P 1 Reissner-Mindlin elemen ts. Theorem 47 ([41 ]) Assume that the energy functional for an elemen t T is spectrally equiv alen t to Z T jr j 2 dx + 1 t 2 + h 2 Z T j r u j 2 dx; (4.24) 92 PAGE 102 where h = diam ( T ), u 2 P 1 ( T ) is the transv ersal displacemen t, and 2 ( P 1 ( T )) 2 is the rotation. Then (4.23) holds. Elemen ts with the energy functional of the form (4.24) include the DKT elemen t as restated in [60 ]. It should be noted that for the related Timoshenk o beam elemen t, the thin limit is exactly the discretization b y cubic splines of the biharmonic equation [35]. 4.5.1 Assumptions Before pro ving a bound on the condition n um ber of the generalized FETI method, w e need to in troduce some specic assumptions. W e refer to the model problem (2.10) and using the spectral equiv alence, that w e discussed in the previous section, w e consider the biharmonic boundary v alue problem in a v ariational form. Find u 2 H 2 0 (n) suc h that a ( u; v ) = f ( v ) ; 8 v 2 H 2 0 (n) ; where a ( u; v ) = Z n @ 11 u@ 11 v + @ 12 u@ 12 v + @ 22 u@ 22 v; 8 u; v 2 H 2 0 (n) Let the domain n be divided in to non-o v erlapping subdomains n i i = 1 ; :::; N s that are shape regular of diameter O ( H ) (Denition 21). Without loss of generalit y w e suppose that H < 1. F urthermore, w e require that, for all r ; s = 1 ; : : :; N s if @ n r \ @ n s is an edge, then its length is larger than a certain prescribed fraction of the length of @ n r W e note that this condition implies the connectivit y condition as that of the original FETI. It is possible to pro v e the result using the general condition, but, for the sak e of simplicit y w e pro v e it using the stronger assumption. 93 PAGE 103 W e assume that the problem is discretized using reduced HCT elemen ts. The general case of plate bending then follo ws from spectral equiv alence of the local elemen t stiness matrices follo wing Assumption 46. Let V HCT h (n) denote the corresponding nite elemen t space, and h denote the c haracteristic elemen t size. Eac h subdomain n i is assumed to be a union of some of the elemen ts. The degrees of freedom are v alues of the transv ersal displacemen t and its deriv ativ es (rotations) at the nodal poin ts of the discretization. W e dene B as follo ws, cf., Fig. 4.3. F or eac h node x on subdomain in terface and eac h pair ( r ; s ) suc h that x 2 n r \ @ n s and @ n r and @ n s share an edge (i.e, do not meet only at a crosspoin t), Bw includes one sub v ector of three elemen ts, ( Bw ) rs ( x ) = rs ( w r ( x ) w s ( x )) ; where w r ( x ) and w s ( x ) denote sub v ectors of w con taining the three degrees of freedom associated with node x and subdomain n r and n s respectiv ely and rs is c hosen to be either 1 or 1. W e poin t out that rs does not depend on x 2 @ n r \ @ n s i.e., coecien ts rs are uniform along the edge bet w een n r and n s Note that, unlik e in the case described in Section 4.3.1, this denition implies redundan t constrain ts at subdomain crosspoin ts. Let us mak e the denition of the matrix C from Section 3.2.3 more precise. Let m = dim and f k j ; j = 1 ; : : :; n g be the the complete list of indices of Lagrange m ultipliers corresponding to conditions (expressed through B ) enforcing con tin uit y of the solution (but not the con tin uit y of the deriv ativ es) at crosspoin ts. Then C = [ c ij ] is the m n matrix satisfying c ij = 1 if i = k j 94 PAGE 104 c ij = 0 otherwise : 4.5.2 Discrete Norm Bounds Using Lemma 33 for subdomain stiness matrices K i for eac h subdomain and summing o v er all subdomains, w e obtain c 1 N s X i =1 jr I HCT R i w j 2 1 2 ; 2 ;@ n i k w k 2 S c 2 N s X i =1 jr I HCT R i w j 2 1 2 ; 2 ;@ n i 8 w 2 W ; (4.25) where the positiv e constan ts c 1 ; c 2 are independen t of the c haracteristic mesh size h and the subdomain diameter H The constan ts c 1 ; c 2 ma y depend on the regularit y of the shape of subdomains. F or denitions of Sobolev spaces, see, for example [1]. Throughout this section, c; C; c 1 ; c 2 ; c 3 ; c 4 ; c 5 and c 6 denote positiv e constan ts independen t of H and h The next t w o lemmas con tain the principal tec hnical estimates. Lemma 48 F or all 2 V 0 and all w 2 W suc h that Bw 2 V there exists a ~ w 2 W suc h that h ; B ~ w i = h ; Bw i and k ~ w k 2 S C (1 + log H=h ) k B T Bw k 2 S : where = 1, and = 0 if 1 2 B B T = I whic h happens when there are no nodes shared b y more than t w o subdomains. Pr o of. Let us rst pro v e that, in the general case, w e obtain 1. Let w 2 W and Bw 2 V That is Z T B T Bw = 0 and C T Bw = 0. W e dene ~ w = B T ( B B T ) + Bw Then, B ~ w = Bw By the triangle inequalit y w e ma y write k ~ w k S k 1 2 B T Bw k S + k 1 2 B T ( I ( 1 2 B B T ) + ) Bw k S : (4.26) Denote z = 1 2 B T ( I ( 1 2 B B T ) + ) Bw F rom the denition of B z is zero at all nodes that belong to at most t w o subdomains. The remaining 95 PAGE 105 nodes lie on subdomain crosspoin ts. A t ev ery suc h node, I HCT R i z ( x ) is a linear com bination of the en tries of B T Bw that correspond to the same node x and the coecien ts of the linear com binations are bounded only in terms of the n um ber of subdomains to whic h the node belongs. In addition, since C T Bw = 0, the transv ersal displacemen t componen ts of z at crosspoin ts are zero. Using (4.25) and Lemma 24 for the subdomain crosspoin t v ertices, one subdomain at a time, w e obtain k z k 2 S C (1 + log H h ) N s X i =1 1 H kr I HCT R i B T Bw k 2 0 ; 2 ;@ n i + jr I HCT R i B T Bw j 2 1 2 ; 2 ;@ n i This together with the P oincar e inequalit y Lemma 34, and (4.26) yields the result. If 1 2 B B T = I w e simply c hoose ~ w = 1 2 B T Bw = w No w w e deriv e the con v erse bound. Lemma 49 F or all 2 V 0 and w 2 W; w ? Ker S there is a ~ w 2 W suc h that B ~ w 2 V h ; Bw i = h ; B ~ w i ; and k B T B ~ w k 2 S C 1 + log H h 2 k w k 2 S : Pr o of. Let 2 V 0 w 2 W; w ? Ker S and consider B ~ w = P Bw + PFPC Since 2 V 0 it is easily v eried that h ; Bw i = h ; B ~ w i V ector can be found from the condition C T B ~ w = 0, whic h can be rewritten as h PFPC; C ~ i = h P Bw; C ~ i 8 ~ : By the denition of F w e can rewrite this as h S + B T PC; B T PC ~ i = h w; B T PC ~ i 8 ~ : 96 PAGE 106 Since S + is positiv e semidenite, the equation yields suc h that h S + 1 = 2 B T PC; S + 1 = 2 B T PC i = h S 1 = 2 w; S +1 = 2 B T PC i and, from the Cauc h y inequalit y k S + 1 = 2 B T PC k k S 1 = 2 w k Therefore, k S + B T PC k S k w k S : (4.27) W e need to estimate k B T B ~ w k S Let B ij be the matrix constructed from B b y zeroing out all the ro ws that do not correspond to the in terface conditions bet w een n i and n j Then, B T B ~ w = N s X i;j =1 ;i PAGE 107 Denote u = w + S + B T PC Then, B ~ w = P Bu The triangle inequalit y and (4.27) yield k u k S k w k S + k S + B T PC k S 2 k w k S : So no w w e only ha v e to pro v e that X k 1 H kr I HCT P Bu k 2 0 ; 2 ; k + jr I HCT P Bu j 2 1 2 ; 2 ; k c 3 k u k 2 S : W e rst use the triangle inequalit y to get jr I HCT P Bu j 1 2 ; 2 ; k jr I HCT Bu j 1 2 ; 2 ; k + jr I HCT ( I P ) Bu j 1 2 ; 2 ; k : F rom the denition of B eac h en try reduces to a linear com bination with coecien ts 1,constan t along ev ery edge. Since ( I P ) Bu 2 Im G I HCT ( I P ) Bu is a restriction of a linear function on ev ery edge, jr I HCT ( I P ) Bu j 1 2 ; 2 ; k = 0. F urthermore, X k jr I HCT Bu j 2 1 2 ; 2 ; k 2 N s X i =1 jr I HCT u j 2 1 2 ; 2 ;@ n i : (4.28) T o estimate the L 2 terms, w e again rst use the triangle inequalit y: kr I HCT P Bu k 0 ; 2 ; k 2 kr I HCT u k 0 ; 2 ; k + kr I HCT ( I P ) Bu k 0 ; 2 ; k : Since r I HCT ( I P ) Bu is a linear function on ev ery edge and P is an orthogonal projection, w e ha v e X k kr I HCT ( I P ) Bu k 2 0 ; 2 ; k c 4 h k ( I P ) Bu k 2 c 5 h k Bu k 2 2 c 5 h k u k 2 : F urthermore, Lemma 32 sho ws that h k u k 2 c 6 (1 + H 2 ) N s X i =1 kr I HCT u k 2 0 ; 2 ;@ n i Finally since u ? Ker S the P oincar e inequalit y (Lemma 34) and the equiv alence of the norms (4.25) conclude the estimate. 98 PAGE 108 4.5.3 Condition Num ber Estimate W e ha v e no w ev erything ready to pro v e the nal result. Theorem 50 Under the assumptions of Section 4.5.1, the condition n um ber of the generalized FETI method with the Diric hlet preconditioner (3.39) satises = max ( QDPF ) min ( QDPF ) C 1 + log H h r with r = 3, and r = 2 if there are no crosspoin ts bet w een more than t w o subdomains. Pr o of. Since V and V 0 can be replaced b y the factorspaces ~ V and ~ V 0 Lemmas 48 and 49 v erify the assumptions (i) and (ii), of Theorem 38, with C 1 = C (1 + log H=h ) = 0 or 1, and C 2 = C (1 + log H=h ) 2 respectiv ely 99 PAGE 109 Figure 4.3. Denition of B for plate bending t t t t t t t t t t t t t t t t t r r r r r r r r r r n 1 n 2 n 3 n 4 @ n 1 @ n 2 @ n 3 @ n 4 w 4 ( x A ) w 4 ( x B ) w 4 ( x C ) w 3 ( x C ) w 3 ( x D ) w 3 ( x E ) w 1 ( x A ) w 1 ( x B ) w 1 ( x C ) w 2 ( x C ) w 2 ( x D ) w 2 ( x E ) w 1 ( x F ) w 2 ( x F ) w 4 ( x G ) w 3 ( x G ) w 1 ( x A ) w 4 ( x A ) w 1 ( x B ) w 4 ( x B ) w 1 ( x C ) w 4 ( x C ) w 2 ( x C ) w 3 ( x C ) w 2 ( x D ) w 3 ( x D ) w 2 ( x E ) w 3 ( x E ) w 1 ( x F ) w 2 ( x F ) w 1 ( x C ) w 2 ( x C ) w 4 ( x C ) w 3 ( x C ) w 4 ( x G ) w 3 ( x G ) 100 PAGE 110 5. Computational Results W e refer to the paper b y F arhat et al. [32] for n umerical results for second order problems conrming our theoretical results. The paper also deals with parallel implemen tation issues and performance. Here w e summarize results for plate bending problems based on the paper b y Mandel, T ezaur and F arhat [53 ], due to F arhat. W e consider the plate bending problem on a unit square. The plate is discretized b y a uniform mesh of three-node triangular DKT plate elemen ts and subjected to a uniform pressure load. The thic kness of the plate is 10 3 the Y oung modulus E = 10 6 and the P oisson coecien t = 0 : 3. The FETI and generalized FETI methods are preconditioned b y the Diric hlet preconditioner (Section 3.3) [32, 34], and the follo wing stopping criterion is used k z k 1 k k f k = 10 3 (5.1) where z k 1 is the preconditioned residual in Algorithm 20. This condition is a good estimator of the global error as sho wn in [34]. W e construct sev eral meshes with dieren t mesh size h and sev eral mesh partitions with dieren t subdomain size H and report the performance of the original FETI and the generalized FETI. Three series of computational experimen ts are reported. The condition n um ber is denoted b y and the n um ber of iterations b y n it in the tables. The rst series of experimen ts sho ws the performance of the method 101 PAGE 111 for three dieren t n um bers of subdomains ( H = 1 = 2 ; 1 = 4 ; 1 = 8) and three differen t meshes corresponding to h = 1 = 10, h = 1 = 20, and h = 1 = 40. Results in T able 5.1) demonstrate that for a giv en H the condition n um ber of the original FETI gro ws fast with the mesh size h while that of the new FETI method is m uc h smaller and gro ws only w eakly with h F or large n um ber of subdomains ( N s = 64) the new FETI method con v erges about 7 times faster than the original one. In the second series of experimen ts (T able 5.2), the mesh size is xed to h = 1 = 120 (28800 elemen ts and 86400 degrees of freedom), and H is v aried bet w een H = 1 = 2 (4 subdomains) and H = 1 = 12 (144 subdomains). In this case, the condition n um bers of both FETI methods are sho wn to decrease with the n um ber of subdomains. This is an expected result because when h is xed and H is decreased, the size of the coarse problem increases for both algorithms. The generalized FETI sho ws only a small dependence on the n um ber of subdomains as predicted b y the theory F or large n um ber of subdomains, the new FETI method clearly outperforms the original FETI. Finally the subdomain problem size is xed to h=H = 1 = 15, and the n um ber of subdomains as w ell as the size of the global problem are increased. The performance results reported in T able 5.3 sho w that in this case too, the new FETI method outperforms signican tly the original one. Since solution time is ultimately the most importan t criterion for assessing performance, Both FETI methods w ere also benc hmark ed for the same plate bending problem described abo v e with 960000 degrees of freedom and 64 subdomains. The performance results obtained on a 64-processor IBM SP2 are summarized in T able 5.4. They sho w that ev en though the new FETI 102 PAGE 112 T able 5.1. Fixed n um ber of subdomains, series of rened meshes FETI Generalized FETI h H n it n it 2x2 subdomains, H= 1 2 1 10 18 2578 12 7.6 1 20 22 30101 15 12.6 1 40 26 409987 17 18.6 4x4 subdomains, H= 1 4 1 10 61 6795 21 11.5 1 20 86 84199 27 17 1 40 119 1038120 36 24.4 8x8 subdomains, H= 1 8 1 10 172 21707 25 13 1 20 247 275004 34 19.4 1 40 323 3920613 42 27.6 T able 5.2. Fixed global mesh 120x120, h= 1 120 series of rened mesh partitions Decomposition FETI Generalized FETI H h H n it n it 1 2 1 60 27 2079032 18 23.2 1 3 1 40 64 839240 29 22.4 1 4 1 30 104 391470 32 21 1 5 1 24 135 234504 33 19.9 1 6 1 20 164 160173 32 18.6 1 8 1 15 222 94285 31 16.6 1 10 1 12 255 63896 29 14.9 1 12 1 10 245 46921 27 13.6 103 PAGE 113 T able 5.3. Fixed local mesh (15x15, h H = 1 15 ), series of rened meshes and mesh partitions Decomposition FETI Generalized FETI H h n it n it 1 2 1 30 20 11088 13 10 1 3 1 45 49 19004 21 13.4 1 4 1 60 74 29041 25 14.6 1 5 1 75 109 40120 28 15.4 1 6 1 90 145 55068 29 15.9 1 8 1 120 222 94285 31 16.6 1 10 1 150 318 144556 32 16.9 104 PAGE 114 T able 5.4. P erformance results for 960000 degrees of freedom and 64 subdomains FETI Generalized FETI n it T otal time Time per iter. n it T otal time Time per iter. 314 265 s 0.8 s 45 105 s 1.1 s method consumes an amoun t of CPU time (55.5 s.) equiv alen t to that of 50 of its iterations to set up and preprocess the coarse problem (3.33), and ev en though eac h of its iterations is 1.3 times more expensiv e than an iteration of the original FETI method, the new FETI method is 2.5 times faster than the original one at solving the system of 960000 plate bending equations on a 64-processor IBM SP2. 105 PAGE 115 REFERENCES [1] R. A. Adams. Sob olev sp ac es Academic press, New Y ork, 1975. [2] V. I. Agoshk o v. P oincar e{Steklo v's operators and domain decomposition methods in nite dimensional spaces. In R. Glo winski, G. H. Golub, G. A. Meuran t, and J. P eriaux, editors, First International Symp osium on Domain De c omp osition Metho ds for Partial Dier ential Equations pages 73{112, Philadelphia, 1988. SIAM. [3] Iv o Babu sk a, Alan W. Craig, Jan Mandel, and Juhanni Pitk aran ta. Efcien t preconditioning for the p -v ersion nite elemen t method in t w o dimensions. SIAM J. Numer. A nal. 28:624{662, 1991. [4] N. S. Bakh v alo v and A. V. Kn y azev. Preconditioned iterativ e methods in a subspace for linear algebraic equations with large jumps in the coecien ts. In Domain De c omp osition Metho ds in Scientic and Engine ering Computing: Pr o c e e dings of the Seventh International Confer enc e on Domain De c omp osition v olume 180 of Contemp or ary Mathematics pages 157{162, Pro vidence, Rhode Island, 1994. American Mathematical Societ y [5] J. L. Batoz. An explicit form ulation for an ecien t triangular plate bending elemen t. Int. J. Numer. Meth. Engr g. 18:1077{1081, 1982. [6] P E. Bjrstad and J. Mandel. Spectra of sums of orthogonal projections and applications to parallel computing. BIT 31:76{88, 1991. [7] J. F. Bourgat, R. Glo winski, P LeT allec, and M. Vidrascu. V ariational form ulation and algorithm for trace operator in domain decomposition calculations. In T. F. Chan, R. Glo winski, J. P eriaux, and O. B. Widlund, editors, Domain De c omp osition Metho ds pages 3{16. SIAM, Philadelphia, 1989. [8] J. H. Bram ble, J. E. P asciak, and A. H. Sc hatz. An iterativ e method for elliptic problems on regions partitioned in to substructures. Math. Comp. 46:361{369, 1986. 106 PAGE 116 [9] J. H. Bram ble, J. E. P asciak, J. W ang, and J. Xu. Con v ergence estimates for product iterativ e methods with applications to domain decomposition. Math. Comp. 57:1{21, 1991. [10] James H. Bram ble, Joseph E. P asciak, and Alfred H. Sc hatz. The construction of preconditioners for elliptic problems b y substructuring, I. Math. Comp. 47(175):103{134, 1986. [11] James H. Bram ble, Joseph E. P asciak, and Alfred H. Sc hatz. The construction of preconditioners for elliptic problems b y substructuring, IV. Math. Comp. 53:1{24, 1989. [12] Sue Brenner. Tw o-lev el additiv e Sc h w arz preconditioners for plate bending problems. T o appear in Numer. Math. [13] T. F. Chan and T. P Mathew. The in terface probing tec hnique in domain decomposition. SIAM J. Matrix A nal. Applic. 13:212{238, 1992. [14] P G. Ciarlet. Basic error estimates for elliptic problems. In P .G. Ciarlet and J. L. Lions, editors, Handb o ok of Numeric al A nalysis v olume II, pages 17{352. North-Holland, Amsterdam, 1989. [15] R. W. Clough and J. L. T oc her. Finite elemen t stiness matrices for analysis of plate bending. In Pr o c. 1965 Conf. Matrix Metho ds Struct. Me ch., Wright-Patterson AFB, Ohio, AFFDL-TR-66-80 pages 515{546, 1966. [16] La wrence Co wsar. P ersonal comm unication, 1992. [17] Q. V. Dinh, R. Glo winski, J. He, V. Kw oc k, T. W. P an, and J. P eriaux. Lagrange m ultiplier approac h to ctitious domain methods: application to ruid dynamics and electro{magnetics. In D. E. Key es, T. F. Chan, G. A. Meuran t, J. S. Scroggs, and R. G. V oigt, editors, Fifth International Symp osium on Domain De c omp osition Metho ds for Partial Dier ential Equations pages 151{194, Philadelphia, 1992. SIAM. [18] Q. V. Dinh, R. Glo winski, and J. P eriaux. Solving elliptic problems b y domain decomposition methods with applications. In G. Birkho and A. Sc hoenstadt, editors, Elliptic Pr oblem Solvers II Academic Press, 1984. [19] M. R. Dorr. On the discretization of in terdomain coupling in elliptic boundary{v alue problems. In T. F. Chan, R. Glo winski, J. P eriaux, and 107 PAGE 117 O. B. Widlund, editors, Domain De c omp osition Metho ds pages 17{37. SIAM, Philadelphia, 1989. [20] M. Dryja. An additiv e Sc h w artz algorithm for t w o{and three{ dimensional nite elemen t elliptic problems. In T. F. Chan, R. Glo winski, J. P eriaux, and O. B. Widlund, editors, Domain De c omp osition Metho ds pages 168{ 172. SIAM, Philadelphia, 1989. [21] M. Dryja and O. B. Widlund. Domain decomposition algorithms with small o v erlap. SIAM J. Sci. Comput. 15:604{620, 1994. [22] Maksymilian Dryja. A method of domain decomposition for 3-D nite elemen t problems. In Roland Glo winski, Gene H. Golub, G erard A. Meuran t, and Jacques P eriaux, editors, First International Symp osium on Domain De c omp osition Metho ds for Partial Dier ential Equations Philadelphia, P A, 1988. SIAM. [23] Maksymilian Dryja, Barry F. Smith, and Olof B. Widlund. Sc h w arz analysis of iterativ e substructuring algorithms for elliptic problems in three dimensions. SIAM J. Numer. A nal. 31(6):1662{1694, Decem ber 1994. [24] Maksymilian Dryja and Olof B. Widlund. T o w ards a unied theory of domain decomposition algorithms for elliptic problems. In T on y Chan, Roland Glo winski, Jacques P eriaux, and Olof Widlund, editors, Thir d International Symp osium on Domain De c omp osition Metho ds for Partial Dier ential Equations, held in Houston, T exas, Mar ch 20-22, 1989 SIAM, Philadelphia, P A, 1990. [25] Maksymilian Dryja and Olof B. Widlund. Sc h w arz methods of NeumannNeumann t ype for three-dimensional elliptic nite elemen t problems. Comm. Pur e Appl. Math 48:121{155, 1995. [26] C. F arhat, P S. Chen, and J. Mandel. Scalable Lagrange m ultiplier based domain decomposition method for time-dependen t problems. Int. J. Numer. Meth. Engr g. 38:3831{3853, 1995. [27] C. F arhat and F.-X. Roux. A method of nite elemen t tearing and interconnecting and its parallel solution algorithm. Int. J. Numer. Meth. Engng. 32:1205{1227, 1991. [28] Charbel F arhat. Lagrange m ultiplier based divide and conquer nite elemen t algorithm. J. Comput. Sys. Engr g. 2:149{156, 1991. 108 PAGE 118 [29] Charbel F arhat, P o-Sh u Chen, Jan Mandel, and F rancois-Xa vier Roux. The t w o-lev el FETI method P art II: Extension to shell problems, parallel implemen tation and performance results. Comp. Meth. Appl. Mec h. Engrg, submitted, 1995. [30] Charbel F arhat, P o-Sh u Chen, F rancois Roux, and Jan Mandel. The t w o-lev el FETI method-part ii: Extension to shell problems, parallel implemen tation, and performance results. CU-CAS Report 96-10, Cen ter for Aerospace Structures, Univ ersit y of Colorado, 1996. [31] Charbel F arhat and Jan Mandel. The t w o-lev el FETI method for static and dynamic plate problems P art I: An optimal iterativ e solv er for biharmonic systems. T ec hnical Report UCB/CAS Report CU-CAS-95-23, Cen ter for Aerospace Structures, Univ ersit y of Colorado at Boulder, 1995. Comp. Meth. Appl. Mec h. Engrg, submitted. [32] Charbel F arhat, Jan Mandel, and F rancois-Xa vier Roux. Optimal conv ergence properties of the FETI domain decomposition method. Comput. Metho ds Appl. Me ch. Engr g. 115:367{388, 1994. [33] Charbel F arhat and F rancois-Xa vier Roux. An uncon v en tional domain decomposition method for an ecien t parallel solution of large-scale nite elemen t systems. SIAM J. Sci. Stat. Comput. 13:379{396, 1992. [34] Charbel F arhat and F rancois-Xa vier Roux. Implicit parallel processing in structural mec hanics. Comput. Me ch. A dvanc es 2:1{124, 1994. [35] Carlos A. F elippa. A surv ey of parametrized v ariational principles and application to computational mec hanics. Comput. Metho ds Appl. Me ch. Engr g 113:109{139, 1994. [36] R. Glo winski and M. F. Wheeler. Domain decomposition and mixed nite elemen t methods for elliptic problems. In R. Glo winski, G. H. Golub, G. A. Meuran t, and J. P eriaux, editors, First International Symp osium on Domain De c omp osition Metho ds for Partial Dier ential Equations pages 144{172, Philadelphia, 1988. SIAM. [37] Roland Glo winski and Mary F. Wheeler. Domain decomposition and mixed nite elemen t methods for elliptic problems. In Roland Glo winski, Gene H. Golub, G erard A. Meuran t, and Jacques P eriaux, editors, First International Symp osium on Domain De c omp osition Metho ds for Partial Dier ential Equations Philadelphia, P A, 1988. SIAM. 109 PAGE 119 [38] Gene H. Golub and Charles F. V an Loan. Matrix Computations Johns Hopkins Univ. Press, 1989. Second Edition. [39] Gene H. Golub and Charles F. V an Loan. Matrix Computations John Hopkins Univ ersit y Press, second edition, 1989. [40] Manoel R. Justino, K. C. P ark, and Carlos A. F elippa. An algebraically partitioned FETI method for parallel structural analysis: Implemen tation and n umerical performance ev aluation. International Journal of Numeric al Metho ds in Engine ering. 40:2739{2758, 1997. [41] P atric k Le T allec, Jan Mandel, and Marina Vidrascu. A NeumannNeumann domain decomposition algorithm for solving plate and shell problems. SIAM J. Numer. Anal., to appear. [42] P atric k Le T allec, Jan Mandel, and Marina Vidrascu. Balancing domain decomposition for plates. Contemp or ary Mathematics 180:515{524, 1994. Proceedings of the 7th In ternational Symposium on Domain Decomposition Methods, P enn State, No v em ber 1993. [43] P LeT allec, Y.-H. de Roec k, and M. Vidrascu. Domain{decomposition methods for large linearly elliptic three dimensional problems. J. Comput. Appl. Math. 34:93{117, 1991. [44] Pierre Louis Lions. On the Sc h w arz alternating method. I. In Roland Glo winski, Gene H. Golub, G erard A. Meuran t, and Jacques P eriaux, editors, First International Symp osium on Domain De c omp osition Metho ds for Partial Dier ential Equations Philadelphia, P A, 1988. SIAM. [45] J. Mandel. Balancing domain decomposition. Comm. Numer. Meth. Engr g. 9:233{241, 1993. [46] Jan Mandel. Balancing domain decomposition. Comm. in Numeric al Metho ds in Engr g. 9:233{241, 1993. [47] Jan Mandel. Hybrid domain decomposition with unstructured subdomains. In Y uri A. Kuznetso v, Jacques P eriaux, Alo Quarteroni, and Olof B. Widlund, editors, Domain De c omp osition Metho ds in Scienc e and Engine ering: The Sixth International Confer enc e on Domain De c omp osition v olume 157. AMS, 1994. Held in Como, Italy June 15{19,1992. [48] Jan Mandel and Marian Brezina. Balancing domain decomposition: Theory and computations in t w o and three dimensions. UCD/CCM Report 2, 110 PAGE 120 Cen ter for Computational Mathematics, Univ ersit y of Colorado at Den v er, 1993. [49] Jan Mandel and Marian Brezina. Balancing domain decomposition for problems with large jumps in coecien ts. Mathematics of Computation 65:1387{1401, 1996. [50] Jan Mandel and Radek T ezaur. Con v ergence of a substructuring method with Lagrange m ultipliers. Numerische Mathematik 73:473{487, 1996. [51] Jan Mandel, Radek T ezaur, and Charbel F arhat. Optimal Lagrange m ultiplier based domain decomposition method for plate bending problems. UCD/CCM Report 61, Cen ter for Computational Mathematics, Univ ersit y of Colorado at Den v er, 1995. [52] Jan Mandel, Radek T ezaur, and Charbel F arhat. Optimal Lagrange m ultiplier based domain decomposition method for plate bending problems. UCD/CCM Report 61, Cen ter for Computational Mathematics, Univ ersit y of Colorado at Den v er, 1995. [53] Jan Mandel, Radek T ezaur, and Charbel F arhat. A scalable substructuring method b y Lagrange m ultipliers for plate bending problems. UCD/CCM Report 84, Cen ter for Computational Mathematics, Univ ersit y of Colorado at Den v er, 1996. Submitted to SIAM J. Sci. Comput. [54] A. M. Matsokin and S. V. Nepomn y asc hikh. A Sc h w arz alternating method in subspaces. In. V uzov 10:61{66, 1985. Also in So viet Mathematics, 10 (1985), pp. 78{84. [55] S. V. Nepomn y asc hikh. Domain de c omp osition and the Schwarz metho d in subsp ac e for appr oximate solution of elliptic b oundary pr oblems PhD thesis, Computing Cen ter of the Siberian Branc h of the USSR Academ y of Sciences, No v osibirsk, Russia, 1986. [56] Jind ric h Ne cas. L es m etho des dir e ctes en th eorie des equations elliptiques Academia, Prague, 1967. [57] J. T. Oden, Abani P atra, and Y usheng F eng. P arallel domain decomposition solv er for adaptiv e hp nite elemen t methods. T ec hnical Report TICIAM 94-11, Univ ersit y of T exas at Austin, 1994. [58] K. C. P ark and Carlos A. F elippa. A v ariational framew ork for solution method dev elopmen ts in structural mec hanics. CU-CAS Report 96-21, 111 PAGE 121 Cen ter for Aerospace Structures, Univ ersit y of Colorado, 1996. T o appear in Journal of Applied Mathematics. [59] K. C. P ark, Manoel R. Justino, and Carlos A. F elippa. An algebraically partitioned FETI method for parallel structural analysis: Algorithm description. International Journal of Numeric al Metho ds in Engine ering. 40:2717{2737, 1997. [60] Juhanni Pitk aran ta. On a simple nite elemen t method for plate bending problems. In W olfgang Hac kbusc h and Kristian Witsc h, editors, Numeric al T e chniques in Continuum Me chanics n um ber 16 in Notes on Numerical Fluid Mec hanics, pages 84{101. View eg, Braunsc h w eig/Wiesbaden, 1987. Proceedings of 2nd GAMM-Seminar, Kiel, Jan uary 1986. [61] Y.-H. de Roec k and P LeT allec. Analysis and test of a local domain decomposition. In R. Glo winski, Y u. A. Kuznetso v, G. A. Meuran t, J. P eriaux, and O. B. Widlund, editors, F ourth International Symp osium on Domain De c omp osition Metho ds for Partial Dier ential Equations pages 112{128, Philadelphia, 1991. SIAM. [62] F.-X. Roux. Acceleration of the outer conjugate gradien t b y reorthogonalization for a domain decomposition method with Lagrange m ultiplier. In T. F. Chan, R. Glo winski, J. P eriaux, and O. B. Widlund, editors, Thir d International Symp osium on Domain De c omp osition Metho ds for Partial Dier ential Equations pages 314{321, Philadelphia, 1990. SIAM. [63] H. A. Sc h w arz. Uber einige Abbildungsaufgaben. Ges. Math. A bh. 11:65{ 83, 1869. [64] Radek T ezaur, P etr V an ek, and Marian Brezina. Tw o-lev el method for solids on unstructured meshes. UCD/CCM Report 73, Cen ter for Computational Mathematics, Univ ersit y of Colorado at Den v er, 1996. Submitted to SIAM J. Sci. Comput. [65] P V an ek, J. Mandel, and M. Brezina. Algebraic m ultigrid b y smoothed aggregation for second and fourth order elliptic problems. T o appear in Computing. [66] P etr V an ek and Jitk a K r zk o v a. Tw o{lev el method on unstructured meshes with con v ergence rate independen t of the coarse-space size. T ec hnical Report 35, Cen ter of Computational Mathematics/UCD, 1995. 112 PAGE 122 [67] P etr V an ek, Jan Mandel, Marian Brezina, and Radek T ezaur. Tw o-lev el substructuring tec hnique for problems with coecien t jumps within subdomains. In preparation, 1996. [68] P etr V an ek, Radek T ezaur, Marian Brezina, and Jitk a K r zk o v a. Tw o{ lev el method with coarse space size independen t con v ergence. In R. Glo winski, J. P eriaux, Z. Shi, and O. Widlund, editors, Domain Dec omp osition Metho ds in Scienc es and Engine ering John Wiley & Sons Ltd., New Y ork, N.Y., 1995. T o appear. [69] Olof B. Widlund. An extension theorem for nite elemen t spaces with three applications. In W olfgang Hac kbusc h and Kristian Witsc h, editors, Numeric al T e chniques in Continuum Me chanics pages 110{122, Braunsc h w eig/Wiesbaden, 1987. Notes on Numerical Fluid Mec hanics, v. 16, F riedr. View eg und Sohn. Proceedings of the Second GAMM-Seminar, Kiel, Jan uary 1986. [70] J. Xu. Iterativ e methods b y space decomposition and subspace correction: A unifying approac h. SIAM R eview 34:581{613, 1992. [71] Xuejun Zhang. Domain decomposition algorithms for the biharmonic Diric hlet problem. In T on y F. Chan, Da vid E. Key es, G erard A. Meuran t, Jerey S. Scroggs, and Robert G. V oigt, editors, Fifth International Symp osium on Domain De c omp osition Metho ds for Partial Dier ential Equations Philadelphia, P A, 1992. SIAM. 113 |