140
5 Systèmes linéaires
3. en fixant la valeur 1 pour les n éléments diagonaux de L, (5.12)
devient un système déterminé qui peut être résolu avec l’algorithme
de Gauss : posons A
(1) = A i.e. a
(1)
ij = a ij pour i, j = 1, . . ., n ;
pour k = 1, . . ., n − 1
pour i = k + 1, . . ., n
l ik =
a
(k)
ik
a
(k)
kk
,
pour j = k + 1, . . . , n
a
(k+1)
ij
= a
(k)
ij − l ik a
(k)
kj
(5.13)
Les termes a
(k)
kk , appelés pivots, doivent être tous non nuls. Pour
k = 1, . . . , n − 1 la matrice A
(k+1) = (a
(k+1)
ij
) a n − k lignes et colonnes.
Remarque 5.1 Il n’est pas nécessaire de stocker toutes les matrices A
(k)
dans l’algorithme (5.13) ; on peut en effet écraser les (n − k) × (n − k) derniers
éléments de la matrice originale A avec les (n−k)×(n−k) éléments de A
(k+1) .
De plus, puisqu’à l’étape k, les éléments sous-diagonaux de la k-ème colonne
n’ont aucun impact sur la matrice finale U, ils peuvent être remplacés par les
coefficients de la k-ème colonne de L (les multiplicateurs). C’est ce qui est fait
dans le Programme 5.1. A l’étape k de l’algorithme, les éléments stockés à la
place des coefficients originaux de A sont
⎡
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎣
a
(1)
11 a
(1)
12 . . .
l21 a
(2)
22
. . .
. . .
. . .
. . . . . . a
(1)
1n
a
(2)
2n
. . .
lk1 . . . lk,k−1
. . .
. . .
ln1 . . . ln,k−1
a
(k)
kk . . . a
(k)
kn
. . .
. . .
a
(k)
nk . . . a
(k)
nn
⎤
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎦
,
où la matrice encadrée est A
(k) .
Ainsi, on peut implémenter l’algorithme en ne stockant qu’un seule matrice,
qu’on initialise avec A et qu’on modifie à chaque itération k ≥ 2 en écrasant les
nouveaux termes a
(k)
ij , pour i, j ≥ k + 1, ainsi que les multiplicateurs lik, pour
i ≥ k + 1. Noter qu’il n’est pas nécessaire de stocker les éléments diagonaux
lii puisqu’on sait qu’ils valent tous 1.
A la fin de cette procédure, les éléments de la matrice triangulaire U
sont donnés par u ij = a
(i)
ij pour i = 1, . . . , n et j = i, . . . , n, tandis que
ceux de L sont donnés par les coefficients l ij calculés par l’algorithme.
Dans (5.13), les termes diagonaux de L ne sont pas considérés, puisque
leur valeur a été fixée à 1.
5 Systèmes linéaires
3. en fixant la valeur 1 pour les n éléments diagonaux de L, (5.12)
devient un système déterminé qui peut être résolu avec l’algorithme
de Gauss : posons A
(1) = A i.e. a
(1)
ij = a ij pour i, j = 1, . . ., n ;
pour k = 1, . . ., n − 1
pour i = k + 1, . . ., n
l ik =
a
(k)
ik
a
(k)
kk
,
pour j = k + 1, . . . , n
a
(k+1)
ij
= a
(k)
ij − l ik a
(k)
kj
(5.13)
Les termes a
(k)
kk , appelés pivots, doivent être tous non nuls. Pour
k = 1, . . . , n − 1 la matrice A
(k+1) = (a
(k+1)
ij
) a n − k lignes et colonnes.
Remarque 5.1 Il n’est pas nécessaire de stocker toutes les matrices A
(k)
dans l’algorithme (5.13) ; on peut en effet écraser les (n − k) × (n − k) derniers
éléments de la matrice originale A avec les (n−k)×(n−k) éléments de A
(k+1) .
De plus, puisqu’à l’étape k, les éléments sous-diagonaux de la k-ème colonne
n’ont aucun impact sur la matrice finale U, ils peuvent être remplacés par les
coefficients de la k-ème colonne de L (les multiplicateurs). C’est ce qui est fait
dans le Programme 5.1. A l’étape k de l’algorithme, les éléments stockés à la
place des coefficients originaux de A sont
⎡
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎣
a
(1)
11 a
(1)
12 . . .
l21 a
(2)
22
. . .
. . .
. . .
. . . . . . a
(1)
1n
a
(2)
2n
. . .
lk1 . . . lk,k−1
. . .
. . .
ln1 . . . ln,k−1
a
(k)
kk . . . a
(k)
kn
. . .
. . .
a
(k)
nk . . . a
(k)
nn
⎤
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎦
,
où la matrice encadrée est A
(k) .
Ainsi, on peut implémenter l’algorithme en ne stockant qu’un seule matrice,
qu’on initialise avec A et qu’on modifie à chaque itération k ≥ 2 en écrasant les
nouveaux termes a
(k)
ij , pour i, j ≥ k + 1, ainsi que les multiplicateurs lik, pour
i ≥ k + 1. Noter qu’il n’est pas nécessaire de stocker les éléments diagonaux
lii puisqu’on sait qu’ils valent tous 1.
A la fin de cette procédure, les éléments de la matrice triangulaire U
sont donnés par u ij = a
(i)
ij pour i = 1, . . . , n et j = i, . . . , n, tandis que
ceux de L sont donnés par les coefficients l ij calculés par l’algorithme.
Dans (5.13), les termes diagonaux de L ne sont pas considérés, puisque
leur valeur a été fixée à 1.