Engineering problems frequently require solving a sequence of dual linear systems. This paper introduces recycling BiCG, that recycles the Krylov subspace from one pair of linear systems to the next pair. Augmented bi-Lanczos algorithm and modified two-term recurrence are developed for using the recycle space. Recycle space is built from the approximate invariant subspace corresponding to eigenvalues close to the origin. Recycling approach is extended to the CGS and the BiCGSTAB algorithms. Experiments on a convection-diffusion problem give promising results. / Master of Science

Krylov subspace iterative methods provide an effective tool for reducing the solution of large linear systems to a size for which a direct solver may be applied. However, the problems of limited storage and speed are still a concern. Therefore, in this dissertation work, we present iterative Krylov subspace algorithms for non-Hermitian systems which do have fixed memory requirements and have favorable convergence characteristics. This dissertation describes three projects. The first project concerns short-term recurrence Krylov subspace methods for nearly-Hermitian linear systems. In 2008, Beckermann and Reichel introduced a short-term recurrence progressive GMRES algorithm for nearly-Hermitian linear systems. However, we have found this method to be unstable. We document the instabilities and introduce a different fixed-memory algorithm to treat nearly-Hermitian problems. We present numerical experiments demonstrating that the performance of this algorithm is competitive. The other two projects involve extending a strategy called Krylov subspace recycling, introduced by Parks and colleagues in 2005. This method requires more overhead than other subspace augmentation methods but offers the ability to recycle subspace information between cycles for a single linear system and recycle information between related linear systems. In the first project, we extend subspace recycling to the block Krylov subspace setting. A block Krylov subspace is a generalization of Krylov subspace where a single starting vector is replaced with a block of linearly independent starting vectors. We then apply our method to a sequence of matrices arising in a Newton iteration applied to fluid density functional theory and present some numerical experiments. In the second project, we extend the methods of subspace recycling to a family of linear systems differing only by multiples of the identity. These problems arise in the theory of quantum chromodynamics, a theory of the behavior of subatomic particles. We wish to build on the class of Krylov methods which allow the simultaneous solution of all shifted linear systems while generating only one subspace. However, the mechanics of subspace recycling complicates this situation and interferes with our ability to simultaneously solve all systems using these techniques. Therefore, we introduce an algorithm which avoids this complication and present some numerical experiments demonstrating its effectiveness. / Mathematics

