User:Compsonheir/Incomplete LU factorization

Given a matrix A, one can solve a linear system Ax = b by computing the LU factorization A = LU, then forward- and back-solving to obtain x. However, if A is sparse, the LU-factors can be prohibitively more dense. The incomplete LU factorization avoids this problem by computing sparse but inexact LU factors and instead using them to precondition an iterative method to solve the linear system.

Definition
The ordinary LU factorization of a matrix A yields a unit lower-triangular matrix L and an upper-triangular matrix U with the following entries:


 * $$ U_{ij} = A_{ij}-\sum_{ki $$
 * $$ L_{ij} = \left(A_{ij}-\sum_{kj $$

In the incomplete LU factorization, one instead chooses a set P of pairs (i,j) where we allow non-zero entries in either of the triangular factors. The last two equations defining L and U still hold but with the proviso that Lij = 0 if (i,j) is not in P, and likewise for U.

Of course, we could take P to be all pairs (i,j), in which case we would recover the usual LU-factorization. The ILU(0) factorization is obtained by taking P to be the set of all entries where A is not zero, or


 * $$ P = \{(i,j) | A_{ij} \neq 0\}.$$

The union of the sparsity patterns for L and U is the same as the sparsity pattern for A.

Other variants
The LU factors obtained from ILU(0) have the same sparsity pattern as A, but their product L*U is more dense. We can then define the ILU(1) factorization to be the static-pattern ILU obtained from the sparsity pattern of L(0)*U(0). Moreover, one can continue in this fashion recursively, defining ILU(p) as the ILU matrices found using the sparsity pattern of L(p-1)*U(p-1).

However, this approach depends solely on the non-zero structure of A with no concern whatsoever for the magnitude of the dropped terms.