KSPTSIRM#

Implements the two-stage iteration with least-squares residual minimization method [CKG16]

Options Database Keys#

  • -ksp_ksp_type - the type of the inner solver (GMRES or any of its variants for instance)

  • -ksp_pc_type - the type of the preconditioner applied to the inner solver

  • -ksp_ksp_max_it - the maximum number of inner iterations (iterations of the inner solver)

  • -ksp_ksp_rtol - sets the relative convergence tolerance of the inner solver

  • -ksp_tsirm_cgls - if 1 use CGLS solver in the minimization step, otherwise use LSQR solver

  • -ksp_tsirm_max_it_ls - the maximum number of iterations for the least-squares minimization solver

  • -ksp_tsirm_tol_ls - sets the convergence tolerance of the least-squares minimization solver

  • -ksp_tsirm_size_ls - the number of residuals for the least-squares minimization step

Notes#

KSPTSIRM is a two-stage iteration method for solving large sparse linear systems of the form \(Ax=b\). The main idea behind this new method is the use a least-squares residual minimization to improve the convergence of Krylov based iterative methods, typically those of GMRES variants. The principle of TSIRM algorithm is to build an outer iteration over a Krylov method, called the inner solver, and to frequently store the current residual computed by the given Krylov method in a matrix of residuals S. After a few outer iterations, a least-squares minimization step is applied on the matrix composed by the saved residuals, in order to compute a better solution and to make new iterations if required. The minimization step consists in solving the least-squares problem \(\min||b-ASa||\) to find ‘a’ which minimizes the residuals \((b-AS)\). The minimization step is performed using two solvers of linear least-squares problems: KSPCGLS or KSPLSQR. A new solution x with a minimal residual is computed with \(x=Sa\).

Defaults to 30 iterations for the inner solve, use option -ksp_ksp_max_it <it> to change it.

Contributed by#

Lilia Ziane Khodja

References#

[CKG16]

Raphaël Couturier, Lilia Ziane Khodja, and Christophe Guyeux. TSIRM: a two-stage iteration with least-squares residual minimization algorithm to solve large sparse linear and nonlinear systems. Journal of Computational Science, 17:535–546, 2016.

See Also#

KSP: Linear System Solvers, KSPCreate(), KSPSetType(), KSPType, KSP, KSPFGMRES, KSPLGMRES, KSPGMRESSetRestart(), KSPGMRESSetHapTol(), KSPGMRESSetPreAllocateVectors(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESModifiedGramSchmidtOrthogonalization(), KSPGMRESCGSRefinementType, KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSPGMRESMonitorKrylov(), KSPSetPCSide()

Level#

advanced

Location#

src/ksp/ksp/impls/tsirm/tsirm.c


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