Sequences of matrices arise in many applications in science and engineering. In this thesis we consider matrices that are closely related (or closely related in groups), and we take advantage of the small differences between them to efficiently solve sequences of linear systems and eigenproblems. Recycling techniques, such as recycling preconditioners or subspaces, are popular approaches for reducing computational cost. In this thesis, we introduce two novel approaches for recycling previously computed information for a subsequent system or eigenproblem, and demonstrate good results for sequences arising in several applications. Preconditioners are often essential for fast convergence of iterative methods. However, computing a good preconditioner can be very expensive, and when solving a sequence of linear systems, we want to avoid computing a new preconditioner too often. Instead, we can recycle a previously computed preconditioner, for which we have good convergence behavior of the preconditioned system. We propose an update technique we call the sparse approximate map, or SAM update, that approximately maps one matrix to another matrix in our sequence. SAM updates are very cheap to compute and apply, preserve good convergence properties of a previously computed preconditioner, and help to amortize the cost of that preconditioner over many linear solves. When solving a sequence of eigenproblems, we can reduce the computational cost of constructing the Krylov space starting with a single vector by warm-starting the eigensolver with a subspace instead. We propose an algorithm to warm-start the Krylov-Schur method using a previously computed approximate invariant subspace. We first compute the approximate Krylov decomposition for a matrix with minimal residual, and use this space to warm-start the eigensolver. We account for the residual matrix when expanding, truncating, and deflating the decomposition and show that the norm of the residual monotonically decreases. This method is effective in reducing the total number of matrix-vector products, and computes an approximate invariant subspace that is as accurate as the one computed with standard Krylov-Schur. In applications where the matrix-vector products require an implicit linear solve, we incorporate Krylov subspace recycling. Finally, in many applications, sequences of matrices take the special form of the sum of the identity matrix, a very low-rank matrix, and a small-in-norm matrix. We consider convergence rates for GMRES applied to these matrices by identifying the sources of sensitivity. / Doctor of Philosophy / Problems in science and engineering often require the solution to many linear systems, or a sequence of systems, that model the behavior of physical phenomena. In order to construct highly accurate mathematical models to describe this behavior, the resulting matrices can be very large, and therefore the linear system can be very expensive to solve. To efficiently solve a sequence of large linear systems, we often use iterative methods, which can require preconditioning techniques to achieve fast convergence. The preconditioners themselves can be very expensive to compute. So, we propose a cheap update technique that approximately maps one matrix to another in the sequence for which we already have a good preconditioner. We then combine the preconditioner and the map and use the updated preconditioner for the current system. Sequences of eigenvalue problems also arise in many scientific applications, such as those modeling disk brake squeal in a motor vehicle. To accurately represent this physical system, large eigenvalue problems must be solved. The behavior of certain eigenvalues can reveal instability in the physical system but to identify these eigenvalues, we must solve a sequence of very large eigenproblems. The eigensolvers used to solve eigenproblems generally begin with a single vector, and instead, we propose starting the method with several vectors, or a subspace. This allows us to reduce the total number of iterations required by the eigensolver while still producing an accurate solution. We demonstrate good results for both of these approaches using sequences of linear systems and eigenvalue problems arising in several real-world applications. Finally, in many applications, sequences of matrices take the special form of the sum of the identity matrix, a very low-rank matrix, and a small-in-norm matrix. We examine the convergence behavior of the iterative method GMRES when solving such a sequence of matrices.

Science and engineering problems frequently require solving a sequence of single linear systems or a sequence of dual linear systems. We develop algorithms that recycle Krylov subspaces and preconditioners from one system (or pair of systems) in the sequence to the next, leading to efficient solutions. Besides the benefit of only having to store few Lanczos vectors, using BiConjugate Gradients (BiCG) to solve dual linear systems may have application-specific advantages. For example, using BiCG to solve the dual linear systems arising in interpolatory model reduction provides a backward error formulation in the model reduction framework. Using BiCG to evaluate bilinear forms -- for example, in the variational Monte Carlo (VMC) algorithm for electronic structure calculations -- leads to a quadratic error bound. Since one of our focus areas is sequences of dual linear systems, we introduce recycling BiCG, a BiCG method that recycles two Krylov subspaces from one pair of dual linear systems to the next pair. The derivation of recycling BiCG also builds the foundation for developing recycling variants of other bi-Lanczos based methods like CGS, BiCGSTAB, BiCGSTAB2, BiCGSTAB(l), QMR, and TFQMR. We develop a generalized bi-Lanczos algorithm, where the two matrices of the bi-Lanczos procedure are not each other's conjugate transpose but satisfy this relation over the generated Krylov subspaces. This is sufficient for a short term recurrence. Next, we derive an augmented bi-Lanczos algorithm with recycling and show that this algorithm is a special case of generalized bi-Lanczos. The Petrov-Galerkin approximation that includes recycling in the iteration leads to modified two-term recurrences for the solution and residual updates. We generalize and extend the framework of our recycling BiCG to CGS, BiCGSTAB and BiCGSTAB2. We perform extensive numerical experiments and analyze the generated recycle space. We test all of our recycling algorithms on a discretized partial differential equation (PDE) of convection-diffusion type. This PDE problem provides well-known test cases that are easy to analyze further. We use recycling BiCG in the Iterative Rational Krylov Algorithm (IRKA) for interpolatory model reduction and in the VMC algorithm. For a model reduction problem, we show up to 70% savings in iterations, and we also demonstrate that solving the problem without recycling leads to (about) a 50% increase in runtime. Experiments with recycling BiCG for VMC gives promising results. We also present an algorithm that recycles preconditioners, leading to a dramatic reduction in the cost of VMC for large(r) systems. The main cost of the VMC method is in constructing a sequence of Slater matrices and computing the ratios of determinants for successive Slater matrices. Recent work has improved the scaling of constructing Slater matrices for insulators, so that the cost of constructing Slater matrices in these systems is now linear in the number of particles. However, the cost of computing determinant ratios remains cubic in the number of particles. With the long term aim of simulating much larger systems, we improve the scaling of computing determinant ratios in the VMC method for simulating insulators by using preconditioned iterative solvers. The main contribution here is the development of a method to efficiently compute for the Slater matrices a sequence of preconditioners that make the iterative solver converge rapidly. This involves cheap preconditioner updates, an effective reordering strategy, and a cheap method to monitor instability of ILUTP preconditioners. Using the resulting preconditioned iterative solvers to compute determinant ratios of consecutive Slater matrices reduces the scaling of the VMC algorithm from O(n^3) per sweep to roughly O(n^2), where n is the number of particles, and a sweep is a sequence of n steps, each attempting to move a distinct particle. We demonstrate experimentally that we can achieve the improved scaling without increasing statistical errors. / Ph. D.

