` ILUFactor`
-- incomplete factorization of quadratic matrices

#include <laspack/factor.h> QMatrix *ILUFactor(QMatrix *Q);

Q = (D + L) D^{-1} (D + U) + R

is performed, where ` D`,
` U`,
and
` L`
are certain diagonal, upper,
and lower triangular matrices, respectively.
The remainder matrix ` R`
contains fill elements,
which have been ignored during the factorization process.
For symmetric matrices, the incomplete Cholesky factorization
with ` L = U^T`
is applied.
In order to be able to perform the incomplete factorization
for singular matrices too,
regularization be means of increasing of diagonal entries is applied.

In the current implementation, matrices ` L`
and
` U`
are laid down
by the same position of non-zero elements as in matrix ` Q.`

Matrices are stored ` D`,
` U`,
and
` L`
in connection with ` Q`
as matrix ` D + L + U`
of type ` QMatrix`.
This one is also the one returned by the procedure ` ILUFactor`.
How to extract the particular matrices from it,
is shown in the example below.

The incomplete factorization is comprehensively describe and analyzed e.g. in:

W. Hackbusch: Iterative Solution of Large Sparse Systems of Equations, Springer-Verlag, Berlin, 1994.

` factor.h ... ` header file

` factor.c ... ` source file

The following example shown the usage of the procedure ` ILUFactor`
in implementation of an ILU preconditioner.
This have to solve the system of equations

W y = c

which arises during the solution of preconditioned systems

W^{-1} A x = W^{-1} b

at which

W = (D + L) D^{-1} (D + U).

The corresponding ` LASPack`
routine could be build as follows:

Vector *ILUPrecond(QMatrix *A, Vector *y, Vector *c, double Omega) { Q_Lock(A); V_Lock(y); V_Lock(c); Asgn_VV(y, MulInv_QV(Add_QQ(Diag_Q(ILUFactor(A)), Upper_Q(ILUFactor(A))), Mul_QV(Diag_Q(ILUFactor(A)), MulInv_QV(Add_QQ(Diag_Q(ILUFactor(A)), Lower_Q(ILUFactor(A))), c)))); Q_Unlock(A); V_Unlock(y); V_Unlock(c); return(y); }

qmatrix(3LAS), operats(3LAS), errhandl(3LAS)

In the current implementation,
it is assumed that the matrix ` Q`
has
symmetric structure with regard to non-zero elements.
Because during discretization of differential equations even such matrices
arise,
this restriction is for many applications not grave.
For matrices which do not have suffice to the above condition,
the error ` LASILUStructErr` is raised.

Tomas Skalicky (skalicky@msmfs1.mw.tu-dresden.de)