Additive Average Schwarz Method for Elliptic Mortar Finite Element Problems with Highly Heterogeneous Coefficients
Ali Khademi, Leszek Marcinkowski, Sanjib Kumar Acharya, Talal Rahman
NNoname manuscript No. (will be inserted by the editor)
Additive Average Schwarz Method for EllipticMortar Finite Element Problems with HighlyHeterogeneous Coefficients
Ali Khademi · Leszek Marcinkowski · Sanjib Kumar Acharya · Talal Rahman
Received: date / Accepted: date
Abstract
In this paper, we extend the additive average Schwarz method tosolve second order elliptic boundary value problems with heterogeneouscoefficients inside the subdomains and across their interfaces by the mortartechnique, where the mortar finite element discretization is on nonmatchingmeshes. In this two-level method, we enrich the coarse space in two differentways, i.e., by adding eigenfunctions of two variants of the generalizedeigenvalue problems. We prove that the condition numbers of the systems ofalgebraic equations resulting from the extended additive average Schwarzmethod, corresponding to both coarse spaces, are of the order O ( H/h ) andindependent of jumps in the coefficients, where H and h are the meshparameters. Keywords
Domain decomposition method · Additive Schwarz method · Generalized eigenvalue problem · Mortar finite elements
Ali KhademiDepartment of Computer science, Electrical engineering and Mathematical sciences,Western Norway University of Applied Sciences, P.O. Box 7030, Bergen, NorwayE-mail: [email protected], [email protected] MarcinkowskiDepartment of Mathematics, Warsaw University, PolandE-mail: [email protected] Kumar AcharyaInstitute of Chemical Technology Mumbai, Indian Oil Campus, Odisha, Bhubaneswar, IndiaE-mail: [email protected] RahmanDepartment of Computer science, Electrical engineering and Mathematical sciences,Western Norway University of Applied Sciences, P.O. Box 7030, Bergen, NorwayE-mail: [email protected] a r X i v : . [ m a t h . NA ] F e b Ali Khademi et al.
Domain decomposition methods are efficient and powerful iterative methodsto solve large algebraic systems arising from a finite element discretization ofelliptic boundary value problems [32,35,37]. They can also be regarded as aprocedure of producing preconditioners for other iterative methods, such as theconjugate gradient method, for achieving fast convergence. In both approaches,to solve an original problem defined on a bounded Lipschitz domain Ω = ∪ Ni =1 Ω i , it is equivalent to solve many subproblems defined locally on thesubdomains Ω i in parallel. To obtain fast convergence, Dryja and Widlund [12],and Matsonki and Nepomnyaschikh [27] proposed to add one global problemand introduced the additive Schwarz methods to solve the global and localproblems in parallel.Corresponding to the global and local problems, different coarse and finespaces can be constructed, for instance, see [5,7,11,29]. It is worth mentioningthat coarse spaces are more important than fine spaces since they have a keyrole in the central of scalability for domain decomposition methods. Therefore,in [7] one of the simplest and efficient ways to construct coarse space, calledthe average coarse space, was proposed. Consequently, the two-level additiveaverage Schwarz method was introduced and developed to solve many kindsof elliptic problems with the continuous and discontinuous coefficients, see[8,14,17,24]. Hence, an extension of the additive average Schwarz method tosolve elliptic model problems arising from many applications such as compositematerials with highly heterogeneous coefficients is of particular interest in thispaper because it has a high-level performance.In general, in terms of distributions of the coefficients, we can classifyelliptic problems with the heterogeneous coefficients into two classes, i.e.,elliptic problems with jumps in the coefficients only inside the subdomainsand ones with jumps in the coefficients both inside the subdomains and onthe subdomain interfaces. For the first class, where the heterogeneouscoefficients are piecewise constants for subdomains Ω i , i = 1 , . . . , N , theadditive Schwarz method was developed and analyzed in [7,9,30,36] andreferences therein. For the second class with the large variation in theheterogeneous coefficients, the classical coarse spaces lead the conditionnumbers of the preconditioned systems to blow up. Consequently, theconvergence rates of iterative methods will deteriorate. To alleviate thisdifficulty, the coarse spaces can be enriched by combining their structurewith spectral spaces.The idea of the coarse spectral space was introduced in [6] and extendedas the spectral algebraic multigrid method in [15]. The spectral constructionof this new space is achieved by solving the generalized eigenvalue problemslocally. Due to its crucial role in preventing the impact of large jumps in thecoefficients on bounds of the condition numbers of preconditioners, the newcoarse space has received considerable attention. As a result, it is extensivelydeveloped for overlapping Schwarz methods [16,18,19,34], nonoverlappingadditive Schwarz method [26], balancing domain decomposition methods [20, itle Suppressed Due to Excessive Length 3 i , is the range of a linearoperator defined on Ω i such that it is either the nodal values of a function u ∈ V h inside Ω i or the average of nodal values of the function u on themortar and nonmortar sides of Ω i , where V h is a finite space of P1conforming elements defined on a fine triangulation of Ω and vanishing on ∂Ω . The second subspace has a particular spectral structure. To obtain thebasis functions for this subspace, we solve the generalized eigenvalueproblems restricted to each subdomain as in [26].To define the proper generalized eigenvalues problems, we requiredetermining minimum values of the coefficients over each subdomain’striangulation to estimate the condition number of additive average Schwarzpreconditioners independent of the large eigenvalues caused by the largejumps in the coefficients. Motivated by the ideas from [26], we consider twodifferent types of the generalized eigenvalue problems based on eitherminimum of the coefficients over the whole subdomain or just minimum ofones over the layer connected to the boundary of the subdomain with onevertex or with one edge of the triangles inside of the subdomain. Solvingthese generalized eigenvalue problems lead to finding orthogonal basisfunctions enriching the coarse spaces. With these new coarse spaces, weprove that the condition numbers of the produced preconditioners are of theorder O ( H/h ) and independent of the number of subdomains. Due to thedefinition of the second type layer, it has faster performance than the firstone through the implementations with numerical-software packages. See alsoSection 6 for numerical results.The outline of this paper is as follows. In Section 2, after introducing adiscrete problem, we define the mortar condition and the space of basis
Ali Khademi et al. functions satisfying that condition. Furthermore, several figures related tothose functions are also given there. Section 3 is devoted to introducing theadditive average Schwarz method, where the average interpolation operatorhas two different types. This operator consists of the natural extension of thestandard average interpolation operator for mortar case and orthogonaloperators defined in the next section. In Section 4, the generalized eigenvalueproblems in terms of how the minimum values of jumps in the coefficientsover the subdomains’ triangulations can be defined are introduced. InSection 5, following the standard additive Schwarz framework [32], we drivean optimal estimate of the condition number of the produced preconditionerswith the aid of removing bad eigenvalues, which are influenced by the largejumps in the coefficients. Finally, to verify our theoretical results’ validity, inSection 6 some numerical experiments are reported.
Let Ω ⊂ R be a bounded Lipschitz domain with a nonoverlapping partition { Ω i } Ni =1 of polygonal subdomains such that Ω = ∪ Ni =1 Ω i . We consider anelliptic model problem defined on Ω : Find u ∗ ∈ H ( Ω ) such that a ( u ∗ , v ) = f ( v ) , v ∈ H ( Ω ) , (1)where a ( u, v ) = N (cid:88) i =1 a i ( u, v ) = N (cid:88) i =1 ( α i ∇ u, ∇ v ) L ( Ω i ) and f ( v ) = (cid:90) Ω f v dx = N (cid:88) i =1 (cid:90) Ω i f v dx. Here, f ∈ L ( Ω ), α ∈ L ∞ ( Ω ) and α i ( x ) is the restriction of α ( x ) over Ω i . Further, we assume that there exists a positive constant α such that α ( x ) > α .The partition { Ω i } Ni =1 forms a coarse triangulation of Ω with the meshparameter H = max { H i , i = 1 , . . . , N } , where H i are diameters of Ω i . Weassume this partition to be geometrically conforming, i.e., ∂Ω i ∩ ∂Ω j ( i (cid:54) = j ) isa vertex or a whole edge of both subdomains Ω i and Ω j or is empty. Further,we denote the set of all vertices of such coarse triangulation, except thosebelonging to ∂Ω , by N H . We also denote the triangulation of the subdomain Ω i by T i which consists of triangles satisfying the shape regular property [10]inside of Ω i , and quasi-uniform triangles touching ∂Ω i with the mesh size h i .Further, we denote the set of all internal nodes of the fine triangulation T i by ı i . We define the product space X h on the computational domain Ω by X h ( Ω ) = X ( Ω ) × . . . × X N ( Ω N ) , itle Suppressed Due to Excessive Length 5 Ω i Ω j γ m ( i ) δ m ( j ) Fig. 1
The mortar and nonmortar sides of subdomains Ω i and Ω j which are denoted by γ m ( i ) and δ m ( j ) , respectively. where X i ( Ω i ) , i = 1 , . . . , N are the finite element spaces of the continuouspiecewise linear functions defined on T i and vanishing on ∂Ω ∩ ∂Ω i . We denoteall nodal points on (cid:83) Ni =1 T i except those on ∂Ω by N Ω . Further, we denote theset of basis functions associated with the set of nodal points N Ω by { φ l } l ∈N Ω .Due to independent triangulations inside of each subdomain, on each side ofthe interface Γ ij = Ω i ∩ Ω j we may have different discretization (cf. Figure 1).We select one side of Γ ij as the mortar side, and the other side as the nonmortarside denoted by γ m ( i ) and δ m ( j ) , respectively. It is obvious that Γ ij = γ m ( i ) = δ m ( j ) . Further, we denote the nodes on the mortar and nonmortar sides by m , m , . . . , m n m +1 and s , s , . . . , s n s +1 , respectively.Let W h i ( Γ ij ) and W h j ( Γ ij ) denote the restrictions of X i ( Ω i ) and X j ( Ω j )onto Γ ij , respectively. Now, the nonmatching grids on the subdomain interfacesimpose discontinuities for the functions belong to X h . Therefore, we need todefine a week continuity condition. To this end, we first define the projection Π m ( u i , Tr v j ) : L ( δ m ( j ) ) → W h j ( δ m ( j ) ) by (cid:90) δ m ( j ) (cid:16) u i | γm ( i ) − Π m ( u i , Tr v j ) (cid:17) ψ dx = 0 , γ m ( i ) = Γ ij = δ m ( j ) , (2)for all functions ψ ∈ M h j ( δ m ( j ) )) = span { ϕ l } n s +1 l =0 and Π m ( u i , Tr v j ) | ∂δm ( j ) = v j | ∂δm ( j ) , where u i , Tr v j , and M h j ( δ m ( j ) ) are the restriction of u into each subdomains Ω i , i = 1 , . . . , N , the trace of v j , and a subspace of W h j ( δ m ( j ) ) with constantvalues on the elements touching ∂δ m ( j ) , more precisely, at two end points of δ m ( j ) , respectively. Now, we say a function u h = { u i } Ni =1 ∈ X h satisfies themortar condition on δ m ( j ) , if (2) holds. We denote the mortar space by V h interms of the mortar condition, i.e., V h = { u h , v h ∈ X h | Π m ( u h , Tr v h ) = 0 , ∀ Γ ij , i, j = 1 , . . . , N } . Ali Khademi et al. Ω i Ω j x k Π m (I) Ω i Ω j x k Π n Π m γ n ( i ) γ m ( i ) (II) Ω i Ω j Π n γ n ( p ) x k γ m ( i ) Π m (III) Ω i Ω j x k (IV) γ m ( i ) Π m γ n ( j ) Π n Fig. 2
The node x k is represented by the black and thick dot, where the value of the basisfunction φ k ( x ) associated with x k is 1. Figures (I)-(IV) illustrate the different positions of x k . To see the structure of the basis functions spanning the mortar space V h ,let φ ( i ) k ( x ) be a basis function defined on Ω i and define V h = span { φ k } , where φ k ( x ) is a basis function associated with a node x k . Due to different positionsof x k , φ k ( x ) takes different forms as follows (cf. Figure 2).1. x k ∈ ı i : φ k ( x ) = φ ( i ) k ( x ) . x k ∈ { m , . . . , m n m } : φ k ( x ) = φ ( i ) k ( x ) , on Ω i ,Π m ( φ ( i ) k ( x ) , x ) , on δ m ( j ) , where γ m ( i ) = δ m ( j ) , , otherwise.3. x k ∈ N H : itle Suppressed Due to Excessive Length 7 (a) x k is a point of intersection between two mortar sides γ m ( i ) and γ n ( i ) : φ k ( x ) = φ ( i ) k ( x ) , on γ m ( i ) and γ n ( i ) ,Π m ( φ ( i ) k , x ) , on δ m ( j ) ,Π n ( φ ( i ) k , x ) , on δ n ( j ) , , otherwise.(b) x k is a point of intersection between mortar side γ m ( i ) and nonmortarside δ n ( i ) : φ k ( x ) = φ ( i ) k ( x ) , on γ m ( i ) ,Π m ( φ ( i ) k , x ) , on δ m ( j ) ,Π n (0 , Tr φ ( i ) k )( x ) , on δ n ( i ) , , otherwise.(c) x k is a point of intersection between two nonmortar sides δ m ( i ) and δ n ( i ) : φ k ( x ) = φ ( i ) k ( x ) , on γ m ( j ) and γ n ( j ) ,Π m (0 , Tr φ ( i ) k )( x ) , on δ m ( i ) ,Π n (0 , Tr φ ( i ) k )( x ) , on δ n ( i ) , , otherwise.We now express the main problem as in the following form: Find u ∗ h = { u hi } Ni =1 ∈ V h such that a h ( u ∗ h , v h ) = f ( v h ) ∀ v h ∈ V h . (3)In what follows, we consider the following matrix representation of thelinear systems arising from the discretization of (3). Av = f , where v is the vector of all unknown coefficients defined on Ω . Further, weconsider the submatrices A Ω i = R i AR Ti , where R i , i = 1 , . . . , N are therestriction matrices such that v i = R i v are the vectors of coefficients definedon Ω i \ ∂Ω i .Employing the mortar condition, we can compute some coefficients of v in terms of other ones. More precisely, where Γ ij = γ m ( i ) = δ m ( j ) assume ν s = ( v h ( s i )) n s i =1 , ν c = ( v h ( s ) , v h ( s n s +1 )) T , ν m = ( v h ( m i )) n m +1 i =0 and considerthe following matrix representations as in [13]. M γ := (cid:0) ( ϕ s i , φ m j ) L ( γ ) (cid:1) n s n m +1 i =1; j =0 ∈ R n s × ( n m +2) , S γ := (cid:0) ( ϕ s i , φ s j ) L ( γ ) (cid:1) n s n s i =1; j =1 ∈ R n s × n s , C γ := (cid:0) ( ϕ s i , φ m j ) L ( γ ) (cid:1) n s i =1; j =0 ,n s +1 ∈ R n s × . Ali Khademi et al.
Hence ν s = S − γ ( M γ ν m − C γ ν c ) . Consequently, in the next sections, we will focus only on all mortar, corner,and interior nodes.
Let V ( i ) , i = 1 , . . . , N be decompositions of V h that are restrictions of V h to Ω i with zero on ∂Ω i and on the subdomains Ω j , j (cid:54) = i . Further, let V ( type )0 standsfor two different types of the coarse spaces, distinguished by two notations Iand II, such that V h = V ( type )0 + N (cid:88) i =1 V ( i ) , type ∈ { I , II } , where V ( type )0 = I V h + N (cid:88) i =1 W typei , type ∈ { I , II } . We define the operator I and the spaces W typei , type ∈ { I , II } in the nextsection. Here, we use the exact bilinear form for all our subproblems. Thus,we have a ( u h , v h ) = a i ( u hi , v hi ) , i = 1 , . . . , N, where u h = { u hi } Ni =1 ∈ V ( i ) and v h = { v hi } Ni =1 ∈ V ( i ) . We define the projection-like operators T ( i ) : V h → V ( i ) i = 1 , . . . , N and T ( type )0 : V h → V ( type )0 suchthat for u h ∈ V h , T ( i ) u h and T ( type )0 u h are the solutions of a ( T ( i ) u h , v h ) = a ( u h , v h ) , v h ∈ V ( i ) , i = 1 , . . . , N and a ( T ( type )0 u h , v h ) = a ( u h , v h ) , v h ∈ V ( type )0 , type ∈ { I , II } . Now, the additive Schwarz operator T type : V h → V h is T type = T ( type )0 + N (cid:88) i =1 T ( i ) , type ∈ { I , II } and the problem (3) can be written as follow. T type u ∗ h = g type , type ∈ { I , II } , where g type = g ( type )0 + (cid:80) Ni =1 g ( i ) with g ( type )0 = T ( type )0 u ∗ h and g ( i ) = T ( i ) u ∗ h . itle Suppressed Due to Excessive Length 9 In this section, we design two different coarse spaces for the additive averageSchwarz method. To this end, we first denote the sets of nodal nodes of allmortar and nonmortar sides, and also all interior nodes of all subdomains Ω i , i = 1 , . . . , N by N m , N s , and N i , respectively. Now, the average interpolationoperator I : V h → V h for the mortar discretization as in [8] has the followingstructure. I u h ( x ) = (cid:40) u hi ( x ) x ∈ N m ∪ N s ,u hi x ∈ N i , i = 1 , . . . , N, where u hi is the average value of u hi over Ω i , i.e., u hi = 1 µ i ( δ,γ ) (cid:88) γ m ( i ) ⊂ ∂Ω i u hγ m ( i ) + (cid:88) δ m ( i ) ⊂ ∂Ω i ,γ m ( j ) = δ m ( i ) u hγ m ( j ) , (4)where u hγ m ( i ) = 1 | γ m ( i ) | (cid:90) γ m ( i ) u hi ds and µ i ( δ,γ ) is the number of all mortar and nonmortar sides of Ω i , and | γ m ( i ) | is the length of γ m ( i ) . To express the interpolation operator I in terms of thematrix form denoted by R , let N c be a set of all nodal nodes at the endof all nonmortar sides of Ω and let also I ( c ) ∈ R N c × N c and I ( m ) ∈ R N m × N m be identity matrices, where N c = dim( N c ) and N m = dim( N m ). Further,consider H = diag( H , . . . , H N ) ∈ R ( N c + N m ) × N i , where H i = ( u hi ) n m n i l =1; j =1 , N i = dim ( N i ), and n i = dim ( ı i ). Hence R = (cid:20) diag( I ( c ) , I ( m ) ) H (cid:21) ∈ R ( N c + N m ) × N Ω , N Ω = N c + N m + N i . Ω i Ω Bi Fig. 3
The boundary layer Ω Bi ⊂ Ω i is highlighted by the colorful triangles. To define two different coarse spaces, we first introduce the boundary layers Ω Bi , i = 1 , . . . , N , where Ω Bi ⊂ Ω i is the sum of all triangles such as τ ∈ T i such that ∂τ ∩ ∂Ω i (cid:54) = φ (cf. Figure 3). We then define the following local minimums of the coefficients over the subdomains Ω i : α i := min x ∈ Ω i α ( x ) , α i,B := min x ∈ Ω Bi α ( x ) , i = 1 , . . . , N. We now proceed to construct the second part of V ( type )0 , i.e., W typei , type ∈ { I , II } , i = 1 , . . . , N which are the spaces of adaptively choseneigenfunctions of specially constructed generalized eigenvalue problemsdefined locally in each subdomain and extended by zero to the rest of thedomain. Hence, the generalized eigenvalue problems is to find all eigen pairs( λ i,typek , ψ i,typek ) ∈ ( R + , V ( i ) ) such that A Ω i x = λ type ( R i B type R Ti ) x , i = 1 , . . . , N, (5)where B type ( · , · ), type ∈ { I , II } are the matrix representations of the followingbilinear forms. b I ( u, v ) := N (cid:88) i =1 (cid:90) Ω i α i ∇ u · ∇ v dx, u, v ∈ V h ,b II ( u, v ) := N (cid:88) i =1 (cid:32)(cid:90) Ω Bi α i,B ∇ u · ∇ v dx + (cid:90) Ω i \ Ω Bi α ∇ u · ∇ v dx (cid:33) , u, v ∈ V h . We also denote all eigenvalues of (5) by D type = diag( D type , D type , . . . , D typeN ) ∈ R N i × N i , type ∈ { I , II } , where D typei = diag( λ i,type , λ i,type , . . . , λ i,typen i ) ∈ R n i × n i , i = 1 , . . . , N, such that λ i,type ≥ λ i,type ≥ . . . ≥ λ i,typen i > . We define W typei := span { ψ i,typek } n i k =1 , type ∈ { I , II } , i = 1 , . . . , N as the space of the eigenfunctions associated with the eigenvalues λ i,typej . Wealso correspond to these spectral spaces the matrix forms denoted by W typei , i = 1 , . . . , N such that W type = diag( W type , W type , . . . , W typeN ) , type ∈ { I , II } . Further, we use the notation R type as the matrix representation of theoperator T ( type )0 which has the following structure. R type = R OO W type , R type ∈ R N Ω × ( N Ω + N i ) , type ∈ { I , II } . itle Suppressed Due to Excessive Length 11 B typeE , type ∈ { I , II } . The idea is based on expressing thesepreconditioners in terms of a combination of B typeC , type ∈ { I , II } and (cid:80) Ni =1 R Ti A − Ω i R i , where B typeC = ( R type ) T (cid:16) R type A N ( R type ) T (cid:17) − R type and A N = A (11) N A (12) N A (21) N A (22) N , A N ∈ R ( N Ω + N i ) × ( N Ω + N i ) such that A (11) N := square submatrix of A whose rows and columns are correspondingto the nodes belong to N c (cid:83) N m (cid:83) N i , A (12) N := rectangular submatrix of A whose rows are corresponding to allnodes belong to N c (cid:83) N m (cid:83) N i , meanwhile whose columns are correspondingto the nodes belong to N i , A (21) N = (cid:16) A (12) N (cid:17) T , A (22) N := square submatrix of A whose rows and columns are correspondingto the nodes belong to N i . Indeed, A (22) N = diag( A Ω , A Ω , . . . , A Ω N ).Hence R type A N ( R type ) T = R A (11) N R T R A (12) N ( W type ) T (cid:16) R A (12) N ( W type ) T (cid:17) T W type A (22) N ( W type ) T . Since W type A (22) N ( W type ) T = D type , type ∈ { I , II } and it is the nonsingularmatrix, we can easily compute B typeC . To this end and for the sake of simplicity,we first assume G type = R A (12) N ( W type ) T , S type = R A (11) N R T − G type (cid:0) D type (cid:1) − (cid:0) G type (cid:1) T . Hence, the inverse of 2 × R type A N ( R type ) T implies B typeC (11) = R T ( S type ) − R , B typeC (12) = − R T ( S type ) − G type ( D type ) − W type , B typeC (21) = − (cid:16) B typeC (12) (cid:17) T , B typeC (22) = ( W type ) T (cid:16) ( D type ) − + (cid:0) G type ( D type ) − (cid:1) T ( S type ) − G type ( D type ) − (cid:17) W type . Due to different sizes of the submatrices of B typeC , we can use the restrictionmatrices R i , i = 1 , . . . , N to obtain new matrices not only with identical sizes,but also with their previous performances. Let R r = diag( R , R , . . . , R N ), B type = B typeC (11) + R Tr B typeC (21) and B type = B typeC (12) R r + R Tr B typeC (22) R r . Hence, B typeE , type ∈ { I , II } as the enrichment preconditioners can be written explicitlyin the following matrix form. B typeE = B type + B type + N (cid:88) i =1 R Ti A − Ω i R i = B type + (cid:16) R Tr − B type A (12) N (cid:17) ( W type ) T ( D type ) − W type R r + N (cid:88) i =1 R Ti A − Ω i R i . To estimate an upper bound for the condition number of B typeE in the nextsection, we define the b typei ( · , · )-orthogonal projection operator π typei : V ( i ) → V ( i ) as π typei v h = N i (cid:88) k =1 b typei ( v h , ψ i,typek ) ψ i,typek , v h ∈ V h , type ∈ { I , II } . For any u h ∈ V h , we consider the function w h = u h − I u h ∈ V h with zerovalue both on ∂Ω i and on the rest of the domain Ω . In addition, we define I type : V h → V ( type )0 as I type u h = I u h + N (cid:88) i =1 π typei w h , u h ∈ V h , type ∈ { I , II } . This section is devoted to obtaining condition number estimates of the additiveaverage Schwarz preconditioners, which are defined in the previous section for type ∈ { I , II } . The proof is based on the standard additive Schwarz framework,where three assumptions to be held, see [32, p. 155]. Here, we need only toshow that there exists the stable splitting for all u ∈ V h as in Theorem 1. Notethat the other assumptions hold since there is no overlapping, and we also usethe exact bilinear forms. For simplicity, throughout this section, we use C as apositive constant, which is independent of the mesh sizes. We add a notationor an index to C if we want to emphasize a special constant. Furthermore,we use the notation (cid:46) to remove any constants except the mesh sizes in theinequalities.Let M typei be a given number such that 0 ≤ M typei < n i and λ i,typeM typei +1 <λ i,typeM typei . We define (cid:102) W typei := span { ψ i,typek } M typei k =1 , type ∈ { I , II } itle Suppressed Due to Excessive Length 13 and (cid:101) V ( type )0 = I V h + (cid:88) Ni =1 (cid:102) W typei , type ∈ { I , II } . For the analysis of the additive average Schwarz method, we define (cid:101) I type : V h → (cid:101) V ( type )0 as (cid:101) I type u h = I u h + (cid:88) Ni =1 (cid:101) π typei w h , u h ∈ V h , type ∈ { I , II } , where (cid:101) π typei v h = (cid:88) M typei k =1 b typei ( v h , ψ i,typek ) ψ i,typek , type ∈ { I , II } . To prove our main result, we need the following lemma.
Lemma 1.
For all u ∈ V ( i ) , i = 1 , . . . , N , | u − (cid:101) π I i u | H ( Ω i ) ,α ≤ Cλ i, I M I i +1 (cid:107) α / i ∇ u (cid:107) L ( Ω i ) , | u − (cid:101) π II i u | H ( Ω i ) ,α ≤ Cλ i, II M II i +1 (cid:16) (cid:107) α / iB ∇ u (cid:107) L ( Ω Bi ) + (cid:107) α / i ∇ u (cid:107) L ( Ω i \ Ω Bi ) (cid:17) , where | · | H ( Ω i ) ,α = a i ( · , · ) .Proof. We first express any u ∈ V ( i ) uniquely in terms of the eigenfunctions,i.e., u = (cid:80) n i k =1 b typei ( u, ψ i,typek ) ψ i,typek . Hence u − (cid:101) π typei u = (cid:88) n i k = M typei +1 b typei ( u, ψ i,typek ) ψ i,typek , type ∈ { I , II } and we have | u − (cid:101) π typei u | H ( Ω i ) ,α ≤ (cid:88) n i k = M typei +1 (cid:16) b typei ( u, ψ i,typek ) (cid:17) | ψ i,typek | H ( Ω i ) ,α . Furthermore, ( W typei ) T A Ω i W typei = D typei is equivalent to |∇ ψ i,typek | H ( Ω i ) ,α = λ i,typek , k = 1 , . . . , n i i = 1 , . . . , N, type ∈ { I , II } . Therefore | u − (cid:101) π typei u | H ( Ω i ) ,α ≤ (cid:88) n i k = M typei +1 (cid:16) b typei ( u, ψ i,typek ) (cid:17) λ i,typek , type ∈ { I , II } . (6)Using the Schwarz inequality for type = I, yields | u − (cid:101) π I i u | H ( Ω i ) ,α ≤ (cid:107) α / i ∇ u (cid:107) L ( Ω i ) (cid:88) n i k = M I i +1 (cid:107) α / i ∇ ψ i, I k (cid:107) L ( Ω i ) λ i, I k . Since ( W typei ) T ( R i B type R Ti ) W typei = I , type ∈ { I , II } and i = 1 , . . . , N we have equivalently (cid:107) α / i ∇ ψ i, I k (cid:107) L ( Ω i ) = 1 , (cid:107) α / iB ∇ ψ i, II k (cid:107) L ( Ω Bi ) + (cid:107) α / i ∇ ψ i, II k (cid:107) L ( Ω i \ Ω Bi ) = 1 . (7) Then | u − (cid:101) π I i u | H ( Ω i ) ,α ≤ Cλ i, I M I i +1 (cid:107) α / i ∇ u (cid:107) L ( Ω i ) . For type = II, from (6) we get | u − (cid:101) π II i u | H ( Ω i ) ,α ≤ C (cid:88) n i k = M II i +1 (cid:16) (cid:107) α / iB ∇ u (cid:107) L ( Ω Bi ) (cid:107) α / iB ∇ ψ i, II k (cid:107) L ( Ω Bi ) + (cid:107) α / i ∇ u (cid:107) L ( Ω i \ Ω Bi ) (cid:107) α / i ∇ ψ i, II k (cid:107) L ( Ω i \ Ω Bi ) (cid:17) λ i, II k . To complete the proof, it is sufficient to use (7).
Theorem 1.
For all u h ∈ V h the following results hold: a ( (cid:101) I type u h , (cid:101) I type u h ) (cid:46) max i λ i,typeM typei +1 Hh a ( u h , u h ) , type ∈ { I , II } , i = 1 , . . . , N. Proof.
We first consider the following splitting. u h = (cid:101) I type u h + N (cid:88) i =1 u i , where u i = (0 , . . . , , w i , , . . . , ∈ V ( i ) and w i = ( u h − (cid:101) I type u h ) | Ω i .Consequently, we get a ( (cid:101) I type u h , (cid:101) I type u h ) (cid:22) a ( u h , u h ) + a ( u h − (cid:101) I type u h , u h − (cid:101) I type u h ) , where a ( u h − (cid:101) I type u h , u h − (cid:101) I type u h ) = N (cid:88) i =1 a ( u i , u i ) = N (cid:88) i =1 a i ( u h − (cid:101) I type u h , u h − (cid:101) I type u h )= N (cid:88) i =1 a i ( w − (cid:101) π typei w, w − (cid:101) π typei w ) , type ∈ { I , II } . We first use Lemma 1 for the type = I. Then, we have a i ( w − (cid:101) π I i w, w − (cid:101) π I i w ) ≤ Cλ i, I M I i +1 α i ||∇ ( u h − I u h ) || L ( Ω i ) ≤ Cλ i, I M I i +1 α i (cid:18) Hh ||∇ u h || L ( Ω i ) + ||∇ I u h || L ( Ω i ) (cid:19) . Further, from [8, pp. 8-9] we get similarly α i ||∇ I u h || L ( Ω i ) ≤ C Hh a i ( u h , u h ) . (8)For type = I, the proof is completed by summing (8) over i = 1 , . . . , N .We now proceed to prove a similar result for type = II. Due to definition of I u hi , we use this fact that ∇ I u hi is equal to zero on each triangle τ (cid:54)∈ Ω Bi , i = 1 , . . . , N . Hence, we use Lemma 1 to get the first inequality as follow. itle Suppressed Due to Excessive Length 15 a i ( w − (cid:101) π II i w, w − (cid:101) π II i w ) ≤ Cλ i, II M II i +1 (cid:16) || α / i ∇ ( u h − I u h ) || L ( Ω i \ Ω Bi ) + α i,B ||∇ ( u h − I u h ) || L ( Ω Bi ) (cid:17) ≤ Cλ i, II M II i +1 (cid:16) || α / i ∇ u h || L ( Ω i \ Ω Bi ) + α i,B ||∇ ( u h − I u h ) || L ( Ω Bi ) (cid:17) ≤ Cλ i, II M II i +1 (cid:18) || α / i ∇ u h || L ( Ω i ) + α i,B (cid:16) ||∇ u h || L ( Ω Bi ) + ||∇ I u h || L ( Ω Bi ) (cid:17) (cid:19) ≤ Cλ i, II M II i +1 (cid:16) a i ( u h , u h ) + || α / i ∇ u h || L ( Ω i ) + α i,B ||∇ I u h || L ( Ω Bi ) (cid:17) ≤ Cλ i, II M II i +1 (cid:16) a i ( u h , u h ) + α i,B ||∇ I u h || L ( Ω Bi ) (cid:17) . Using the inverse inequality and the definition of the operator I , implies ||∇ I u h || L ( Ω Bi ) = (cid:88) τ ∈ Ω Bi ||∇ I u h || L ( τ ) = (cid:88) τ ∈ Ω Bi ||∇ ( I u h − u h ) || L ( τ ) ≤ C (cid:88) τ ∈ Ω Bi h − τ || I u h − u h || L ( τ ) ≤ C (cid:88) x ∈ ∂Ω ih ( u hi ( x ) − u hi ) . (9)In what follow, we need to use the following Poincar´e inequality [36]. || g || L ( Ω i ) ≤ C | g | H ( Ω i ) + C (cid:18)(cid:90) Ω i g dx (cid:19) , g ∈ H ( Ω i ) , (10)where C and C are two positive constants independent of the mesh size of Ω . Consider the function g = u hi − c , where c = Ω i ) (cid:82) Ω i u hi dx . Now, (10)implies || u hi − c || L ( Ω i ) ≤ C | u hi | H ( Ω i ) . (11)Furthermore, we use (4) and this property that all boundary elements ofsubdomains Ω i , i = 1 , . . . , N are quasi-uniform which implies that the numberof all nodes belonging to ∂Ω i , denoted by M i , is of the order H i /h i , to get M i (cid:16) u hi − c (cid:17) ≤ CM i µ i ( δ,γ ) (cid:88) Γ ⊂ ∂Ω i Γ || u hi − c || L ( Γ ) ≤ Ch − i || u hi − c || L ( ∂Ω i ) , (12)where Γ is γ m ( i ) and δ m ( i ) . Now, we use g and (12) to estimate an upper bound for (9) as follows. ||∇ I u h || L ( Ω Bi ) ≤ C (cid:88) x ∈ ∂Ω ih ( g ( x ) − g ) ≤ Ch − i || g || L ( ∂Ω i ) ≤ Ch − i H i || ˆ g || L ( ∂ ˆ Ω ) , where ˆ Ω is the reference element of unit diameter. From this inequality wededuce ||∇ I u h || L ( Ω Bi ) ≤ Ch − i H i (cid:16) | ˆ g | H ( ˆ Ω ) + || ˆ g || L ( ˆ Ω ) (cid:17) ≤ Ch − i (cid:16) H i | g | H ( Ω i ) + H − i || g || L ( Ω i ) (cid:17) ≤ C H i h i a i ( u h , u h ) , by the trace inequality, Theorem 3.1.2 in [10] and (11). To complete the proof,it suffices to take a summation over i = 1 , . . . , N .To estimate upper bounds of the condition numbers of B typeE A , type ∈{ I , II } it suffices to use Lemma 3 in [32, pp. 156-158] and Theorem 1. Theorem 2.
The condition numbers of the enriched additive average Schwarzpreconditioners are bounded by κ ( B typeE A ) ≤ C (cid:18) Hh (cid:19) max i λ i,typeM typei +1 , type ∈ { I , II } , i = 1 , . . . , N, where C is a positive constant. In this section, numerical results confirm the additive average Schwarzmethod’s validity and efficiency with adaptive enrichment, where the jumpsin the coefficient α ( x ) in (1) are very large and even change rapidly. Thosejumps might be occurred inside of the subdomains, or even across thesubdomain boundaries. To have a complicated distribution of jumps, we usethe following pattern [26], depicted in Figure 4, which is a periodic patternwhen the number of subdomains are increased (cf. Figure 5 ). We alsoconsider background channels, crossing channels and corner channels denotedby α b , α i and α c , respectively. For different values of α b , α c and α i , our testproblem has the right-hand side function f ( x, y ) = 2 π sin( πx ) sin( πy )defined in Ω = [0 , × [0 , H ∈ { / , / } , h ∈ { / , / } , and h ∗ ∈ { / , / } (cf. Figure 5). We use the similar itle Suppressed Due to Excessive Length 17 Fig. 4
The decomposition of Ω = [0 , × [0 ,
1] into 3 × α ( x ) consisting of α b , α c and α i in the white, red and green areas,respectively. Further, all mortar sides are denoted by the thick and black line segments. Fig. 5
The distribution of all jumps in α ( x ) following the extended pattern in Figure 4,where Ω = [0 , × [0 ,
1] is divided into 6 × H = 1 / h = 1 /
36 and h ∗ = 1 /
54. Further, all mortar sides are the coarsemortar. locations for α c and α i as the periodic patterns and the similar mortar sides,as depicted in Figure 4, for N ≥
9. For instance, see Figure 5, where N = 36.Further, in what follows, we use two terms, i.e., coarse mortar and finemortar to distinguish between the coarse and fine triangulations connected to the mortar sides. To withdraw numerical results, we use the additiveaverage Schwarz method to produce the enrichment preconditioners. Weestimate the condition number of such preconditioners for type ∈ { I , II } .Moreover, the iteration numbers in all tables come from the preconditionedconjugate gradient method, with the tolerance 5 e −
6, based on the producedpreconditioners.
ADD
Different Values for αb , αc and αi
1, 1e3, 1e4 1, 1e4, 1e6
Coarse Mortar Fine Mortar Coarse Mortar Fine Mortar N u m b e r o f Subd o m a i n s × . . . . × . . . . Table 1
The condition number of B typeE A and the number of iterations of thepreconditioned conjugate gradient method (in parentheses) for type = II with differentvalues for α b , α c , and α i . For 6 × × h and h ∗ belong to { / , / } and { / , / } , respectively. In addition, selecting the number of eigenfunctions for eachsubdomain to construct the enrichment coarse space is based on the adaptive enrichment,where the threshold is 50. In Table 1, we set different values for α b , α c , and α i for different numberof subdomains, for instance, see Figure 5 including 6 × type = II, where the threshold is 50,i.e., all eigenfunctions associated with all eigenvalues larger than 50 areconsidered. Note that the condition number of the non-enrichmentpreconditioner is very large. For instance, in the case 6 × . H/h andindependent of the number of subdomains. We also observe that if weconsider the same number of eigenfunctions for each subdomain based on aproper threshold, consequently, the condition number estimates areindependent of values of α b , α c and α i .Table 2 demonstrates only the condition number estimates and iterationnumbers of the preconditioned conjugate gradient method (in parentheses) inthe cases type = II, coarse mortar, 6 × × α b = 1, α c = 1e4, and α i = 1e6. For the other cases and also different values of jumpsin the coefficient α ( x ), we have similar results. Also, the construction of theenrichment coarse space is relied on the fixed number of eigenfunctions foreach subdomain varying from 0 to 7.It is conclusive that both approaches to enrich the standard coarse space,i.e., using the given threshold and imposing the fixed numbers of itle Suppressed Due to Excessive Length 19ADD Number of Eigenfunctions in Each Subdomain N u m b e r o f Subd o m a i n s × . . . . . . . . × . . . . . . Table 2
The implementation of the additive average Schwarz method, type = II, with thefixed number of eigenfunctions for each subdomain to estimate the condition number of B typeE A and the number of iterations of the preconditioned conjugate gradient method (inparentheses) for 6 × h = 1 / h ∗ = 1 /
54) and 9 × h = 1 / h ∗ = 1 /
81) subdomainsand the coarse mortar case. Further, the distribution of jumps in α ( x ) are α b = 1, α c = 1e4and α i = 1e6. ADD Different Values for αb , αc and αi
1, 1e3, 1e4 1, 1e4, 1e6
Coarse Mortar Fine Mortar Coarse Mortar Fine Mortar N u m b e r o f Subd o m a i n s × t y p e I 625 525 635 534II 75 80 75 809 × t y p e I 1988 1680 1992 1681II 240 233 240 233
Table 3
The total numbers of the eigenfunctions associated with the eigenvalues greaterthan 50 to enrich the coarse spaces used in the additive average Schwarz method for bothtypes I and II. eigenfunctions for each subdomain, lead to similar results. This fact can beviewed by comparing, for instance, the third column of Table 1 with the lastcolumn of Table 2.Table 3 gives the total number of required eigenfunctions for the adaptiveenrichment coarse spaces in the cases 6 × × type = Iand II and different values for α ( x ). As we can see from this table, solving thesecond type of the generalized eigenvalue problem (5) is more efficient thansolving the first type. To analyze the distribution of eigenfunctions, Figure 6contains the polar histograms for both coarse and fine mortars in the case 6 × type ∈ { I , II } . It clearly shows that we need only considera few eigenfunctions for each subdomain in the case type = II compared to type = I. a ) Fine mortar and type = I . ( c ) Coarse mortar and type = I . ( d ) Coarse mortar and type = II . ( b ) Fine mortar and type = II . Fig. 6
The histograms of the number of eigenfunctions associated with the eigenvaluesgreater than 50 corresponding to the partition of Ω = [0 , × [0 ,
1] into 6 × α b = 1, α c = 1e4 and α i = 1e6. The number of subdomains are ordered anticlockwisearound the polar histograms from 1 to 36. For type = I and II, the number of eigenfunctionsare grouped into the ranges 0 to 70 and 0 to 10, respectively. The largest numbers of theeigenfunctions in ( a ) − ( d ) are 39, 6, 63 and 5, respectively.itle Suppressed Due to Excessive Length 21 In this paper, we have employed the additive average Schwarz method withthe enrichment coarse spaces to solve second order elliptic boundary problemcoupled with very large jumps in the function α ( x ), where the finite elementdiscretization of that problem has been based on the nonmatchingtriangulations across the subdomains interfaces. We have proved that thecondition number estimates for the produced preconditioners by the additiveaverage Schwarz methods for type ∈ { I , II } are proportional to the ratio H/h and independent of the number of subdomains. Besides, we have comparedthe numerical results to conclude that, in practice, type = II has much betterperformance than type = I.
Acknowledgements
The work of Leszek Marcinkowski was partially supported by Polish ScientificGrant: National Science Center: 2016/21/B/ST1/00350.
References
1. T. Arbogast, L. C. Cowsar, M. F. Wheeler and I. Yotov, Mixed finite element methodson nonmatching multiblock grids,
SIAM J. Numer. Anal.
Multiscale Model. Simul.
SIAM J. Numer. Anal.
Math. Comp.
65, 897–921, (1996)6. M. Brezina, C. Heberton, J. Mandel and P. Vanˇek, An iterative method with convergencerate chosen a priori, Tech. Report 140, Center for Computational Mathematics, Universityof Colorado Denver, Denver, (1999)7. P. E. Bjørstad, M. Dryja and E. Vainikko, Parallel implementation of a Schwarz domaindecomposition algorithms, in Applied Parallel Computing in Industrial Problems andOptimization (J. Wasniewski, J. Dongara, K. Madsen, and D. Olsem, eds.) Lecture Notesin Computer Science, Springer, 1184, 141–157, (1996)8. P. E. Bjørstad, M. Dryja and T. Rahman, Additive Schwarz methods for elliptic mortarfinite element problems,
Numer. Math.
95, 427–457, (2003)9. T. F. Chan and T. P. Mathew, Domain decomposition algorithms,
Acta Numerica.
J. Numer. Math.
Comput. Methods Appl. Math . 10, 164–176, (2010)15. Y. Efendiev, J. Galvis and P. S. Vassilevski, Spectral element agglomerate algebraicmultigrid methods for elliptic problems with high contrast coefficients. In: Y. Huang, R.Kornhuber, O. B. Widlund, J. Xu (eds.) Domain Decomposition Methods in Science andEngineering XIX, Lecture Notes in Computational Science and Engineering, vol. 78, pp.407–414. Springer, Berlin (2011)16. E. Eikeland, L. Marcinkowski and T. Rahman, Overlapping Schwarz methods withadaptive coarse spaces for multiscale problems in 3D,
Numer. Math.
J. Numer. Math.
10, 109–125, (2002)18. M. J. Gander, A. Loneland and T. Rahman, Analysis of a new harmonicallyenriched multiscale coarse space for domain decomposition methods, arXiv preprintarXiv:1512.05285 (2015)19. A. Heinlein, A. Klawonn, J. Knepper and O. Rheinbach, Adaptive GDSW coarse spacesfor overlapping Schwarz methods in three dimensions,
SIAM J. Sci. Comput.
MultiscaleModel. Simul.
13, 571–593, (2015)21. H. H. Kim, E. Chung and J. Wang, BDDC and FETI-DP preconditioners with adaptivecoarse spaces for three-dimensional elliptic problems with oscillatory and high contrastcoefficients,
J. Comput. Phys.
SIAM J. Numer. Anal.
Numer. Math.
Math.Comput. Simulation.
82, 1812–1831, (2012)26. L. Marcinkowski and T. Rahman, Additive average Schwarz with adaptive coarse spaces:scalable algorithms for multiscale problems,
Electron. Trans. Numer. Anal.
49, 28–40.,(2018)27. A. M. Matsokin and S. V. Nepomnyaschikh, A Schwarz alternating method in asubspace.
Sov. Math.
Numer. Math.
Numer. Math. hp convergence results for the mortar finite elementmethod, Math. Comp.
69, 521–546, (2000)32. B. F. Smith, P. E. Bjørstad and W. D. Gropp, Domain Decomposition: ParallelMultilevel Methods for Elliptic Partial Differential Equations, Cambridge University Press,(1996).33. N. Spillane and D. J. Rixen, Automatic spectral coarse spaces for robust finite elementtearing and interconnecting and balanced domain decomposition algorithms,
Int. J.Numer. Methods. Eng.
Numer. Math.
J. Comput. Appl. Math.
34, 93–117, (1991)36. A. Toselli and O. B. Widlund, Domain decomposition methods-algorithms and theory,volume 34 of Springer Series in Computational Mathematics, Springer-Verlag, Berlin,(2005)37. J. Xu and J. Zou, Some nonoverlapping domain decomposition methods,