The main theme of this work is effectiveness and efficiency of Krylov subspace methods and Krylov subspace recycling. While solving long, slowly changing sequences of large linear systems, such as the ones that arise in engineering, there are many issues we need to consider if we want to make the process reliable (converging to a correct solution) and as fast as possible. This thesis is built on three main components. At first, we target bilinear and quadratic form estimation. Bilinear form $c^TA^{-1}b$ is often associated with long sequences of linear systems, especially in optimization problems. Thus, we devise algorithms that adapt cheap bilinear and quadratic form estimates for Krylov subspace recycling. In the second part, we develop a hybrid recycling method that is inspired by a complex CFD application. We aim to make the method robust and cheap at the same time. In the third part of the thesis, we optimize the implementation of Krylov subspace methods on Graphic Processing Units (GPUs). Since preconditioners based on incomplete matrix factorization (ILU, Cholesky) are very slow on the GPUs, we develop a preconditioner that is effective but well suited for GPU implementation. / Ph. D. / In many applications we encounter the repeated solution of a large number of slowly changing large linear systems. The cost of solving these systems typically dominates the computation. This is often the case in medical imaging, or more generally inverse problems, and optimization of designs. Because of the size of the matrices, Gaussian elimination is infeasible. Instead, we find a sufficiently accurate solution using iterative methods, so-called Krylov subspace methods, that improve the solution with every iteration computing a sequence of approximations spanning a Krylov subspace. However, these methods often take many iterations to construct a good solution, and these iterations can be expensive. Hence, we consider methods to reduce the number of iterations while keeping the iterations cheap. One such approach is Krylov subspace recycling, in which we recycle judiciously selected subspaces from previous linear solves to improve the rate of convergence and get a good initial guess. In this thesis, we focus on improving efficiency (runtimes) and effectiveness (number of iterations) of Krylov subspace methods. The thesis has three parts. In the first part, we focus on efficiently estimating sequences of bilinear forms, c^TA⁻¹b. We approximate the bilinear forms using the properties of Krylov subspaces and Krylov subspace solvers. We devise an algorithm that allows us to use Krylov subspace recycling methods to efficiently estimate bilinear forms, and we test our approach on three applications: topology optimization for the optimal design of structures, diffuse optical tomography, and error estimation and grid adaptation in computational fluid dynamics. In the second part, we focus on finding the best strategy for Krylov subspace recycling for two large computational fluid dynamics problems. We also present a new approach, which lets us reduce the computational cost of Krylov subspace recycling. In the third part, we investigate Krylov subspace methods on Graphics Processing Units. We use a lid driven cavity problem from computational fluid dynamics to perform a thorough analysis of how the choice of the Krylov subspace solver and preconditioner influences runtimes. We propose a new preconditioner, which is designed to work well on Graphics Processing Units.