Derivation of Krylov space methods (Conjugate Gradients, GMRES, BiCGstab, et cetera) has long been a dark art. The road from a simple idea to a practical code is long and arduous, and uniqueness and optimality of the endpoint is not guaranteed. We try to remedy this situation by offering a formal treatment that requires a researcher to specify a very high-level predicate describing the operation. After this, the actual derivation will be mechanical and in such a way that the result is proved correct.
In this project we derive Krylov methods such as the Conjugate Gradients method in the FLAME framework. With this, we have two goals in mind. First, we show how FLAME can be used to simplify the process of deriving iterative methods. Second, we make a case that a Flame-based environment for deriving new iterative methods is a distinct, and attractive, possibility.
Traditional expositions of the \ac{CG} method, and ones related to it, posit the basic form of relations between matrices and vectors, and compute the scalar coefficients in them by `lengthy induction arguments'. In the spirit of Householder's derivation, we summarize vector and scalar sequences as matrices; the FLAME framework then allows us to derive in a formal manner the actual iteration from properties on the quantities constructed in the algorithm. The big advantage here is that we can dispense with quantified statements over sequences, and instead consider simple predicates over simple, unindexed, objects. Simultaneously, the inductive arguments that have always been at the heart of traditional expositions are captured in a framework that guides the derivation of the algorithm and the proof of its correctness. As a first demonstration of the potential of this approach, we have derived a CG method for nonsymmetric systems.
The conciseness of our derivation should be contrasted with the lengthy research papers in the classical approach to polynomial iterative methods. Beyond simply presenting an alternative derivation of these methods, we argue that the essential calculations in a FLAME worksheet, contained in the update step, can be derived mechanically from the loop invariant of the algorithm. Coupling this to ongoing projects for automatic code generation from FLAME worksheets, this raises the possibility of automatic generation of numerical libraries. Krylov subspace methods are then merely a proof-of-concept of a much more general idea: the mechanical derivation of algorithms and tuned library software incorporating these algorithms.
Victor Eijkhout Research Scientist
National Science Foundation, NSF Grant 0904907
