A framework for deflated BiCG and related solvers Martin H. - - PowerPoint PPT Presentation

A framework for deflated BiCG and related solvers

Martin H. Gutknecht

Seminar for Applied Mathematics, ETH Zurich

SIAM Conference on Applied Linear Algebra Valencia, 20 June 2012 Joint work with André Gaul, Jörg Liesen, Reinhard Nabben

Augment./Deflat. History CG BiCG Conclusions

Outline

Augmentation and Deflation: Basics History Galerkin: CG, CR, GCR, ... Petrov-Galerkin: BICG, BICR Conclusions

M.H. Gutknecht SIAM-ALA12

• p. 2

Augment./Deflat. History CG BiCG Conclusions

Iterative methods based on (Petrov-)Galerkin conditions

To solve: Ax = b with A ∈ CN×N nonsingular. Construct sequence xn such that rn :≡ b − Axn → o . Choose xn from search space x0 + Sn such that some Galerkin or Petrov-Galerkin condition is satisfied: xn ∈ x0 + Sn , rn = A(x⋆ − xn) ⊥ BH Sn with some (formal) inner product matrix BH , i.e., rn ∈ r0 + ASn , rn ⊥ BH Sn . r0 is approximated from ASn such that “error” rn ⊥ Sn . Two cases: Galerkin: Sn = Sn , Petrov-Galerkin: Sn = Sn

B :≡ I for CG/BICG, B :≡ AH for CR/GCR, B :≡ A for BICR.

M.H. Gutknecht SIAM-ALA12

• p. 3

Augment./Deflat. History CG BiCG Conclusions

Augmentation and deflation

Ass.: know basis U of approximately A-invariant subspace U, i.e., U = R(U), where U ∈ CN×k full rank. Search space Sn and test space Sn are split up: xn ∈ x0 + Sn , rn ⊥ BH Sn , Sn :≡ Kn ⊕ U ,

• Sn :≡

Ln ⊕ U ,

• Kn :≡ Kn(

A, r0) :≡ span { r0, . . . , An−1 r0} . Still to be specified:

• A,
• r0,
• Ln,
• U.

Galerkin case: Ln :≡ Kn,

• U :≡ U.

Petrov-Galerkin case: e.g., Ln :≡ Kn( AH,

• r0),

but other options for Kn and Ln exist.

M.H. Gutknecht SIAM-ALA12

• p. 4

Augment./Deflat. History CG BiCG Conclusions

Rationale of augmentation and deflation

Ideally: columns of U ∈ CN×k span A-invariant subspace U belonging to eigenvalues close to 0 . Let Z :≡ AU , Z :≡ AU = U. Note: images of the restriction A−1

• Z are trivial to compute:

if z = Zc ∈ Z , then A−1z = Uc . Choose projector P such that N(P) = Z . Split up space: CN = R(P) ⊕ Z. Choose A :≡ PA and r0 :≡ Pr0 so that Kn ⊆ R(P). Split up r0 accordingly: r0 =

• r0
• ∈ R(P)

+ r0 − r0 ∈ Z . A−1 (r0 − r0) is trivial;

• A−1

r0 is found with Krylov space solver acting on R(P) .

M.H. Gutknecht SIAM-ALA12

• p. 5

Augment./Deflat. History CG BiCG Conclusions

Since N(P) = Z , we have N( A) = N(PA) = U . So the (absolutely) small eigenvalues of A represented by U that caused trouble are replaced by a k-fold EVal o in A (deflation). Projector P is not fully determined by its nullspace since it may be oblique. Hopefully, A

• R(P) =

A

• R(P) .

(1) If A is Hermitian and U is A-invariant, and if P is chosen such that R(P) = U⊥ or R(P) = Z⊥ , Eq. (1) holds. If R(P) = Z⊥ , P is an orthogonal projector.

M.H. Gutknecht SIAM-ALA12

• p. 6

Augment./Deflat. History CG BiCG Conclusions

How to find an approximately invariant subspace?

It may be known from a theoretical analysis of the problem. It may result from the solution of previous systems with the same A . ( linear system with multiple right-hand sides) It may results from the solution of previous systems with nearby A . It may results from previous cycles of the solution process (if the method is restarted). There are lots of examples in the literature.

M.H. Gutknecht SIAM-ALA12

• p. 7

Augment./Deflat. History CG BiCG Conclusions

Things to distinguish

Augmented bases: xn ∈ x0 + Kn( A, r0) + U , where

