In this book, which focuses on the use of iterative methods for solving large sparse systems of linear equations, templates are introduced to meet the needs of both the traditional user and the high-performance specialist Templates, a description of a general algorithm rather than the executable object or source code more commonly found in a conventional software library, offer whatever degree of customization the user may desire

/pdf/templates-for-the-solution-of-linear-systems-building-blocks-28xlh68wjc.pdf

Templates for the Solution of Linear Systems: Building Blocks for Iterative Methods

In recent years, the SUite of Nonlinear and DIfferential/ALgebraic equation Solvers (SUNDIALS) has been redesigned to better enable the use of application-specific and third-party algebraic solvers and data structures. Throughout this work, we have adhered to specific guiding principles that minimized the impact to current users while providing maximum flexibility for later evolution of solvers and data structures. The redesign was done through creation of new classes for linear and nonlinear solvers, enhancements to the vector class, and the creation of modern Fortran interfaces that leverage interoperability features of the Fortran 2003 standard. The vast majority of this work has been performed "behind-the-scenes," with minimal changes to the user interface and no reduction in solver capabilities or performance. However, these changes now allow advanced users to create highly customized solvers that exploit their problem structure, enabling SUNDIALS use on extreme-scale, heterogeneous computational architectures.

Enabling New Flexibility in the SUNDIALS Suite of Nonlinear and Differential/Algebraic Equation Solvers.

http://www.inf.ufes.br/~luciac/comcie/JFNK-Knoll.pdf

Jacobian-free Newton-Krylov methods: a survey of approaches and applications

http://www.mathcs.emory.edu/%7Ebenzi/Web_papers/survey.pdf

Preconditioning techniques for large linear systems: a survey

The standard model for sea ice dynamics treats the ice pack as a visco‐plastic material that flows plastically under typical stress conditions but behaves as a linear viscous fluid where strain rates are small and the ice becomes nearly rigid. Because of large viscosities in these regions, implicit numerical methods are necessary for time steps larger than a few seconds. Current solution methods for these equations use iterative relaxation methods, which are time consuming, scale poorly with mesh resolution, and are not well adapted to parallel computation. To remedy this, the authors developed and tested two separate methods. First, by demonstrating that the viscous‐plastic rheology can be represented by a symmetric, negative definite matrix operator, the much faster and better behaved preconditioned conjugate gradient method was implemented. Second, realizing that only the response of the ice on timescales associated with wind forcing need be accurately resolved, the model was modified so that it reduces to the viscous‐plastic model at these timescales, whereas at shorter timescales the adjustment process takes place by a numerically more efficient elastic wave mechanism. This modification leads to a fully explicit numerical scheme that further improves the model’s computational efficiency and is a great advantage for implementations on parallel machines. Furthermore, it is observed that the standard viscous‐plastic model has poor dynamic response to forcing on a daily timescale, given the standard time step (1 day) used by the ice modeling community. In contrast, the explicit discretization of the elastic wave mechanism allows the elastic‐viscous‐plastic model to capture the ice response to variations in the imposed stress more accurately. Thus, the elastic‐viscous‐plastic model provides more accurate results for shorter timescales associated with physical forcing, reproduces viscous‐plastic model behavior on longer timescales, and is computationally more efficient overall.

/pdf/an-elastic-viscous-plastic-model-for-sea-ice-dynamics-3dfsoyzprj.pdf

An Elastic–Viscous–Plastic Model for Sea Ice Dynamics

The numerical simulation of groundwater flow through heterogeneous porous media is discussed. The focus is on the performance of a parallel multigrid preconditioner for accelerating convergence of conjugate gradients, which is used to compute the pressure head. The numerical investigation considers the effects of boundary conditions, coarse grid solver strategy, increasing the grid resolution, enlarging the domain, and varying the geostatistical parameters used to define the subsurface realization. Scalability is also examined. The results were obtained using the PARFLOW groundwater flow simulator on the CRAY T3D massively parallel computer.

A parallel multigrid preconditioned conjugate gradient algorithm for groundwater flow simulations

