The Multigrid Workbench: Linear Iterations

Iterative methods for a linear system of equations
Au = f, (1)
like the Gauß-Seidel method, over-relaxation (SOR), or the Gauß-Jacobi method can all be written in the general form
Buk+1 = f - ( A - B ) uk (2)
uk+1 = uk + B-1 ( f - Auk ) (3)
where the upper index indicates the number of the iterate. B is a matrix characteristic for the method.

Clearly, the exact solution u*=A-1f is a fixed point of iteration (2). To check whether the iteration converges uk=u*+ek and uk+1=u*+ek+1 may be inserted into (3). This simplifies to
ek+1 = ( I - B-1A ) ek. (4)
Convergence is obtained, when the error norm is reduced. Thus
||ek+1 || || I - B-1A || || ek ||. (5)
assures convergence when || I - B-1A || < 1. This statement is a precise way of saying that for (3) to converge fast, B must approximate A. Clearly, when B=A this norm vanishes and we trivially get a method which converges in one step. However, this iteration would require the inversion of B, or the solution of a system, where B is the system matrix. This would essentially mean the original system is solved directly, and nothing would be won.

The difficulty is to find the right compromise of having B easily invertible, and being a good approximation to A.

Relaxation methods

Relaxation methods use particularly simple choices of B. For the Gauß-Jacobi method B=diag(A), for the Gauß-Seidel method, B is taken as the lower triangle of A. These B are certainly easily invertible, however, they do not make very good approximations to A, at least when A is large. However, they provide good smoothers for the multigrid method.

Coarse grid correction

Here we now study an alternative method, the so-called coarse grid correction method.

Our idea is to construct B-1 form a similar grid problem as A, however, a simpler one by taking fewer gridlines. Let us assume that AC is constructed just like A, but from a grid with only (N-1)x(N-1) gridlines. The dimension of the matrix is therefore less than ¼ of the original dimension. It is therefore obviously easier to invert AC than A. Unfortunately, however, equation (3) cannot be used any more, because now AC-1 does not have a dimension compatible with f, A, and u. To fix this, we introduce the interpolation (prolongation) operator P, and the restriction operator R. R and P are mappings between these two spaces:
R: |RNxN -> |R(N-1)/2x(N-1)/2 (6)
P: |R(N-1)/2x(N-1)/2 -> |RNxN (7)
and define
BC = P AC-1 R (8)

Now we replace B-1 in (3) by BC. Note that since BC has rank smaller than A, therefore there exists no formulation of the iteration equivalent to (2).

For the same reason, the product BCA is singular, and
|| I - BCA || 1.
Thus the coarse grid correction cannot be a convergent method. However, it is a good iteration with fast convergence for the low frequency solution modes.

Multigrid method

The multigrid method alternates between smoothers and coarse grid corrections. The combination of the two results in a new iterative algorithm, where the different convergence properties of the two basic iteration combine favorably. The two grid algorithm resulting from this idea is extended to a multigrid algorithm by extending this idea recursively to the solution of A-1.

Other iterative methods

Besides the simple iterative schemes many more possibilities exist. We refer to the chapter on numerical linear algebra in the Computational Science Education Project.

Ulrich Ruede , Thu Feb 2 21:06:09 MEZ 1995
Updated by Craig C. Douglas