A domain decomposition method derived from the Primal Hybrid Formulations for 2nd order elliptic problems

XX Congreso de Ecuaciones Diferenciales y Aplicaciones X Congreso de Matem´ atica Aplicada Sevilla, 24-28 septiembre 2007 (pp. 1–8) A domain decomposition method derived from the Primal Hybrid Formulations for 2nd order elliptic problems C. Bernardi1, T . Chac´ on Rebollo2, E. Chac´ on Vera2 1Laboratoire Jacques-Louis Lions, Universit´e Paris VI et CNRS, 175 rue du Chevaleret 75013 Paris. E-mail: [email protected] 2Dpto. E.D.A.N., Universidad de Sevilla, Aptdo. 1160, E-41080 Sevilla. E-mails: [email protected], [email protected] Palabras clave: Domain Decomposition Methods, Lagrange Multipliers Resumen We consider the primal hybrid formulation for second order elliptic problems introduced by Raviart-Thomas [9] and apply the classical iterative method of Uzawa to obtain a non overlapping domain decomposition method that converges geometrically with a mesh independent ratio. The proposed method connects with the Finite Element Tearing and Interconnecting (FETI) method proposed by Farhat-Roux and collaborators [7]-[8]. In this research work we use the detailed work on domains with corners developed by Grisvard [6], which clarifies the situation of cross-points, and the direct computation of the duality H−1/2−H1/2using the H1/2scalar product; therefore no consistency error appears. 1. Introduction The primal hybrid formulation for second order elliptic problems inforce via the Lagrange multipliers the continuity of the approximations across interfaces, and this is expressed via the duality H−1/2−H1/2, see Raviart-Thomas [9]. Usually, for numerical discretizations, this duality is worked out by means of some projection operator onto the L2space on the interfaces, see Ben Belgacem [3]. In our approach we use Riesz representation and replace the duality with the H1/2scalar product that is explicitly computed. As a consequence, we have a formulation in terms of a saddle point problem suitable for iterative techniques, see Bacuta [2]. The method presented is similar to the classical Lagrange Finite Element Tearing and Interconnecting (FETI) method proposed by Farhat-Roux and collaborators [7]-[8]. In this research work we use the detailed work on domains with corners developed by Grisvard 1 C. Bernardi, T . Chac´on Rebollo, E. Chac´on Vera [6] which clarifies the situation of cross-points and the direct computation of the duality H−1/2−H1/2using the H1/2scalar product; therefore no consistency error appears. 2. Formulation with Lagrange Multipliers Let Ω ⊂R2be a poligonal domain and consider f∈L2(Ω). Then, our departure problem looks for u∈H1 0(Ω) such that (∇u, ∇v)Ω+ (u, v)Ω= (f, v)Ω∀v∈H1 0(Ω).(1) Assume that Ω is a polygonal bounded domain in R2with a Lipschitz–continuous boundary and consider a decomposition without overlapping in polygonal subdomains Ω = ∪R r=1Ωrand Ωr∩Ωr0=∅,1≤r < r0≤R(2) where each Ωrhas a Lipschitz–continuous boundary. We describe ∂Ωrin terms of its edges via ∂Ωr= Γr,0∪Γr,1∪... ∪Γr,Jr(3) where Γr,0=∂Ωr∩∂Ω such that ∂Ω = ∪R r=1Γr,0and assume that Ωr∩Ωsis either empty, a single point or a full edge Γr,s. On each Γr,j we consider the classical Hilbert space of traces H1/2 00 (Γr,j) and its dual space H−1/2 00 (Γr,j), see Adams [1]. We call skeleton of Ω, and denote it by E, the set of all interfaces in Ω E=∪I i=1Γi(4) where Γi= Γi,0for i= 1, ..., R describe the boundary ∂Ω, and for i≥R+ 1 we set Γi= Γr,j for some r, j ≥1. Green’s formulae on polygonal domains will be used Lemma 1 (Grisvard [6]) When O ⊂ R2is a poligonal domain and ∂O=∪J j=1Γj, then H2(O)is dense on E={u∈H1(O); ∆u∈L2(Ω)}. The mapping u7→ ∂nju|Γjhas a unique continuous extension from Eto H−1/2 00 (Γj)dual space of H1/2 00 (Γj). Moreover, for each u∈Eand v∈H1(O)such that v|Γj∈H1/2 00 (Γj)we have −(∆u, v)O= (∇u, ∇v)O− J X j=1 < ∂nju, v >−1/2,00,Γj.(5) Also, using that D(O)dis dense on H(div;O), for any ~q ∈H(div;O), we have nj·~q ∈ H−1/2 00 (Γj)and for any v∈H1(Ω) with v|Γj∈H1/2 00 (Γj) (~q, ∇v)O+ (div(~q), v)O= J X j=1 <nj·~q, v >−1/2,00,Γj.(6) 2 A domain decomposition method derived from the Primal Hybrid Formulations Next, on each Ωrwe consider the classical Hilbert space H1 b(Ωr) = {vr∈H1(Ωr); vr= 0 on ∂Ωr∩∂Ω}(7) with scalar product (ur, vr)1,Ωr= (ur, vr)Ωr+ (∇ur,∇vr)Ωrthe dense subspace Wrof H1 b(Ωr) given by Wr={u∈H1 b(Ωr); u|Γr,j ∈H1/2 00 (Γr,j), j = 1, ..., Kr}, needed to apply Green’s formulae, the global, defined on Ω, Hilbert space X={v∈L2(Ω); vr=v|Ωr∈H1 b(Ωr), r = 1, ..., R} ≈ R Y r=1 H1 b(Ωr) (8) with scalar product and norm given by (u, v)X= R X r=1 (ur, vr)1,Ωr,kvk2 X= (v, v)X= R X r=1 kvrk2 1,Ωr,∀v, u ∈X(9) and finally X0≈QR r=1 Wrthat is also a dense subspace of X. We also need the Hilbert space H0(div; Ω) = {~q ∈L2(Ω)d;div(~q)∈L2(Ω),nr,0·~q = 0 en Γr,0,1≤r≤R}(10) where nr,j ·~q ∈H−1/2 00 (Γr,j), j = 1, ..., Kr, consider Tr=QKr j=1 H−1/2 00 (Γr,j) and Mgiven by M={~µ ∈ R Y r=1 Tr;µr,j =nr,j ·~q, for some ~q ∈H0(div; Ω)}.(11) Let b:M×X7→ Rgiven for v∈X0,~ λ∈Mby b(~ λ, v) = I X i=R+1 < λi, vs−vt>−1/2,00,Γi(12) when Ωs∩Ωt= Γiand extended by density to all v∈X. Then H1 0(Ω) = {v∈X;b(~ λ, v) = 0,∀~ λ∈M}. Define the bilinear form a:X×X7→ Rgiven by a(u, v) = R X r=1 {(∇ur,∇vr)Ωr+ (ur, vr)Ωr}= R X r=1 ZΩr {∇ur· ∇vr+urvr}dx. (13) Then, the primal hybrid formulation for Poisson problem (1) consists in looking for a pair (u,~ λ)∈X×Msuch that a(u, v) + b(~ λ, v) = R X r=1 (f, vr)Ωr,∀v∈X(v|Ωr=vr) (14) b(~µ, u)=0,∀~µ ∈M. (15) 3 C. Bernardi, T . Chac´on Rebollo, E. Chac´on Vera Theorem 1 If u∈H1 0(Ω) solves the Dirichlet problem (1) then there exists a unique (u,~ λ)∈X×Mthat solves problem (14)-(15). If (u,~ λ)∈X×Msolves (14)-(15) then u∈H1 0(Ω) and solves the Dirichlet problem (1). Moreover, for i=R+ 1, ..., I λi=−∂niu∈H−1/2 00 (Γi).(16) Next, via Riesz representation we identify H−1/2 00 (Γi) (dual space of H1/2 00 (Γi)) with H1/2 00 (Γi), write the duality in terms of the scalar product in H1/2 00 (Γi), identify Mwith its dual space M0and define b:M×X7→ Rgiven for any v∈X0,~ λ∈Mby b(~ λ, v) = I X i=R+1 (λi, vs−vt)1/2,00,Γi(17) when Ωs∩Ωt= Γiand extended by density to all v∈X. Then, the formulation of Poisson problem (1) that we shall use is: Find a pair (u,~ λ)∈X×Msuch that a(u, v) + b(v,~ λ) = R X r=1 (f, vr)Ωr,∀v∈X, (18) b(u, ~µ) = 0,∀~µ ∈M. (19) Thanks to Theorem 1 this is equivalent to (1) but it also is whithin the saddle point problems framework, see Girault-Raviart [5], which allows the use of different methods for computing the solution. Also, the analysis at the continuous level is reproduced in the discrete version of the saddle point problems as a simple consequence of the finite element extension theorems, see for instance Bernardi-Maday-Rapetti [4]. 3. Domain decomposition methods A rephrasing of the problem in terms of functional operators will clarify what we do. To fix ideas we work with Ω split up slicewise into two subdomains. Let B:X7→ Mgiven by Bv = (v1)|Γ−(v2)|Γ, i.e., the jump of vacross the interface Γ, set R:X07→ Xas the Riesz isomorphism associated with the scalar product a(·,·) on Xand F:X7→ Rgiven by < F, v >=P2 r=1(f, vr)Ωr.Then, our saddle point problem looks for (u, λ)∈X×M such that R−1u+B0λ=Fon X0(20) B u = 0 on M, (21) where B0is the transpose operator to B. Then, u=R(F−B0λ)⇒Bu =BRF −BRB0λ and (using Bu = 0) from here we have the dual problem associated to the saddle point problem (BRB0)λ=BRF on M. (22) Thanks to the inf-sup condition the operator BRB0is symmetric positive definite, see Bacuta [2]. Now the resolution of (22) via an iterative method is possible; we propose 4 A domain decomposition method derived from the Primal Hybrid Formulations the use of the iterative method of Richardson, which amounts to the classical Uzawa’s Method. Given ρ > 0 and λ0∈M, for m= 0,1,2,3, ... set rm=BRF −(BRB0)λm=B um,using (??) (23) λm+1 =λm+ρ rm(24) which unfolds from (20)-(21) as 2 X r=1 (um,r, vr)1,Ωr= 2 X r=1 (f, vr)Ωr−(λm, v1−v2,)1/2,00,Γ,∀v∈X, (25) and update λm+1 =λm+ρ(um,1−um,2).(26) Following standard convergence results, see Bacuta [2] and references therein, we have geometric convergence for this iterative process by simply blocking to zero the test functions alternatively on each subdomain Theorem 2 The iterative process: Given ρ > 0and λ0∈M, find for m≥0um∈Xvia (∇um,1,∇v1)Ω1+ (um,1, v1)Ω1= (f, v1)Ω1−(λm, v1)1/2,00,Γ,∀v1∈X1, (∇wm,2,∇v2)Ω2+ (um,2, v2)Ω2= (f, v2)Ω2+ (λm, v2)1/2,00,Γ,∀v2∈X2, and update λm+1 =λm+ρ(um,1−um,2)on Γ is a non overlapping domain decomposition method geometrically convergent with a ratio of convergence independent of the mesh size. The drawback that this method presents is how to fix the optimal parameter ρ > 0. In the numerical experiments that we present the value of ρhas been tuned easily by hand thanks to the great speed of convergence that the method exhibits. For a method that has no need of fixing any parameter we could use the application of the Conjugate Gradient Method which is the core of the FETI methods. 4. Numerical experiments We compute on a non convex domain with three subdomains. We use a Galerkin approximation with P1Lagrange finite elements on a uniform triangular mesh size hof Ω and its restriction to each of the Ωifor i= 1,2,3. The numerical results show a geometric rate of convergence with a mesh independent ratio as the theory predicts. We set Ω = (−1,1)2\ {(−1,0) ×(−1,0)}and decompose it into three squares so that our interfaces are Γ1={0} × (0,1) and Γ2= (0,1) × {0}. Then we solve −∆u= 1,on Ω, u = 0 on ∂Ω. We take λ0= (0,0) and stop iterating for m=niter(h) such that kum+1 h−um hkX kum hkX =(P3 i=1 RΩi|∇(um+1 i,h −um i,h)|2dx)1/2 (P3 i=1 RΩi|∇um i,h|2dx)1/2≤10−7; (27) 5 C. Bernardi, T . Chac´on Rebollo, E. Chac´on Vera 0 2 4 6 8 10 12 14 16 18 20 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 Error ratio 0 2 4 6 8 10 12 14 16 18 0 0.5 1 1.5 2 Uzawa h=1/4 Uzawa h=1/8 Uzawa h=1/16 Uzawa h=1/32 CG h=1/4 CG h=1/8 CG h=1/16 CG h=1/32 Figura 1: Decrease error ratio given by (29) on the L-shaped domain test using Uzawa’s Method and Conjugate Gradient Method (CG). we compute the errors and their decrease ratio given for m≥0 by euh(h, m) = kuh−um hkX= ( 3 X i=1 ZΩi |∇(uh−um i,h)|2dx)1/2(28) r(h, m) = euh(h, m + 1) euh(h, m), r(h)≈l´ım mr(h, m) (29) where uhis the P1solution computed on the whole domain Ω. For Uzawa’s Method we found by performing few several tests that ρ≈0.12 seems to be the closest value to the optimal one. The results are shown in Table 1. 1/h 4 8 16 32 #iterations 17 17 17 21 ρ≈0.12 ≈0.12 ≈0.12 ≈0.12 r(h)≈0.43.. ≈0.43... ≈0.43... ≈0.44... Table 1: Uzawa’s Method: Number of iterations, values of ρand of r(h) obtained with the inverted L-shape domain for different values of h. We also computed the solution using the Conjugate Gradient Method. In Figure 1, For both of the iterative methods proposed, we show the ratio between consecutive errors. Figure 2 shows the Galerkin solution computed in Ω and Figure 3 shows the solution computed right after the first iteration of Uzawa’s Method. We see the lack of jumps on the two interfaces. Acknowledgments Research work partially supported by Ministerio de Educacion y Ciencia with public funds: proyect MTM2006-01275. 6 A domain decomposition method derived from the Primal Hybrid Formulations 0 .0072936 0.021881 0.036468 0.051055 0.065643 0.08023 0.094817 0.1094 0.12399 0.13858 Iso Surfaces of Aproximated Solution with h=1/8. GiD x y z Figura 2: Approximated solution computed with standard Galerkin P1finite elements on the whole domain and with h= 1/8. 7 C. Bernardi, T . Chac´on Rebollo, E. Chac´on Vera 0 .0073365 0.022009 0.036682 0.051355 0.066028 0.080701 0.095374 0.11005 0.12472 0.13939 Iso Surfaces of Aproximated Solution with ddm. for h=1/8. GiD x y z Figura 3: Computed solution after the first iteration with h= 1/8using Uzawa’s Method. Referencias [1] Adams, R. A., Sobolev Spaces. Pure and Applied Mathematics, Vol. 65. Academic Press, New YorkLondon, 1975. [2] Bacuta, C., A unified approach for Uzawa algorithm SIAM J. Numer. Anal., Vol. 44, No. 6, pp.-26332649, 2006. [3] Ben Belgacem, F., The Mortar finite element method with Lagrange multipliers Numerische Mathematik, 84:173-197, 1999. [4] C. Bernardi, Y. Maday and F. Rapetti, Discr´etisations variationnelles de probl`emes aux limites elliptiques. Math´ematiques & Applications, 45. Springer-Verlag, New-York, 2002. [5] Girault, V. and Raviart, P.-A., Finite Element Methods for Navier-Stokes Equations. Theory and Algorithms. Springer Series in Computational Mathematics, 5. Springer-Verlag, Berlin, 1986. [6] Grisvard, P., Singularities in Boundary value problems. Recherches en Math´ematiques Appliqu´ees 22, Masson, 1992. [7] C. Farhat, F.-X. Roux, A method of finite element tearing and interconnecting and its parallel solution algorithm, Int. J. Numer. Methods Engrg. 32 (1991) 1205-1227. [8] C. Farhat, F.-X. Roux,An unconventional domain decomposition method for an efficient parallel solution of large-scale finite element systems SIAM J. Sci. Statist. Comput., 13, 379-396 (1992) [9] Raviart, P.A. and Thomas, J.-M., Primal Hybrid Finite Element Methods for second order eliptic equations Math Comp, Vol 31, number 138, pp.- 391-413, 1977. 8