In many parts of the world, groundwater resources are under increasing threat from growing demands, wasteful use, and contamination. To face the challenge, good planning and management practices are needed. A key to the management of groundwater is the ability to model the movement of fluids and contaminants. The purpose of this book is to construct conceptual and mathematical models that can provide the information required for making decisions associated with the management of groundwater resources, and the remediation of contaminated aquifers. The basic approach of this book is to accurately describe the underlying physics of groundwater flow and solute transport in heterogeneous porous media, starting at the microscopic level, and to rigorously derive their mathematical representations at the macroscopic levels.

Modeling groundwater flow and contaminant transport

The conjugate gradient method of Hestenes and Stiefel is an effective method for solving large, sparse Hermitian positive definite (hpd) systems of linear equations, $Ax = b$. Generalizations to non-hpd matrices have long been sought. The recent theory of Faber and Manteuffel gives necessary and sufficient conditions for the existence of a CG method. This paper uses these conditions to develop and organize such methods. It is shown that any CG method for $Ax = b$ is characterized by an hpd inner product matrix B and a left preconditioning matrix C. At each step the method minimizes the B-norm of the error over a Krylov subspace. This characterization is then used to classify known and new methods. Finally, it is shown how eigenvalue estimates may be obtained from the iteration parameters, generalizing the well-known connection between CG and Lanczos. Such estimates allow implementation of a stopping criterion based more nearly on the true error.

A taxonomy for conjugate gradient methods

Field-scale modeling of three-dimensional groundwater flow in randomly heterogeneous porous media is considered. The system of linear equations resulting from the finite difference discretization of this problem may easily involve more than a million unknowns. Problems of such magnitude can only be solved practically on supercomputers and require an efficient iterative solution method that is well suited to the particular computer architecture being used. The preconditioned conjugate gradient method is highly efficient for this groundwater flow problem, and by using an appropriate preconditioning matrix the method can be adapted to different computers. The numerical software package CgCode enhances this adaptability; the user of CgCode must only provide subroutines for preconditioning and matrix-vector multiplication. These subroutines should be tailored to the specific problem and machine architecture. Numerical results demonstrating the efficiency of polynomial preconditioning for the conjugate gradient method on both a vector machine (Cray X-MP/48) and a vector-parallel machine (Alliant FX/8) are presented.

A numerical investigation of the conjugate gradient method as applied to three‐dimensional groundwater flow problems in randomly heterogeneous porous media

The solution of a linear system of equations, Ax = b, arises in many scientific applications. If A is large and sparse, an iterative method is required. When A is hermitian positive definite (hpd), the conjugate gradient method of Hestenes and Stiefel is popular. When A is hermitian indefinite (hid), the conjugate residual method may be used. If A is ill-conditioned, these methods may converge slowly, in which case a preconditioner is needed. In this thesis we examine the use of polynomial preconditioning in CG methods for both hermitian positive definite and indefinite matrices. Such preconditioners are easy to employ and well-suited to vector and/or parallel architectures.
We first show that any CG method is characterized by three matrices: an hpd inner product matrix B, a preconditioning matrix C, and the hermitian matrix A. The resulting method, CG(B,C,A), minimizes the B-norm of the error over a Krylov subspace. We next exploit the versatility of polynomial preconditioners to design several new CG methods. To obtain an optimum preconditioner, we solve a constrained minimax approximation problem. The preconditioning polynomial, C($\lambda$), is optimum in that it minimizes a bound on the condition number of the preconditioned matrix, $p\sb{m}$(A). An adaptive procedure for dynamically determining the optimum preconditioner from the CG iteration parameters is also discussed. Finally, in a variety of numerical experiments, conducted on a Cray X-MP/48, we demonstrate the effectiveness of polynomial preconditioning.

S. F. Ashby

Papers

A parallel multigrid preconditioned conjugate gradient algorithm for groundwater flow simulations

Modeling groundwater flow and contaminant transport

A taxonomy for conjugate gradient methods

A numerical investigation of the conjugate gradient method as applied to three‐dimensional groundwater flow problems in randomly heterogeneous porous media

Polynomial preconditioning for conjugate gradient methods