Preconditioning of Elliptic Saddle Point Systems by Substructuring and a Penalty Approach th International Conference on Domain 16 th International Conference on Domain 16 Decomposition Methods Decomposition Methods January 12-15, 2005 Clark R. Dohrmann Structural Dynamics Research Department Sandia National Laboratories Albuquerque, New Mexico Sandia is a multiprogram laboratory operated by Sandia Corporation, a Lockheed Martin Company, for the United States Department of Energy’s National Nuclear Security Administration under contract DE-AC04-94AL85000.
Thanks • DD16 scientific & organizing committees • Sandia – Rich Lehoucq – Pavel Bochev – Kendall Pierson – Garth Reese • University – Jan Mandel
Overview • Saddle Point Preconditioner – related penalized problem – positive real eigenvalues – conjugate gradients • BDDC Preconditioner – overview – modified constraints • Examples – H(grad), H(div), H(curl)
Motivation original system: ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ T u f A B A > 0 on kernel of B , C ≥ 0 ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ = ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ A T = A , C T = C , B full rank − p g B C ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ penalized (regularized) system: Axelsson (1979) ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ T u f A B ~ ~ ~ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ = > = C 0 , C C T ~ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ − p g B C ⎝ ⎠ ⎝ ⎠ ⎝ ⎠
exact solution of penalized system ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ z r T A B ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ u u = ~ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ − z r B C ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ p p is ~ = − z C ( Bz r ) − 1 ~ − p u p = + T 1 S A B C B ~ A = + z S ( r B C r ) − − 1 T 1 u A u p primal rather than dual Schur complement considered ~ − = − + 1 T S ( C BA B ) ~ C
Some Connections penalty method for C = 0 and g = 0: = − + ρ G u Au / 2 u f ( Bu ) ( Bu ) / 2 T T T ⇒ + ρ = ( A B B ) u f T ~ C = ρ matrix same as S A for I / = + augmented Lagrangian: L G p Bu T + ρ ⎛ A B B B ⎞ ⎛ u ⎞ ⎛ f ⎞ T T ⇒ ⎜ ⎟ ⎜ ⎟ = ⎜ ⎟ B 0 p 0 ⎝ ⎠ ⎝ ⎠ ⎝ ⎠
regularized constraint preconditioning: Benzi survey (2005) ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ z r T A B ⎜ ⎟ ⎜ u ⎟ ⎜ u ⎟ = ~ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ z r − B C ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ p p ~ 1 = − z C ( Bz r ) − p u p regularized constraint equation satisfied exactly
Penalty Preconditioner recall exact solution of penalized system ~ = − + − 1 T 1 z S ( r B C r ) ~ − u A u p = + T 1 S A B C B ~ A = − − 1 z C ( Bz r ) p u p ⇒ preconditioner: ~ = − + − 1 T 1 z M ( r B C r ) u u p ~ M preconditioner for S A = − − 1 z C ( Bz r ) p u p
matrix representation: ~ ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ z − r ⎛ ⎞ 1 − I 0 T 1 M 0 I B C ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ u ⎜ ⎟ u ⎜ ⎟ = ~ ~ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ − − z 1 − r − 1 − C B I 0 C 0 I ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ p p � � � � � � � � � � � � � � � � � − 1 M − M 1 question: what about spectrum of A ⎛ ⎞ T A B ⎜ ⎟ = A ⎜ ⎟ − B C ⎝ ⎠
eigenproblem: = λ z z A M introduce: Bramble & Pasciak (1988), Klawonn (1998) − ⎛ ⎞ S A M 0 = ⎜ ⎟ H ~ ⎜ ⎟ − 0 C C ⎝ ⎠ eigenvalues same as those for = λ z z HM − A H 1 � � � � � symmetric can show H > 0 ⇒ HM -1 A > 0
CG Connection (B&P too) original linear system: A w = b equivalent linear system: HM -1 A w = HM -1 b ⇒ solve using pcg with H as preconditioner eigenvalues of preconditioned system same as those of M -1 A − ⎛ ⎞ S A M 0 = ⎜ ⎟ H not available, but not a problem H ~ ⎜ ⎟ − 0 C C ⎝ ⎠
required calculations: ~ ⎛ ⎞ ⎛ ⎞ + d M ( a B C a ) − − 1 T 1 ⎜ ⎟ ⎜ ⎟ = = M a u − 1 u p ~ ⎜ ⎟ ⎜ ⎟ − d C ( Bd a ) − 1 ⎝ ⎠ ⎝ ⎠ p u p ~ ⎛ ⎞ − + S d ( a B C a ) − T 1 ⎜ ⎟ = a HM − 1 A u u p ⎜ ⎟ − − Bd a Cd ⎝ ⎠ u p p recurrences: = + α w w p A − k k 1 k k = − α r r p A − k k 1 k k = − α z z p M A − 1 − k k 1 k k ~ ~ = − α r r p HM A − 1 − k k 1 k k
Theory Theorem (C = 0): Given α 1 > 1, γ 1 > 0, and α ≤ ≤ α ∀ ∈ ℜ u Mu u S u u Mu u T T T n 1 A 2 ~ γ ≤ ≤ γ ∀ ∈ ℜ p BM B p p C p p BM B p p − − T 1 T T T 1 T m 1 2 eigenvalues satisfy − A δ ≤ σ ≤ δ ( ) M 1 1 2 where { } δ = σ α α σ α γ min ( / ), /( ) 1 2 1 2 1 2 2 { } δ = α − σ − σ α γ max 2 , ( 2 / ) / 2 2 2 1 2 1 and σ 1 > 0, σ 2 > 0, σ 1 + σ 2 = 1
~ = ε C C Simplification for 0 recall with α 1 > 1 α ≤ ≤ α ∀ ∈ ℜ u Mu u S u u Mu u T T T n 1 A 2 ~ γ ≤ ≤ γ ∀ ∈ ℜ p BM B p p C p p BM B p p − − T 1 T T T 1 T m 1 2 ~ → as ε → 0 BS B C − 1 T A ⇒ γ 1 bounded below by 1/ α 2 and γ 2 bounded above by 1/ α 1 ⇒ − A ( ) α α ≤ σ ≤ α / / 2 ( ) 3 / 2 M 1 1 2 2
two goals: 1. Ensure α 1 > 1 (i.e. S A – M > 0) 2. Minimize α 2 / α 1 ( M good preconditioner for S A ) potential issues: 1. Scaling preconditioner M to satisfy Goal 1 2. S A becomes very poorly conditioned as ε → 0 conclusion: effective preconditioner for S A is essential
Recap of Penalty Preconditioner • Based on approximate solution of related penalized problem • Conjugate gradients can be used to solve saddle point system • Theory provides conditions for scalability • Effectiveness hinges on preconditioner for primal Schur complement S A Ref: C. R. Dohrmann and R. B. Lehoucq, “A primal based penalty preconditioner for elliptic saddle point systems ,” Sandia National Laboratories, Technical report SAND 2004-5964J.
BDDC Overview • Primal Substructuring Preconditioner – no coarse triangulation needed – recent theory by Mandel et al M7 – multilevel extensions straightforward • Numerical Properties – C (1+log(H/h)) 2 condition number bounds – preconditioned eigenvalues all ≥ 1 (woo hoo!) • eigenvalues identical to FETI-DP – local and global problems only require sparse solver for definite systems
FE discretization of pde elements partitioned into substructures solve for unknowns on sub boundaries (Schur complement)
Building Blocks Additive Schwarz: – Coarse grid correction v 1 − = + + 1 M r v v v – Substructure correction v 2 1 2 3 – Static condensation correction v 3 − 1 FETI-DP counterparts: ↔ λ * T k v F K F 1 I cc I rc rc N ∑ s − ↔ 1 T ↔ λ v s s s k v B K B Dirichlet preconditioner 3 2 r rr r = s 1
Coarse Grid and Substructure Problems coarse problem: Φ ⎛ ⎞ ⎛ ⎞ K Q ⎛ 0 ⎞ T ⎜ ⎟ ⎜ ⎟ = ⎜ ⎟ i i i ⎜ ⎟ ⎜ ⎟ Λ Q 0 I ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ i i ⇒ = Φ Φ T K K ci i i i • K ci coarse element matrix • assemble K ci ⇒ K c • K c positive definite substructure problem: ⎛ ⎞ ⎛ ⎞ K Q v ⎛ r ⎞ T ⎜ ⎟ ⎜ ⎟ = ⎜ ⎟ i i 2 i i ⎜ ⎟ ⎜ ⎟ λ Q 0 0 ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ i 2 i
Mixed Formulation of Linear Elasticity ⎛ ⎞ ⎛ ⎞ ⎛ ⎞ A B u f T ⎜ ⎟ ⎜ ⎟ = ⎜ ⎟ ⇒ + = ( A B C B ) u f − T 1 − B C p 0 ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ • A + B T C -1 B has same sparsity as A (discontinuous pressure p ) • BDD preconditioner with enriched coarse space investigated by Goldfeld (2003) • related work for incompressible problems in Pavarino and Widlund (2002) and Li (2001), see also M7
recall + = ( A B C B ) u f − T 1 Let A = µ A 0 and C = (1/ λ ) M p where ν E E µ = λ = + ν + ν − ν 2 ( 1 ) ( 1 )( 1 2 ) ⇒ + = µ + λ A B C B A B M B − − T 1 T 1 0 p notice as ν → ½ that λ → ∞ ⇒ condition number → ∞
2D Plane Strain Example condition number estimates for 4 substructure problem with H/h = 8 ν κ of preconditioned system 0.3 2.0 0.4 1.8 0.49 2.1 κ ≈ 1/(1 − 2 ν ) for ν near ½ 0.499 2.2 0.4999 7.2 0.49999 63 Q: why does κ → ∞ as ν → ½ for 0.499999 6.2e2 BDDC preconditioner? 0.4999999 6.2e3
• Explanation – Movement of substucture boundary nodes is weighted average of coarse grid and substructure corrections (v 1 and v 2 ) from neighboring subs – If substructure volume changes, then strain energy of substructure → ∞ as ν → ½ – cg step length α → 0 since cg minimizes energy • Solution – Modify “standard” BDDC constraint equations to enable enforcement of zero volume change – Same number of constraints in 2D slightly more in 3D
Recommend
More recommend