134
4 M´ ethodes it´ eratives pour la r´ esolution des syst` emes lin´ eaires
2, . . ., n, on effectue les op´ erations suivantes : si lev ik ≤ p, k = 1, . . ., i−1,
le coefficient m ik de L in et les coefficients a
(k+1)
ij
de U in , j = i + 1, . . ., n,
sont mis ` a jour. De plus, si a
(k+1)
ij
= 0 la valeur lev ij est choisie comme
le minimum entre l’ancienne valeur de lev ij et lev ik + lev kj + 1. La raison
de ce choix est que |a
(k+1)
ij
| = |a
(k)
ij − m ik a
(k)
kj | | |δ
levij
− δ
levik+levkj +1
|,
et qu’on peut donc supposer que |a
(k+1)
ij
| est le maximum entre δ
levij et
δ
levik+levkj +1 .
La proc´ edure de factorisation que l’on vient de d´ ecrire s’appelle ILU(p)
et se r´ ev` ele extr` emement efficace (pour p petit) quand on la couple avec
une renum´ erotation convenable des lignes et des colonnes de la matrice.
Le Programme 18 propose une impl´ ementation de la factorisation ILU(p) ;
il renvoie en sortie les matrices approch´ ees L in et U in (stock´ ees dans la
matrice initiale a), avec des 1 sur la diagonale de L in , et la matrice lev
contenant le niveau de remplissage de chaque coefficient ` a la fin de la
factorisation.
Programme 18 - ilup : Factorisation incompl` ete ILU(p)
function [A,lev]=ilup(A,p)
%ILUP Factorisation incompl` ete ILU(p).
% [Y,LEV]=ILUP(A): U est stock´ ee dans la partie triangulaire sup´ erieure
% de Y et L dans la partie triangulaire inf´ erieure stricte de Y.
% Les matrices L et U ont un niveau de remplissage P.
% LEV contient le niveau de remplissage de chaque terme `
a la fin
% de la factorisation.
[n,m]=size(A);
if n ˜= m, error(’Seulement pour les matrices carr´ ees’); end
lev=Inf*ones(n,n);
i=(A˜=0);
lev(i)=0,
for i=2:n
for k=1:i-1
if lev(i,k) <= p
if A(k,k)==0, error(’Pivot nul’); end
A(i,k)=A(i,k)/A(k,k);
for j=k+1:n
A(i,j)=A(i,j)-A(i,k)*A(k,j);
if A(i,j) ˜= 0
lev(i,j)=min(lev(i,j),lev(i,k)+lev(k,j)+1);
end
end
end
end
for j=1:n, if lev(i,j) > p, A(i,j) = 0; end, end
end
return
4 M´ ethodes it´ eratives pour la r´ esolution des syst` emes lin´ eaires
2, . . ., n, on effectue les op´ erations suivantes : si lev ik ≤ p, k = 1, . . ., i−1,
le coefficient m ik de L in et les coefficients a
(k+1)
ij
de U in , j = i + 1, . . ., n,
sont mis ` a jour. De plus, si a
(k+1)
ij
= 0 la valeur lev ij est choisie comme
le minimum entre l’ancienne valeur de lev ij et lev ik + lev kj + 1. La raison
de ce choix est que |a
(k+1)
ij
| = |a
(k)
ij − m ik a
(k)
kj | | |δ
levij
− δ
levik+levkj +1
|,
et qu’on peut donc supposer que |a
(k+1)
ij
| est le maximum entre δ
levij et
δ
levik+levkj +1 .
La proc´ edure de factorisation que l’on vient de d´ ecrire s’appelle ILU(p)
et se r´ ev` ele extr` emement efficace (pour p petit) quand on la couple avec
une renum´ erotation convenable des lignes et des colonnes de la matrice.
Le Programme 18 propose une impl´ ementation de la factorisation ILU(p) ;
il renvoie en sortie les matrices approch´ ees L in et U in (stock´ ees dans la
matrice initiale a), avec des 1 sur la diagonale de L in , et la matrice lev
contenant le niveau de remplissage de chaque coefficient ` a la fin de la
factorisation.
Programme 18 - ilup : Factorisation incompl` ete ILU(p)
function [A,lev]=ilup(A,p)
%ILUP Factorisation incompl` ete ILU(p).
% [Y,LEV]=ILUP(A): U est stock´ ee dans la partie triangulaire sup´ erieure
% de Y et L dans la partie triangulaire inf´ erieure stricte de Y.
% Les matrices L et U ont un niveau de remplissage P.
% LEV contient le niveau de remplissage de chaque terme `
a la fin
% de la factorisation.
[n,m]=size(A);
if n ˜= m, error(’Seulement pour les matrices carr´ ees’); end
lev=Inf*ones(n,n);
i=(A˜=0);
lev(i)=0,
for i=2:n
for k=1:i-1
if lev(i,k) <= p
if A(k,k)==0, error(’Pivot nul’); end
A(i,k)=A(i,k)/A(k,k);
for j=k+1:n
A(i,j)=A(i,j)-A(i,k)*A(k,j);
if A(i,j) ˜= 0
lev(i,j)=min(lev(i,j),lev(i,k)+lev(k,j)+1);
end
end
end
end
for j=1:n, if lev(i,j) > p, A(i,j) = 0; end, end
end
return