petsc-3.13.6 2020-09-29
Report Typos and Errors

KSPHPDDM

Interface with the HPDDM library. This KSP may be used to further select methods that are currently not implemented natively in PETSc, e.g., GCRODR [2006], a recycled Krylov method which is similar to KSPLGMRES, see [2016] for a comparison. ex75.c shows how to reproduce the results from the aforementioned paper [2006]. A chronological bibliography of relevant publications linked with KSP available in HPDDM through KSPHPDDM, and not available directly in PETSc, may be found below.

Options Database Keys

-ksp_richardson_scale <scale, default=1.0> - see KSPRICHARDSON
-ksp_gmres_restart <restart, default=40> - see KSPGMRES
-ksp_hpddm_type <type, default=gmres> - any of gmres, bgmres, cg, bcg, gcrodr, bgcrodr, or bfbcg
-ksp_hpddm_deflation_tol <eps, default=-1.0> - tolerance when deflating right-hand sides inside block methods (no deflation by default, only relevant with block methods)
-ksp_hpddm_enlarge_krylov_subspace <p, default=1> - split the initial right-hand side into multiple vectors (only relevant with nonblock methods)
-ksp_hpddm_orthogonalization <type, default=cgs> - any of cgs or mgs, see KSPGMRES
-ksp_hpddm_qr <type, default=cholqr> - distributed QR factorizations with any of cholqr, cgs, or mgs (only relevant with block methods)
-ksp_hpddm_variant <type, default=left> - any of left, right, or flexible (this option is superseded by KSPSetPCSide())
-ksp_hpddm_recycle <n, default=0> - number of harmonic Ritz vectors to compute (only relevant with GCRODR or BGCRODR)
-ksp_hpddm_recycle_target <type, default=SM> - criterion to select harmonic Ritz vectors using either SM, LM, SR, LR, SI, or LI (only relevant with GCRODR or BGCRODR)
-ksp_hpddm_recycle_strategy <type, default=A> - generalized eigenvalue problem A or B to solve for recycling (only relevant with flexible GCRODR or BGCRODR)

References

1980 - The Block Conjugate Gradient Algorithm and Related Methods. O'Leary. Linear Algebra and its Applications.
2006 - Recycling Krylov Subspaces for Sequences of Linear Systems. Parks, de Sturler, Mackey, Johnson, and Maiti. SIAM Journal on Scientific Computing
2013 - A Modified Block Flexible GMRES Method with Deflation at Each Iteration for the Solution of Non-Hermitian Linear Systems with Multiple Right-Hand Sides. Calandra, Gratton, Lago, Vasseur, and Carvalho. SIAM Journal on Scientific Computing.
2016 - Block Iterative Methods and Recycling for Improved Scalability of Linear Solvers. Jolivet and Tournier. SC16.
2017 - A breakdown-free block conjugate gradient method. Ji and Li. BIT Numerical Mathematics.

See Also

KSPCreate(), KSPSetType(), KSPType (for list of available types), KSP, KSPGMRES, KSPCG, KSPLGMRES, KSPDGMRES

Level

intermediate

Location

src/ksp/ksp/impls/hpddm/hpddm.cxx

Examples

src/ksp/ksp/tutorials/ex75.c.html
src/ksp/ksp/tutorials/ex77.c.html
src/ksp/ksp/tutorials/ex75f.F90.html
src/ksp/ksp/tutorials/ex77f.F90.html

Index of all KSP routines
Table of Contents for all manual pages
Index of all manual pages