Deflated Hermitian Lanczos Methods for Multiple Right-Hand Sides

For Hermitian positive definite linear systems and eigenvalue problems, the eigCG algorithm is a memory efficient algorithm that solves the linear system and simultaneously computes some of its eigenvalues. The algorithm is based on the Conjugate-Gradient (CG) algorithm, however, it uses only a window of the vectors generated by the CG algorithm to compute approximate eigenvalues. The number and accuracy of the eigenvectors can be increased by solving more right-hand sides. For Hermitian systems with multiple right-hand sides, the computed eigenvectors can be used to speed up the solution of subsequent systems. The algorithm was tested on Lattice QCD problems by solving the normal equations and was shown to give large speed up factors and to remove the critical slowing down as we approach light quark masses. Here, an extension to the non-symmetric case based on the two-sided Lanczos algorithm is given. The new algorithm is tested on Lattice QCD problems and is shown to give promising results. We also study the removal of the critical slowing down and compare results with those of the eigCG algorithm. We also discuss the case when the system is gamma5-Hermitian.

Most calculations in lattice Quantum Chromodynamics (QCD) involve the solution of a series of linear systems of equations with exceedingly large matrices and a large number of right hand sides. Iterative methods for these problems can be sped up significantly if we deflate approximations of appropriate invariant spaces from the initial guesses. Recently we have developed eigCG, a modification of the Conjugate Gradient (CG) method, which while solving a linear system can reuse a window of the CG vectors to compute eigenvectors almost as accurately as the Lanczos method. The number of approximate eigenvectors can increase as more systems are solved. In this paper we review some of the characteristics of eigCG and show how it helps remove the critical slowdown in QCD calculations. Moreover, we study scaling with lattice volume and an extension of the technique to nonsymmetric problems.

This paper focuses on total least squares (TLS) problems $AX\approx B$ with multiple right-hand sides. Existence and uniqueness of a TLS solution for such problems was analyzed in the paper [I. Hnětynková et al., SIAM J. Matrix Anal. Appl., 32, 2011, pp. 748--770]. For TLS problems with single right-hand sides the paper [C. C. Paige and Z. Strakoš, SIAM J. Matrix Anal. Appl., 27, 2006, pp. 861--875] showed how necessary and sufficient information for solving $Ax\approx b$ can be revealed from the original data through the so-called core problem concept. In this paper we present a theoretical study extending this concept to problems with multiple right-hand sides. The data reduction we present here is based on the singular value decomposition of the system matrix $A$. We show minimality of the reduced problem; in this sense the situation is analogous to the single right-hand side case. Some other properties of the core problem, however, cannot be extended to the case of multiple right-hand sides.

We present a seed method for solving large nonsymmetric linear systems with multiple right-hand sides. The method uses a single augmented Krylov subspace corresponding to a seed system as a generator of approximations to the nonseed systems. The residual evaluate of the method is shown, and a new strategy to form a seed system which could supply information shareable among the right-hand sides is given. Numerical experiments indicate that our seed selection strategy is more efficient than two existing strategies and our method has significant time saving compared with the block GMRES method and the GMRES method with a projection process.

Block Krylov subspace methods are efficient iterative methods for solving linearsystems with multiple right-hand sides. The Block BiCGGR method is one of the Block Krylov subspace methods, and can generate highly accurate approximatesolutions. However, when the number of right-hand sides is large, this method becomes numerically unstable. In this paper, the Block BiCGGR2 method is proposed. This method is a numerically stable and highly accurate method than the Block BiCGGR method and the Block BiCGSTAB method. Numerical results illustrate the robustness and the accuracy of the approximate solutions of the Block BiCGGR2 method.

