Transcript ASI 3 Méthodes numériques pour l’ingénieur Résolution de systèmes linéaires par des méthodes directes : LDL’ et Choleski.

ASI 3
Méthodes numériques
pour l’ingénieur
Résolution de systèmes linéaires
par des méthodes directes :
LDL’ et Choleski
Résoudre un système linéaire
2 x1  4 x2
 x  3x
 1
2
Ax  b  
3 x1  x2

 x2
Fonction x = Gauss(A,b)
 2 x3
 x4
 6
 0
 x3
 2 x4
 8
 2 x3
 x4
 6
Fonction x = LU(A,b)
L,U = decompose(A)
U,c = descent(A,b)
y = triang(L,b)
x = triang(U,c)
x = triang(U,y)
Cas particulier : A est symétrique définie positive
Matrice à diagonale dominante
Définition :
une matrice carrée A est dite à diagonale dominante ssi :
n
i  1, n, aii   aij
j 1, j  i
4  3
7 2 0 
 6




 3 5 1  ,  4  2 0 
 0 5  6
3 0
1 



Théorème : Si A est a diagonale strictement dominante,
Alors est elle alors non singulière,
De plus Gauss est stable
et peut fonctionner sans changement de colonne
Éléments de démonstration : Ax=0, par l’absurde  x  0, A singulière
Matrice symétrique définie positive
Symétrique : A’=A
Définie positive : x,  Ax, x  0
( 0 : strictement d.f.)
Exemple :
2
2
min y avec y  Bx  min Bx  min  Bx,Bx 
x
x
x
 min x' ( B' B) x
x
 min x' Ax
avec A  B' B
x,  Ax , x   Bx  0
2
 Ax  0
éléments de démonstrat ion :
2
n  m
 nm
 m

   bij x j       bij x j   bij x j  
i 1  j 1
  i 1  j 1
 j 1

m n m
 
      bik bkj x j  xi 
 
 i 1  j 1k 1
x
Rappelez vous
des moindres carrés
Exemple
 2 -1 0 
la matrice - 1 2 - 1 est elle définie positive ?


 0 - 1 2 
 2 - 1 0   x1 
x' Ax  x1 x2 x3  - 1 2 - 1  x2 

 
 0 - 1 2   x3 
 x1
x2
2 x1 - x2
x3   - x1  2 x2

- x2


- x3 

 2 x3 
 2 x12  2 x1 x2  2 x22  2 x2 x3  2 x32
 x12   x1  x2 2   x2  x3 2  x32  0
Exemple
 2 -1 0 
la matrice - 1 2 - 1 est elle définie positive ?


 0 - 1 2 
 2 - 1 0   x1 
x' Ax  x1 x2 x3  - 1 2 - 1  x2 

 
 0 - 1 2   x3 
 x1
x2
2 x1 - x2
x3   - x1  2 x2

- x2


- x3 

 2 x3 
 2 x12  2 x1 x2  2 x22  2 x2 x3  2 x32
 x12   x1  x2 2   x2  x3 2  x32  0
Propriétés des matrices définies positives
Théorème :
Si
A est une matrice n x n strictement définie positive
Alors :
– A est non singulière
– aii > 0 pour i=1,n
 aij 2  aii a jj
i  j
 max a jk  max aii
k, j
i
A singulière  ker( A)  0
 x  0  ker( A)
or x  ker( A)  Ax  0 et donc x' Ax  0
ce qui est contradict oire avec l' hypothèse
Propriétés des matrices définies positives
Théorème :
Si
A est une matrice n x n strictement définie positive
Alors :
– A est non singulière
– aii > 0 pour i=1,n
 aij 2  aii a jj
i  j
 max a jk  max aii
k, j
i
 (i )
x, x' Ax  0, en particulier pour x


or x (i ) Ax (i )  aii  0 par hypothèse
0 

 
 1  i ème position

 
0
Propriétés des matrices définies positives
Théorème :
Si
A est une matrice n x n symétrique strictement définie positive
Alors :
 max a  max a
jk
k, j
i
ii
0 si i  j et i  j
 ( jk ) 
x, x' Ax  0, en particulier pour x
  1 si i  j
 - 1 si i  k
 ( jk )  ( jk )

or x Ax
 a jj  akk  a jk  akj  0
(1) A symétrique  2a jk  a jj  akk
 ( jk )  ( jk )
et x
Ax
 aii  akk  a jk  akj  0
(2) A symétrique  2a jk  a jj  akk
(1)et (2)  a jk 
a jj  akk
2
 max aii
i
et donc : max a jk  max aii
j ,k
i
Propriétés des matrices définies positives
Théorème :
Si
A est une matrice n x n symétrique strictement définie positive
Alors :
 a 2  a a