• A = A
• r

spec( A) ⊂ spec(A) ∪ {0} (Spectral) deflation: A A :≡ PA s.t. small EVals 0 EVal translation: A A :≡ AP s.t. small EVals |λ

Krylov space recycling: choice of U based on prev. cycles Flexible KSS: adaptation of P at each restart Two basic approaches: Augmentation of basis with or without spectral deflation. EVal translation by suitable preconditioning (no deflation!).

M.H. Gutknecht SIAM-ALA12

• p. 8

Augment./Deflat. History CG BiCG Conclusions

History

Early contributions (many more papers appeared since): Nicolaides ’85/’87SINUM: deflated 3-term CG (w/augm. basis) Dostál ’87/’88IntJCompMath: deflated 2-term CG (w/augm. basis) Morgan ’93/’95SIMAX: GMRES with augmented basis de Sturler ’93/’96JCAM: inner-outer GMRES/GCR (and, briefly, inner-outer BiCGStab/GCR) with augmented basis Chapman / Saad ’95/’97NLAA GMRES with augmented basis Saad ’95/’97SIMAX Analysis of KSS with augmented basis de Sturler ’97/’99SINUM inner-outer GMRES/GCR w/truncation Vuik / Segal / Meijerink ’98/’99JCP 2-term CG w/augm. basis Bristeau / Erhel ’98/’98NumAlg CG with augmented basis Erhel / Guyomarc’h ’97/’00SIMAX defl. 2-term CG w/augm. basis Saad / Yeung / Erhel / Guyomarc’h ’98/’00SISC the same

M.H. Gutknecht SIAM-ALA12

• p. 9

Augment./Deflat. History CG BiCG Conclusions

Galerkin: CG, CR, GCR, ... : some details

Given: A, B ∈ CN×k and U ∈ CN×k Most relevant cases: B = I for CG, B = AH for CR, GCR E :≡ UHBAU ∈ Ck×k , assumed nonsingular, M :≡ UE−1UH , P :≡ I − AMB , projector onto (BH U)⊥ along Z , Q :≡ I − MBA , projector onto (AHBH U)⊥ along U ,

• A :≡ PA = AQ = PAQ

with N( A) = U , R( A) = (BHU)⊥,

• rn :≡ P(b − A

xn) = Pb − A xn ∈ (BHU)⊥ ,

• Kn :≡ K(

A, r0) ⊆ (BHU)⊥.

M.H. Gutknecht SIAM-ALA12

• p. 10

Augment./Deflat. History CG BiCG Conclusions

An equivalence theorem

THEOREM

For n ≥ 1 the two pairs of conditions, xn ∈ x0 + Kn + U , rn ⊥ B Kn + B U , (2) and

• xn ∈ x0 +

Kn ,

• rn ⊥ B

Kn . (3) are equivalent in the sense that xn = Q xn + MBHb and rn = rn . (4)

• DEF. The direct deflation approach is given by (2),

the indirect deflation approach is given by (3)–(4).

M.H. Gutknecht SIAM-ALA12

• p. 11

Augment./Deflat. History CG BiCG Conclusions

Computing xn satisfying (3) means solving the singular linear system

• A

x = Pb with a Krylov space solver characterized by (3). What are the properties of A ? N( A) = U , R( A) = (BHU)⊥. If U is A-invariant, the corresp. EVals become 0. What can we say about the others? Consider partitioned Jordan decomposition of A A = SJS−1 =

• S1

S2 J1 J2 SH

1

• SH

2

• ,

where S1, S1 ∈ CN×k, S2, S2 ∈ CN×(N−k) and either R(S1) = U or R( S1) = BHU.

M.H. Gutknecht SIAM-ALA12

• p. 12

Augment./Deflat. History CG BiCG Conclusions

Theorem

(1) If U = R(S1), if U ∈ CN×k is any matrix satisfying R(U) = U, and if UHBAU is nonsingular, then

• A = PA =
• U

PS2 J2 U PS2 −1 (5) with

• U

PS2 −1 =

• BHU(UHBHU)−1
• S2

H . (2) If BHU = R( S1), if U ∈ CN×k is any matrix satisfying R(U) = U, and if UHBAU is nonsingular, then

• A = PA =
• U

S2 J2 U S2 −1 (6) with

• U

S2 −1 =

• BHU(UHBHU)−1

QH S2 H . In particular, in both cases the spectrum Λ( A) of A is given