Contemporary large-scale Partial Differential Equation (PDE) simulations usually require the solution of large and sparse linear systems. Moreover, it is often needed to solve these linear systems with different or multiple Right-Hand Sides (RHSs). In this paper, various strategies will be presented to extend the scalability of existing multigrid or domain decomposition linear solvers using appropriate recycling strategies or block methods—i.e., by treating multiple right-hand sides simultaneously. The scalability of this work is assessed by performing simulations on up to 8,192 cores for solving linear systems arising from various physical phenomena modeled by Poisson's equation, the system of linear elasticity, or Maxwell's equation. This work is shipped as part of on open-source software, readily available and usable in any C/C++, Python, or Fortran code. In particular, some simulations are performed on top of a well-established library, PETSc, and it is shown how our approaches can be used to decrease time to solution down by 30%.

A new resolution method for solving linear systems resulting from the discrete adjoint approach in aerodynamic shape optimization with a single objective and multiple constraint functions is presented. The steady flow is governed by the three-dimensional compressible Favre-averaged Navier–Stokes equations combined with a two-equation turbulence model. Block flexible variants of the generalized minimal residual method are Krylov methods designed for the solution of linear systems with multiple right-hand sides allowing variable preconditioning. Then, they can be applied to the discrete adjoint equations of Favre-averaged Navier–Stokes equations for all objective and constraint functions as a right-hand side. The paper focuses on the development of a new algorithm combining the three following features: flexibility, multiple right-hand-side concurrent processing, and spectral deflation. This latter property enables to recycle approximate spectral information from one restart to the next one in order to avoid generalized minimal residual stagnation when short restart recurrences are performed. This novel block method is named block flexible generalized minimal residual with deflated restarting. It is evaluated on two aerodynamic configurations showing its capabilities to reduce computational times over standard Krylov methods.

In this paper, we introduce two new methods for solving large sparse nonsymmetric linear systems with several right-hand sides. These methods are the global Hessenberg and global CMRH methods. Using the global Hessenberg process, these methods are less expensive than the global FOM and global GMRES methods [9]. Theoretical results about the new methods are given, and experimental results that show good performances of these new methods are presented.

The block conjugate orthogonal conjugate gradient method (BCOCG) is recognized as a common method to solve complex symmetric linear systems with multiple right-hand sides. However, breakdown always occurs if the right-hand sides are rank deficient. In this paper, based on the orthogonality conditions, we present a breakdown-free BCOCG algorithm with new parameter matrices to handle rank deficiency. To improve the spectral properties of coefficient matrix A, a precondition version of the breakdown-free BCOCG is proposed in detail. We also give the relative algorithms for the block conjugate A-orthogonal conjugate residual method. Numerical results illustrate that when breakdown occurs, the breakdown-free algorithms yield faster convergence than the non-breakdown-free algorithms.

The controlled Cholesky factorization (CCF) has been shown to be a robust preconditioner for the conjugate gradient (CG) method. In this scheme the amount of fill-in is defined in terms of a parameter $\eta$, the number of extra nonzero elements allowed per column. This parameter $\eta$ can be chosen in such a way that the preconditioning factor will require more or less storage than the coefficient matrix. The selection of nonzero elements kept for the preconditioner is made by value instead of by position, so the sparsity pattern of the original matrix is never considered by the proposed factorization. This preconditioned method is termed the CCCG($\eta$) method and is a generalization of a method developed by Jones and Plassmann [\textit{ACM Trans.\ Math.\ Software}, 21 (1995), pp.\ 5--17]. It is demonstrated how an optimum value of $\eta$ can be automatically determined when solving systems with multiple right-hand sides. Typically such systems occur when using time-stepping procedures. Numerical examples are presented for some time-dependent parabolic partial differential equations discretized using linear triangular finite elements on unstructured grids. Comparisons made between CCCG($\eta$), standard ICCG, and the Jones and Plassmann method show CCCG($\eta$) to be an efficient solver.

