Krylov subspace iteration

In the simulation of continuous events, such as the ow of a uid through a pipe, or the ow of air around an aircraft, one usually imposes a grid over the area of interest and one restricts oneself to the computation of relevant parameters, for instance the pressure or the velocity of the ow or the temperature, in the gridpoints. Physical laws lead to approximate relations between these parameters in neighbouring gridpoints, and together with the prescribed behavior at the boundary gridpoints and with given sources, this leads eventually to very large linear systems of equations Ax = b. The vector x stands for the unknown values of the parameter(s) in the gridpoints, b represents the given input, and the matrix A describes the relations between parameters in the gridpoints. Because these relations are often restricted to nearby gridpoints, most of the elements of the matrix are zero. The model becomes more accurate when the grid is reened, that is when the distance between gridpoints is decreased, and in 3-dimensional simulation this easily leads to very large systems of equations. Even only a few hundreds gridpoints in each coordinate direction leads already to systems with millions of unknowns. Many other problems also lead to very large systems: electric circuit simulation, the computation of magnetic elds, weather prediction, chemical processes, semiconductor device simulation, nuclear reactor safety problems, stress in mechanical structures, etcetera. The standard numerical solution methods for these linear systems are based on clever implementations of Gaussian elimination, where the sparse structure of the linear systems is exploited as much as possible in order to avoid computations with zero elements and in order to avoid storage of zero elements. But for very large systems these methods are often too expensive even on todays fastest supercomputers, except for cases where A has a very special structure. One can show mathematically, that for many of the problems listed above, the standard solution methods will not lead to solutions in any reasonable amount of time. For that reason, it has since long been tried to approximate the solution x in some iterative manner. One starts with a good guess for the solution, for instance by solving a much easier nearby (idealized) problem and one attempts to improve this guess by reducing the error with a convenient cheap approximation for A. This can be done in a repetitive way and in that case we …