On modern large-scale parallel computers, the performance of Krylov subspace iterative methods is limited by global synchronization. This has inspired the development of @math -step Krylov subspace method variants, in which iterations are computed in blocks of @math , which can reduce the number of global synchronizations per iteration by a factor of @math . Although the @math -step variants are mathematically equivalent to their classical counterparts, they can behave quite differently in finite precision depending on the parameter @math . If @math is chosen too large, the @math -step method can suffer a convergence delay and a decrease in attainable accuracy relative to the classical method. This makes it difficult for a potential user of such methods - the @math value that minimizes the time per iteration may not be the best @math for minimizing the overall time-to-solution, and further may cause an unacceptable decrease in accuracy. Towards improving the reliability and usability of @math -step Krylov subspace methods, in this work we derive the , a variable @math -step CG method where in block @math , the parameter @math is determined automatically such that a user-specified accuracy is attainable. The method for determining @math is based on a bound on growth of the residual gap within block @math , from which we derive a constraint on the condition numbers of the computed @math -dimensional Krylov subspace bases. The computations required for determining the block size @math can be performed without increasing the number of global synchronizations per block. Our numerical experiments demonstrate that the adaptive @math -step CG method is able to attain up to the same accuracy as classical CG while still significantly reducing the total number of global synchronizations.
In both @math -step and classical variants of Krylov subspace methods, finite precision roundoff error in updates to the approximate solution @math and the residual @math in each iteration can cause the updated residual @math and the true residual @math to grow further and further apart as the iterations proceed. If this deviation grows large, it can limit the , i.e., the accuracy with which we can solve @math on a computer with unit round-off @math . Analyses of maximum attainable accuracy in CG and other classical KSMs are given by Greenbaum @cite_14 , van der Vorst and Ye @cite_40 , Sleijpen, van der Vorst, and Fokkema @cite_12 , Sleijpen, van der Vorst, and Modersitzki @cite_16 , Bj "o rck, Elfving, and Strako s @cite_3 , and Gutknecht and Strako s @cite_33 . One important result of these analyses is the insight that loss of accuracy can be caused at a very early stage of the computation, which can not be corrected in later iterations. Analyses of the maximum attainable accuracy in @math -step CG and the @math -step biconjugate gradient method (BICG) can be found in @cite_19 @cite_6 .
@cite_14 @cite_33 @cite_3 @cite_6 @cite_19 @cite_40 @cite_16 @cite_12
2060710554 1560263223 2033180095 2002054778
It has been widely observed that Krylov space solvers based on two three-term recurrences can give significantly less accurate residuals than mathematically equivalent solvers implemented with three two-term recurrences. In this paper we attempt to clarify and justify this difference theoretically by analyzing the gaps between recursively and explicitly computed residuals. It is shown that, in contrast with the two-term recurrences analyzed by Sleijpen, van der Vorst, and Fokkema [ Numer. Algorithms, 7 (1994), pp. 75--109] and Greenbaum [SIAM J. Matrix Anal. Appl., 18 (1997), pp. 535--551], in the two three-term recurrences the contributions of the local roundoff errors to the analyzed gaps may be dramatically amplified while propagating through the algorithm. This result explains, for example, the well-known behavior of three-term-based versions of the biconjugate gradient method, where large gaps between recursively and explicitly computed residuals are not uncommon. For the conjugate gradient method, however, such a devastating behavior---although possible---is not observed frequently in practical computations, and the difference between two-term and three-term implementations is usually moderate or small. This can also be explained by our results. List of Algorithms Preface 1. Introduction. Brief Overview of the State of the Art Notation Review of Relevant Linear Algebra Part I. Krylov Subspace Approximations. 2. Some Iteration Methods. Simple Iteration Orthomin(1) and Steepest Descent Orthomin(2) and CG Orthodir, MINRES, and GMRES Derivation of MINRES and CG from the Lanczos Algorithm 3. Error Bounds for CG, MINRES, and GMRES. Hermitian Problems-CG and MINRES Non-Hermitian Problems-GMRES 4. Effects of Finite Precision Arithmetic. Some Numerical Examples The Lanczos Algorithm A Hypothetical MINRES CG Implementation A Matrix Completion Problem Orthogonal Polynomials 5. BiCG and Related Methods. The Two-Sided Lanczos Algorithm The Biconjugate Gradient Algorithm The Quasi-Minimal Residual Algorithm Relation Between BiCG and QMR The Conjugate Gradient Squared Algorithm The BiCGSTAB Algorithm Which Method Should I Use? 6. Is There A Short Recurrence for a Near-Optimal Approximation? The Faber and Manteuffel Result Implications 7. Miscellaneous Issues. Symmetrizing the Problem Error Estimation and Stopping Criteria Attainable Accuracy Multiple Right-Hand Sides and Block Methods Computer Implementation Part II. Preconditioners. 8. Overview and Preconditioned Algorithms. 9. Two Example Problems. The Diffusion Equation The Transport Equation 10. Comparison of Preconditioners. Jacobi, Gauss--Seidel, SOR The Perron--Frobenius Theorem Comparison of Regular Splittings Regular Splittings Used with the CG Algorithm Optimal Diagonal and Block Diagonal Preconditioners 11. Incomplete Decompositions. Incomplete Cholesky Decomposition Modified Incomplete Cholesky Decomposition 12. Multigrid and Domain Decomposition Methods. Multigrid Methods Basic Ideas of Domain Decomposition Methods. The conjugate gradient method applied to the normal equations ATAx=ATb (CGLS) is often used for solving large sparse linear least squares problems. The mathematically equivalent algorithm LSQR based on the Lanczos bidiagonalization process is an often recommended alternative. In this paper, the achievable accuracy of different conjgate gradient and Lanczos methods in finite precision is studied. It is shown that an implementation of algorithm CGLS in which the residual sk=AT(b-Axk) of the normal equations is recurred will not in general achieve accurate solutions. The same conclusion holds for the method based on Lanczos bidiagonalization with starting vector ATb. For the preferred implementation of CGLS we bound the error ||r-rk|| of the computed residual rk. Numerical tests are given that confirm a conjecture of backward stability. The achievable accuracy of LSQR is shown to be similar. The analysis essentially also covers the preconditioned case. Krylov subspace methods are a popular class of iterative methods for solving linear systems with large, sparse matrices. On modern computer architectures, both sequential and parallel performance of classical Krylov methods is limited by costly data movement, or communication, required to update the approximate solution in each iteration. These motivated communication-avoiding Krylov methods, based on @math -step formulations, reduce data movement by a factor of @math by reordering the computations in classical Krylov methods to exploit locality. Studies on the finite precision behavior of communication-avoiding Krylov methods in the literature have thus far been empirical in nature; in this work, we provide the first quantitative analysis of the maximum attainable accuracy of communication-avoiding Krylov subspace methods in finite precision. Following the analysis for classical Krylov methods, we derive a bound on the deviation of the true and updated residuals in communication-avoiding conjugate gradient...