i  j
ij
ii
jj
0 si i  j et i  j


x, x' Ax  0, en particulier pour x ( jk )   1 si i  j
  si i  k
 ( jk )  ( jk )

or x Ax
 a jj   2 akk  2a jk  0
l' équation quadratiqu e n' adettant pas de racine réelle,
sont déterminan t est négatif
  4a 2jk  4a jj akk  0  a 2jk  a jj akk
Autres propriétés
Définition : une sous matrice principale d’une matrice A
est une matrice carrée de la forme A(1:i,1:i) quelque soit i
Théorème : Une matrice symétrique est définie positive
ssi chacune de ses sous matrice principales à un déterminant positif
 2 -1 0 
A  - 1 2 - 1


 0 - 1 2 
det A1  2  0;
 2  1
det A2  det 
 4 1  3  0

 1 2 
det A2  det A  2 det A2  1 det ....  4
Autres propriétés
Définition : une sous matrice principale d’une matrice A
est une matrice carrée de la forme A(1:i,1:i) quelque soit i
Théorème : Une matrice symétrique est définie positive
ssi chacune de ses sous matrice principales à un déterminant positif
 2 -1 0 
A  - 1 2 - 1


 0 - 1 2 
det A1  2  0;
 2  1
det A2  det 
 4 1  3  0

 1 2 
det A2  det A  2 det A2  1 det ....  4
Autres propriétés
Théorème :
Une matrice est symétrique définie positive si la méthode de Gauss
être appliquée sans permutations n’admet que des pivots positifs
De plus, le résultat est stable par rapport aux erreurs d’arrondi
Corollaire
si A est une matrice symétrique non singulière,
Alors il existe une matrice diagonale D
et une matrice triangulaire avec des 1 sur la diagonale L
telles que :
A = LDL’
éléments de démonstration : A non singulière => A=LU =LDV
et A’=V’DL’ que l’on identifie car A est symétrique
Factorisation LDL’
V=DL’
0
0
i
A
vj=S lijdj
aii=di+S lijvj
aij=dilij+S likvk
L
Exemple :
1
1   1
0
0  4 0 0  1  0.25 0.25 
4

 



0  0 4 0  0
1
0.75 
  1 4,25 2,75     0.25 1
 1 2,75 3,5   0.25 0.75 1  0 0 1  0

0
1

 



La factorisation LDL’
Fonction L,D = décomposeLDL(A)
pour i  1 jusqu' à n
somme  0
pour j  1 jusqu' à i  1
v j   ij d j
somme  somme   ij v j
fait
d i  aii  somme
pour j  i  1 jusqu' à n
somme  0
pour k  1 jusqu' à i  1
somme  somme   jk vk
fait
a ji  somme
l ji 
di
fait
fait
Choleski : LL’
A  LDL'
et
~~
A  LL'
~
LL D
D doit être positif !
Théorème :
toute matrice A symétrique définie positive
admet une décomposition unique sous la forme
A=LL’
ou L est une matrice triangulaire inférieure dont tous les éléments
diagonaux sont positifs
Choleski :
l’algorithme
Fonction L = Choleski(A)
 11  a11
pour j  2 jusqu' à n
a j1
 j1 
l11
fait
pour i  2 jusqu' à n  1
somme  0
pour k  1 jusqu' à i  1
somme  somme   ik 2
fait
 ii  aii  somme
pour j  i  1 jusqu' à n
somme  0
pour k  1 jusqu' à i  1
somme  somme   jk  ik
fait a  somme
ji
 ji 
lii
fait
 nn  ann  somme
Comparaison : temps de calcul
Gauss
Choleski
Additions
(n3-n)/3 n3/6
Multiplications
(n3-n)/3 n3/6
division
n(n-1)/2 n /2
32
n
racines
Total
2n3/3
n3/3
Plus stable !
Logiciels
– Cas général : PA=LU
– Matrice symétrique définie positive : A=LL’ (Choleski)
– Matrice symétrique : A=LDL’
– Matrice tridiagonale (heisenberg) : LU par bande (cf TD)
– Matrice triangulaire : « remontée en n2 »
Matlab : x =A\b ;
si A triangulaire : x=trisup(A,b)
sinon si A symétrique définie positive ; L=chol(A)
sinon : (* cas général*)
[L,U,P]=lu(A); z=L\(P*b); x=U\z;
LAPACK = BLAS (basic linear algebra subprograms)
- blaise, octave, matlab (interprétés calcul)
- scilab, maple (interprétés formels)
- IMSL, NAG (bibliothèques)