M.H. Gutknecht SIAM-ALA12

• p. 13

Augment./Deflat. History CG BiCG Conclusions

Petrov-Galerkin: BICG, BICR, ...

Generalized BICG (GENBICG) [G. ’90(CopperMtn), ’97ActaNum] with formal inner product matrix B requires A and B to commute (to maintain short recurrences). For deflated solvers based on Petrov-Galerkin condition we need projectors and operators for creating split dual spaces: Sn :≡ Kn ⊕ U ,

• Sn :≡

Ln ⊕ U . May consider solving two dual systems at once [Ahuja ’09Diss]: Ax = b , AH x = b such that xn ∈ x0 + Sn ,

• xn ∈

x0 + Sn , rn ⊥ BH Sn ,

• rn ⊥ BSn .

M.H. Gutknecht SIAM-ALA12

• p. 14

Augment./Deflat. History CG BiCG Conclusions

Petrov-Galerkin: Projectors and other operators

definition null space range further properties E

• UHBAU

{0} Ck assumed to be nonsingular M UE−1 UH

• U⊥

U rank M = k P I − AMB Z (BH U)⊥ P2 = P Q I − MBA U (AHBH U)⊥ Q2 = Q

• A

PA U (BH U)⊥

• A = PA = AQ = PAQ
• E

UHBHAH U {0} Ck assumed to be nonsingular

• M
• U

E−1UH U⊥

• U

rank M = k

• P

I − AH MBH

• Z

(BU)⊥

• P2 =

P

• Q

I − MBHAH

• U

(ABU)⊥

• Q2 =

Q

• A
• PAH
• U

(BU)⊥

• A =

PAH = AH Q = PAH Q

Krylov spaces used:

• Kn :≡ Kn(

A, r0) ,

• Ln :≡ Kn(

A, r0) .

M.H. Gutknecht SIAM-ALA12

• p. 15

Augment./Deflat. History CG BiCG Conclusions

THEOREM

For n ≥ 1 the two sets of conditions xn ∈ x0 + Kn + U , rn ⊥ BH( Ln + U) ,

• xn ∈

x0 + Ln + U ,

• rn ⊥ B(

Kn + U) and

• xn ∈ x0 +

Kn ,

• rn⊥ BH

Ln ,

• xn ∈

x0 + Ln ,

• rn ⊥ B

Kn are equivalent in the sense that xn = Q xn + MBb and rn = rn , (7)

• xn =

Q

• xn + MHBH

b and

• rn =
• rn .

(8)

M.H. Gutknecht SIAM-ALA12

• p. 16

Augment./Deflat. History CG BiCG Conclusions

Deflated BICG

For solving a single system Ax = b with deflated BICG the indirect approach requires to solve A x = Pb such that

• xn ∈ x0 +

Kn ,

• rn ⊥

Ln . This works with BICG if A = AH . Fortunately, we have:

THEOREM

If B = I, then PH = Q and QH = P , and therefore AH = A . It is not clear if one can obtain an equally efficient deflated BICR (where B = A). But BICGSTAB, IDR(s), etc. can be done.

M.H. Gutknecht SIAM-ALA12

• p. 17

Augment./Deflat. History CG BiCG Conclusions

Conclusions

We have set up a framework for deflated Krylov space solvers based on a Galerkin condition. We are developing an analogue framework for those based

• n a Petrov-Galerkin condition.

By our indirect approach we obtain a deflated BICG with coupled two-term recurrences (or three-term recurrences) which reduces to deflated CG when applied to an Hpd problem:

• If B = I, BICG applied to the projected problem
• A

x = Pb yields xn, rn, and shadow residuals

• rn.
• The transformation (7) yields xn.

M.H. Gutknecht SIAM-ALA12

• p. 18

Augment./Deflat. History CG BiCG Conclusions

References

• 1. A. Gaul, M.H.G., J. Liesen, R. Nabben, A framework for

deflated and augmented Krylov subspace methods. Submitted (Jan. 2011 / June 2012).

• 2. A. Gaul, M.H.G., J. Liesen, R. Nabben, Deflated and

Augmented Krylov Subspace Methods: A Framework for Deflated BiCG and Related Methods. In preparation.

• 3. M.H.G., Spectral Deflation in Krylov Solvers: a Theory of

Coordinate Space Based Methods. Electronic Trans.

• Numer. Anal. 39, 156–185 (2012).

M.H. Gutknecht SIAM-ALA12

• p. 19