Objective. Real-time reconstruction of magnetic particle imaging (MPI) shows promising clinical applications. However, prevalent reconstruction methods are mainly based on serial iteration, which causes large delay in real-time reconstruction. In order to achieve lower latency in real-time MPI reconstruction, we propose a parallel method for accelerating the speed of reconstruction methods. Approach. The proposed method, named adaptive multi-frame parallel iterative method (AMPIM), enables the processing of multi-frame signals to multi-frame MPI images in parallel. To facilitate parallel computing, we further propose an acceleration strategy for parallel computation to improve the computational efficiency of our AMPIM. Main results. OpenMPIData was used to evaluate our AMPIM, and the results show that our AMPIM improves the reconstruction frame rate per second of real-time MPI reconstruction by two orders of magnitude compared to prevalent iterative algorithms including the Kaczmarz algorithm, the conjugate gradient normal residual algorithm, and the alternating direction method of multipliers algorithm. The reconstructed image using AMPIM has high contrast-to-noise with reducing artifacts. Significance. The AMPIM can parallelly optimize least squares problems with multiple right-hand sides by exploiting the dimension of the right-hand side. AMPIM has great potential for application in real-time MPI imaging with high imaging frame rate.

Many scientific applications require one to solve successively linear systems Ax=b with different right-hand sides b and a symmetric positive definite matrix A. The conjugate gradient method applied to the first system generates a Krylov subspace which can be efficiently recycled thanks to orthogonal projections in subsequent systems. A modified conjugate gradient method is then applied with a specific initial guess and initial descent direction and a modified descent direction during the iterations. This paper gives new theoretical results for this method and proposes a new version. Numerical experiments show the efficacy of our method even for quite different right-hand sides.

This paper focuses on saddle point problems with a 2-by-2 block coefficient matrix. When the number of columns in the upper-right block and the number of rows in the lower-left block of the coefficient matrix is large, the convergence behavior of Krylov subspace methods for the saddle point problems tends to be poor even if the upper-left block is a well-conditioned matrix. In this paper, an efficient approach for solving the saddle point problems using block structure of the problems is proposed. The most time-consuming part of our proposed approach is the solution of a linear system with multiple right-hand sides. To solve the linear system with multiple right-hand sides efficiently, we propose to apply Block Krylov subspace methods to this linear system. Numerical experiments show that the proposed approach with Block Krylov subspace methods can solve the saddle point problems more efficiently than the conventional approach in terms of the number of iterations and the computation time.

In this paper, we consider using conjugate gradient (CG) methods for solving multiple linear systems $A^{(i)} x^{(i)} = b^{(i)},$ for $1 \le i \le s,$ where the coefficient matrices $A^{(i)}$ and the right-hand sides $b^{(i)}$ are different in general. In particular, we focus on the seed projection method which generates a Krylov subspace from a set of direction vectors obtained by solving one of the systems, called the seed system, by the CG method and then projects the residuals of other systems onto the generated Krylov subspace to get the approximate solutions. The whole process is repeated until all the systems are solved. Most papers in the literature [T. F. Chan and W. L. Wan, {\it SIAM J. Sci. Comput.}, 18 (1997), 1698--1721; B. Parlett {\it Linear Algebra Appl.}, 29 (1980), 323--346; Y. Saad, {\it Math. Comp.}, 48 (1987), 651--662; V. Simoncini and E. Gallopoulos, {\it SIAM J. Sci. Comput.}, 16 (1995), 917--933; C. Smith, A. Peterson, and R. Mittra, {\it IEEE Trans. Antennas and Propagation}, 37 (1989), pp. 1490--1493] considered only the case where the coefficient matrices $A^{(i)}$ are the same but the right-hand sides are different. We extend and analyze the method to solve multiple linear systems with varying coefficient matrices and right-hand sides. A theoretical error bound is given for the approximation obtained from a projection process onto a Krylov subspace generated from solving a previous linear system. Finally, numerical results for multiple linear systems arising from image restorations and recursive least squares computations are reported to illustrate the effectiveness of the method.

We describe a new technique for solvinga sparse linear system Ax = b as a block system AX = B, where multiple startingv ectors and right-hand sides are chosen so as to accelerate convergence. Efficiency is gained by reusing the matrix A in block operations with X and B. Techniques for reducingthe cost of the extra matrix-vector operations are presented.KeywordsKrylov SubspaceApproximation SpaceLinear SolverMemory HierarchyKrylov Subspace MethodThese keywords were added by machine and not by the authors. This process is experimental and the keywords may be updated as the learning algorithm improves.