<< spilu_ilu0 Spilu spilu_iludp >>

Spilu >> Spilu > spilu_ilud

spilu_ilud

sparse ILU factorization with single dropping and diagonal compensation

Calling Sequence

[L,U]=spilu_ilud(A)
[L,U]=spilu_ilud(A,alph)
[L,U]=spilu_ilud(A,alph,drop)

Parameters

A

a n-by-n sparse real matrix of doubles

alph

a 1-by-1 matrix of doubles, positive, diagonal compensation parameter (default: 0.5)

drop

a 1-by-1 matrix of doubles, positive, threshold parameter for dropping small terms in the factorization (default: drop=0.001*max(abs(A)))

L

a n-by-n sparse real matrix of doubles, lower triangular

U

a n-by-n sparse real matrix of doubles, upper triangular

Description

Builds an incomplete LU factorization of the sparse matrix A, that is, computes a lower triangular matrix L and an upper triangular matrix U such that

        A ≈ L*U
      

All diagonal elements of the input matrix must be nonzero.

Any optional input argument equal to the empty matrix [] is replaced by its default value.

This routine computes the ILU factorization with standard threshold dropping.

There is no control on memory size required for the factors as is done in ILUT. This routines computes also various diagonal compensation ILU's such MILU. These are defined through the parameter alph.

The parameter alph is the diagonal compensation parameter. The term

alph*(sum of all dropped out elements in a given row)
      

is added to the diagonal element of U of the factorization. Thus:

The drop parameter is the threshold parameter for dropping small terms in the factorization. During the elimination, a term A(i,j) is dropped whenever

drop * (weighted norm of A(i,:)) > abs(A(i,j))
      

Here weighted norm = 1-norm / (number of nnz elements in the row).

Hence, a smaller value of drop (e.g. drop=0.) leads to a better decomposition of the matrix, but also leads to a larger number of nonzeros in the matrices L and U.

Examples

// Incomplete factorization of a 5-by-5 sparse matrix
// nnz(A)=11
A=[
    18.    0.    -6.     0.     3.  
    0.     16.    0.     8.     0.  
   -6.     0.     36.    0.     0.  
    0.     8.     0.     16.    0.  
    3.     0.     0.     0.     6.  
];
A=sparse(A);
// With default options
[L,U]=spilu_ilud(A); 
// Check the quality of the decomposition
[norm(A-L*U,"inf")/norm(A,"inf") nnz(L)+nnz(U)]
// Configure alpha
[L,U]=spilu_ilud(A,1.); 
[norm(A-L*U,"inf")/norm(A,"inf") nnz(L)+nnz(U)]
// Configure alpha and drop
[L,U]=spilu_ilud(A,1.,0.); 
[norm(A-L*U,"inf")/norm(A,"inf") nnz(L)+nnz(U)]

// Incomplete factorization of a 960-by-960 sparse matrix
// nnz(A)=15844
path = spilu_getpath (  );
filename = fullfile(path,"tests","matrices","nos3.mtx");
A=mmread(filename);
n=size(A,1);
b=ones(n,1);
alph=0;
drop=0;
[L,U]=spilu_ilud(A,alph,drop);
x=U\(L\b);
norm(A*x-b)

Authors

Copyright (C) 2011 - DIGITEO - Michael Baudin

Copyright (C) 2011 - NII - Benoit Goepfert

Copyright (C) 2005 - INRIA - Sage Group (IRISA)

Copyright (C) 1993 - Univ. of Tennessee and Oak Ridge National Laboratory - Youcef Saad


Report an issue
<< spilu_ilu0 Spilu spilu_iludp >>