For sequences of related linear systems, the computation of a preconditioner for every system can be expensive. Often a fixed preconditioner is used, but this may not be effective as the matrix changes. This research examines the benefits of both reusing and recycling preconditioners, with special focus on ILUTP and factorized sparse approximate inverses and proposes an update that we refer to as a sparse approximate map or SAM update. Analysis of the residual and eigenvalues of the map will be provided. Applications include the Quantum Monte Carlo method, model reduction, oscillatory hydraulic tomography, diffuse optical tomography, and Helmholtz-type problems. / Master of Science

Please download this dissertation to see the abstract.

The Method of Moments is a numerical technique for solving electromagnetic problems with integral equations. The method discretizes a surface in three dimensions, which reduces the dimension of the problem with one. A drawback of the method is that it yields a dense system of linear equations. This effectively prohibits the solution of large scale problems. Papers I-III describe the Fast Multipole Method. It reduces the cost of computing a dense matrix vector multiplication. This implies that large scale problems can be solved on personal computers. In Paper I the error introduced by the Fast Multipole Method is analyzed. Paper II and Paper III describe the implementation of the Fast Multipole Method. The problem of computing monostatic Radar Cross Section involves many right hand sides. Since the Fast Multipole Method computes a matrix times a vector, iterative techniques are used to solve the linear systems. It is important that the solution time for each system is as low as possible. Otherwise the total solution time becomes too large. Different techniques for reducing the work in the iterative solver are described in Paper IV-VI. Paper IV describes a block Quasi Minimal Residual method for several right hand sides and Sparse Approximate Inverse preconditioner that reduce the number of iterations significantly. In Paper V and Paper VI a method based on linear algebra called the Minimal Residual Interpolation method is described. It reduces the work in an iterative solver by accurately computing an initial guess for the iterative method. In Paper VII a hybrid method between Physical Optics and the Fast Multipole Method is described. It can handle large problems that are out of reach for the Fast Multipole Method.

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.

Robust preconditioners on block-triangular and block-factorized form for three types of linear systems of two-by-two block form are studied in this thesis. The first type of linear systems, which are dense, arise from a boundary element type of discretization of crack propagation problems. Numerical experiment show that simple algebraic preconditioning strategies results in iterative schemes that are highly competitive with a direct solution method. The second type of algebraic systems, which are sparse, indefinite and nonsymmetric, arise from a finite element (FE) discretization of the partial differential equations (PDE) that describe (visco)elastic glacial isostatic adjustment (GIA). The Schur complement approximation in the block preconditioners is constructed by assembly of local, exactly computed Schur matrices. The quality of the approximation is verified in numerical experiments. When the block preconditioners for the indefinite problem are combined with an inner iterative scheme preconditioned by a (nearly) optimal multilevel preconditioner, the resulting preconditioner is (nearly) optimal and robust with respect to problem size, material parameters, number of space dimensions, and coefficient jumps. Two approaches to mathematically formulate the PDEs for GIA are compared. In the first approach the equations are formulated in their full complexity, whereas in the second their formulation is confined to the features and restrictions of the employed FE package. Different solution methods for the algebraic problem are used in the two approaches. Analysis and numerical experiments reveal that the first strategy is more accurate and efficient than the latter. The block structure in the third type of algebraic systems is due to a fine-coarse splitting of the unknowns. The inverse of the pivot block is approximated by a sparse matrix which is assembled from local, exactly inverted matrices. Numerical experiments and analysis of the approximation show that it is robust with respect to problem size and coefficient jumps.

No description available.

No description available.