ENSEM 1A-ISN Analyse Numérique

Travaux Dirigés d’Analyse Numérique

TD 07 : Résolution des systèmes linéaires - corrigé

Généralités

Exercice 1. Normes vectorielles

Sur Rn (n ∈ N∗ ), on définit les normes

n n
! 12
X X
∥x∥1 = |xk | (norme 1), ∥x∥2 = |xk |2 (norme 2), ∥x∥∞ = max |xk | (norme infinie)
k∈J1,nK
k=1 k=1

De manière générale, étant donné p ∈ N∗ , on définit la norme p d’un vecteur de Rn par

n
! p1
X
p
∥x∥p = |xk | .
k=1

Montrer que pour tout vecteur x ∈ Rn on a :

lim ∥x∥p = ∥x∥∞

p→+∞

ce qui justifie le terme ”norme infinie”.

Correction 1.
Soit x ∈ Rn . Il existe alors i ∈ J1, nK (pas forcément unique) tel que

|xi | = max |xk | = ∥x∥∞ .

k∈J1,nK

Nous avons alors

n
X
|xi |p ≤ |xk |p ≤ n|xi |p , ∀ p ∈ N∗ ,
k=1

d’où 1
∥x∥∞ ≤ ∥x∥p ≤ n p ∥x∥∞ , ∀ p ∈ N∗ .
1
Maintenant, comme lim n p = 1, le théorème des encadrements assure que
p→∞

lim ∥x∥p = ∥x∥∞ .

p→+∞

Exercice 2. Normes matricielles

On se place dans l’espace vectoriel des matrices à coefficients réels (matrices carrées ou rectangulaires)

1
ENSEM 1A-ISN Analyse Numérique

1. À partir des normes vectorielles 1 et infinie, on définit les normes matricielles subordonnées 1 et infinie par

∥Ax∥1 ∥Ax∥∞
∥A∥1 = sup , ∥A∥∞ = sup .
x ∈ Cn ∥x∥1 x ∈ Cn ∥x∥∞
x ̸= 0 x ̸= 0

Montrer que X X
∥A∥1 = max |ai,j |, ∥A∥∞ = max |ai,j |,
j i
i j

2. Application : calculer la norme infinie (et/ou 1 !) de la matrice A suivante ainsi que celle de son inverse
   
10 7 8 7 25 −41 10 −6
 7 5 6 5  −1
 −41 68 −17 10 
A=  8 6 10 9  A =  10 −17
  
5 −3 
7 5 9 10 −6 10 −3 2

Correction 2.
Pour éviter d’introduire trop de notations la correction, est faite dans le cas d’une matrice A = (ai,j ) carrée
d’ordre n. Vous pouvez essayer de reproduire ce corrigé pour des matrices à m lignes et n colonnes.
X
1. (a) ∥A∥1 = max |ai,j |
j
i
Soit x ∈ Rn . Nous avons alors
 
n X
 n n X
n n n
X  X X X

∥Ax∥1 = 
 ai,j j  ≤
x  |ai,j ||xj | ≤ |xj | |ai,j |
i=1  j=1  i=1 j=1 j=1 i=1
 
Xn n
X n
X
≤  |xj | max |ai,j | ≤ ∥x∥1 max |ai,j |.
j∈J1,nK j∈J1,nK
j=1 i=1 i=1

Par conséquent
n
X
∥A∥1 ≤ max |ai,j |.
j∈J1,nK
i=1

Maintenant, il existe p ∈ J1, nK tel que

n
X n
X
|ai,p | = max |ai,j |.
j∈J1,nK
i=1 i=1

Soit y le vecteur de Rn défini par yi = 0 si i ̸= p et yp = 1. Par construction, ∥y∥1 = 1. De plus,

(Ay)i = ai,p , ∀ i ∈ J1, nK,

par conséquent
n
X n
X n
X
∥Ay∥1 = |(Ay)i | = |ai,p | = max |ai,j |.
j∈J1,nK
i=1 i=1 i=1

On a donc égalité pour ce vecteur.

2
ENSEM 1A-ISN Analyse Numérique

X
(b) ∥A∥∞ = max |ai,j |
i
j
Soit x ∈ Rn . Nous avons alors pour tout i ∈ J1, nK
     
 n n n n
X  X X X

|(Ax)i | =  ai,j xj  ≤ |ai,j ||xj | ≤  |ai,j | ∥x∥∞ ≤  max |ai,j | ∥x∥∞ .
 j=1  j=1 i∈J1,nK
j=1 j=1

D’où  
n
X
∥Ax∥∞ = max |(Ax)i | ≤  max |ai,j | ∥x∥∞ .
i∈J1,nK i∈J1,nK
j=1

Par conséquent
n
X
∥A∥∞ ≤ max |ai,j |.
i∈J1,nK
j=1

Maintenant, il existe p ∈ J1, nK tel que

n
X n
X
|ap,j | = max |ai,j |.
i∈J1,nK
j=1 j=1

ap,j
Soit y le vecteur de Rn défini par yj = si ap,j ̸= 0 et yj = 1 si ap,j = 0. Par construction, ∥y∥∞ = 1.
|ap,j |
Nous avons alors  
Xn
Ay =  ai,j yj 
j=1
i∈J1,nK

et en particulier
n n n n n
X X X a2p,j X X
(Ay)p = ap,j yj = ap,j yj = = |ap,j | = |ap,j |.
j=1 j=1 j=1
|ap,j | j=1 j=1
ap,j ̸=0 ap,j ̸=0 ap,j ̸=0

Comme nous avons aussi

|(Ay)p | ≤ ∥Ay∥∞ ≤ ∥A∥∞ ∥y∥∞ ,
nous obtenons que
n
X
max |ai,j | ≤ ∥A∥∞
i∈J1,nK
j=1

et le résultat est démontré avec les deux inégalités.

2. Des calculs simples donnent

∥A∥∞ = ∥A∥1 = 33 et ∥A1 ∥∞ = ∥A−1 ∥1 = 136.

(remarque si A est symétrique alors ∥A∥∞ = ∥A∥1 !)

Exercice 3. conditionnement d’une matrice

On reprend la matrice A de l’exercice précédent.

3
ENSEM 1A-ISN Analyse Numérique

1. Résoudre les systèmes linéaires

Ax = b, et A(x + δx) = b + δb
avec   
32 32, 1
 23   22, 9 
b= 
 33  et b + δb =  
 33, 1 
31 30, 9
∥δx∥∞ ∥δb∥∞
2. Comparer les erreurs relatives ∥x∥∞ et ∥b∥∞ .

3. Que peut-on dire du conditionnement de la matrice A ? À l’aide de l’exercice précédent donner le condi-
tionnement en norme infinie de la matrice A.

Correction 3.

1. Après calculs on obtient que

  
1 9, 2
 1   −12, 6 
x= 
 1  et x + δx = 
 4, 5 


1 −1, 1

On remarque qu’une petite perturbation du second membre donne une modification assez conséquente de
la solution obtenue. D’un point de vue des erreurs relatives, nous avons

∥δx∥∞ 13, 6 ∥δb∥∞ 0, 1

= = 13, 6 et = .
∥x∥∞ 1 ∥b∥∞ 33

2. Le conditionnement de la matrice A doit donc être mauvais. Effectivement, nous avons

cond∞ (A) = ∥A∥∞ ∥A−1 ∥∞ = 136 × 33 = 4488!

Remarque à faire :
Nous avons vu en cours le résultat suivant :
∥δx∥ ∥δb∥
≤ cond(A) .
∥x∥ ∥b∥

Cet exemple montre que cette estimation est optimale puisque nous avons ici

∥δx∥∞ 0, 1 ∥δb∥∞
= 13, 6 = 136 × 33 × = cond∞ (A) .
∥x∥∞ 33 ∥b∥∞

Exercice 4. Matrice strictement diagonalement dominante

Soit A = (ai,j ) une matrice carrée d’ordre n.
1. Montrer que  
n 
[ X 
sp(A) ⊂ z ∈ C; |z − ai,i | ≤ |ai,j | .
 
i=1 j̸=i

4
ENSEM 1A-ISN Analyse Numérique

2. Montrer que si la matrice A est strictement diagonalement dominante, c’est-à-dire si

X
|aii | > |aij |, ∀ i ∈ J1, nK,
j̸=i

alors elle est inversible.

3. Soit A une matrice strictement diagonalement dominante et soit δ > 0 défini par
X
δ = min{|aii | − |aij |}.
i
j̸=i

Montrer que ∥A−1 ∥∞ ≤ δ −1 .

Correction 4.
1. Soit λ Sp(A). Il existe alors une vecteur non nul x tel que Ax = λx. Nous avons alors pour tout i ∈ J1, nK
n
X
ai,j xj = λxi ,
j=1

soit encore
n
X
(λ − ai,i )xi = ai,j xi .
j=1
j̸=i

Maintenant, il existe p ∈ J1, nK tel que |xp | = maxi∈J1,nK |xi | = ∥x∥∞ . Remarquons que puisque x ̸= 0,
alors |xp | =
̸ 0 et que pour cet indice p nous avons alors
 
n
X 
|λ − ap,p ||xp | ≤ 
 |ai,j |
 |xp |,
j=1
j̸=i

d’où
n
X
|λ − ap,p | ≤ |ai,j |.
j=1
j̸=i

Nous venons de montrer que

   
 X  [n  X 
λ ∈ z ∈ C; |z − ap,p | ≤ |ap,j | ⊂ z ∈ C; |z − ai,i | ≤ |ai,j | .
   
j̸=i i=1 j̸=i

Le résultat est donc démontré.

2. Si la matrice A est strictement diagonalement dominante, il est alors facile de montrer à l’aide de la question
précédente que 0 ∈
/ Sp(A). Par conséquent la matrice est inversible.
3. Soit un vecteur non nul y et soit p ∈ J1, nK tel que |yp | = ∥y∥∞ . Nous avons alors
    
      
 n n n
X      
 X   X 
∥Ay∥∞ ≥ |(Ay)p | =  ap,j yj  = 
  ap,j yj + ap,p yp  ≥ |ap,p ||yp | − 
   ap,j yj 
 j=1   j=1    j=1 
j̸=p   j̸=p 

5
ENSEM 1A-ISN Analyse Numérique

Or,    
 
 n n
X 
 X 
≤ |ap,j |
ap,j yj  |yp | < |ap,p ||yp |.

  
 j=1  j=1
j̸=p  j̸=p

Par conséquent,    
 
 n n
X  X
  
∥Ay∥∞ ≥ |ap,p ||yp | −  ap,j yj  ≥ 
 |ap,p | −  |yp |.
|ap,j |
 j=1  j=1
j̸=p  j̸=p

En conclusion, nous avons montré que

∀ y ̸= 0, ∥Ay∥∞ ≥ δ∥y∥∞

soit encore
∥y∥∞ 1
∀ y ̸= 0, ≤ .
∥Ay∥∞ δ
Ceci peut maintenant se réécrire en posant x = Ay soit y = A−1 x

∥A−1 x∥∞ 1
∀ x ̸= 0, ≤
∥x∥∞ δ

d’où ∥A−1 ∥∞ ≤ δ −1 .
Remarque : pour une matrice strictement diagonalement dominante, il est donc très facile d’obtenir une
majoration du conditionnement en norme infinie (sans calculer A−1 ).

Méthodes itératives

Exercice 5.
Rappel, on note J la matrice de Jacobi et L1 la matrice de Gauss-Seidel. Pour chacune des deux matrices suivantes
calculer les rayons spectraux ρ(J) et ρ(L1 ) et en déduire s’il y a convergence de la méthode de Jacobi et de la
méthode de Gauss-Seidel :    
1 2 −2 2 −1 1
A =  1 1 1 , B =  2 2 2 .
2 2 1 −1 −1 2

Correction 5.

1. Nous avons A = D − E − F avec

     
1 0 0 0 0 0 0 −2 2
D= 0 1 0 , E =  −1 0 0 , F = 0 0 −1  .
0 0 1 −2 −2 0 0 0 0

Par conséquent  
0 −2 2
J = D−1 (E + F ) =  −1 0 −1 
−2 −2 0

6
ENSEM 1A-ISN Analyse Numérique

    
1 0 0 0 −2 2 0 −2 2
L1 = (D − E)−1 F =  −1 1 0  0 0 −1  =  0 2 −3  .
0 −2 1 0 0 0 0 0 2
Nous pouvons maintenant calculer les valeurs propres : Pour Jacobi
 
 −λ −2 2 
−λ −1  = −λ3 .

 −1

 −2 −2 −λ 

Par conséquent ρ(J) = 0 et la méthode de Jacobi converge !

Pour Gauss-Seidel  
 −λ −2 2
−3  = −λ(2 − λ)2

 0 2−λ

 0 0 2−λ
Par conséquent ρ(L1 ) = 2 et la méthode de Gauss-Seidel diverge.
Remarque :
Il n’est pas nécessaire de calculer J et L1 pour déterminer leur rayon spectral.
En effet λ ∈ Sp(J) si et seulement si il existe x ̸= 0 tel que

Jx = λx ⇔ D−1 (E + F )x = λx ⇔ (E + F )x = λDx ⇔ (λD − E − F )x = 0 ⇔ det(λD − E − F ) = 0

et la matrice λD − E − F est très facile à construire.

De même, λ ∈ Sp(L1 ) si et seulement si il existe x ̸= 0 tel que

L1 x = λx ⇔ (D − E)−1 F x = λx ⇔ F x = λ(D − E)x ⇔ (λD − λE − F )x = 0 ⇔ det(λD − λE − F ) = 0

et encore une fois la matrice λD − λE − F est très facile à construire !

2. En utilisant la remarque λ ∈ Sp(J) si et seulement si det(λD − E − F ) = 0, nous avons alors
 
 2λ −1 1 
2λ 2  = 2λ(4λ2 + 5).

det(λD − E − F ) =  2
 −1 −1 2λ 
√
Nous en déduisons que ρ(J) = 25 et la méthode de Jacobi diverge.
Comme λ ∈ Sp(L1 ) si et seulement si det(λD − λE − F ) = 0, nous avons
 
 2λ −1 1 
det(λD − λE − F ) =  2λ 2λ 2  = 2λ(2λ + 1)2 .
 
 −λ −λ 2λ 

Nous en déduisons que ρ(L1 ) = 21 et la méthode de Gauss-Seidel converge.

Pour ceux qui souhaitent déterminer J et L1
1
− 12
 
0 2
−1
J = D (E + F ) =  0 − 21 − 12  .
0 0 − 12

Exercice 6.
On considère le système linéaire Ax = b où A est une matrice triangulaire inférieure inversible.
1. Montrer que le rayon spectral de la matrice de Jacobi J est nul. La méthode de Jacobi est-elle convergente ?

7
ENSEM 1A-ISN Analyse Numérique

2. Montrer qu’il existe p ∈ N tel que J p = 0.

3. Soit (xk )k∈ N une suite obtenue avec la méthode de Jacobi. Montrer que le vecteur x, unique solution du
système, vérifie :
xk+1 − x = J(xk − x), ∀k ≥ 0.
4. Montrer que dans ce cas la méthode de Jacobi est directe, c’est-à-dire que la solution est atteinte en un
nombre fini d’itérations.
5. Calculer tous les termes de la suite (xk ) pour le système :
    
1 0 0 x 1
 −2 1 0   y  =  1 
0 −2 1 z 1

avec x0 = [0, 0, 0]T .

Correction 6.

1. Tout d’abord, remarquons que si la matrice A est une matrice triangulaire inférieure, alors A = D−E −F =
D − E puisque ici F = 0 et la matrice J = D−1 (E + F ) = D−1 E est aussi triangulaire inférieure.
Nous savons aussi que A est inversible, donc det(A) ̸= 0 et comme A est triangulaire inférieure det(A) =
Qn
ai,i ̸= 0.
i=1
Maintenant, nous avons vu dans les exercices précédents que
λ ∈ Sp(J) ⇔ det(λD − E − F ) = det(λD − E) = 0.
Or
n
Y n
Y
det(λD − E) = λai,i = λn ai,i .
i=1 i=1
Donc
λ ∈ Sp(J) ⇔ λ = 0.
Le rayon spectral de la matrice J est donc nul et la méthode de Jacobi est donc toujours convergente pour
une matrice triangulaire inférieure inversible.
2. Il s’agit simplement d’une application du théorème de Cayley Hamilton. Puisque J est triangulaire inférieure,
son polynôme caractéristique est de la forme
n
Y
PJ (X) = det(J − XId ) = (Ji,i − X).
i=1

Or nous avons vu que 0 est l’unique valeur propre de J, donc pour tout i ∈ J0, nK, Ji,i = 0 et PJ (X) = X n .
Le théorème de Cayley Hamilton assure alors que
PJ (J) = J n = 0.

3. Nous avons pour tout k ∈ N

xk+1 = Jxk + D−1 b
et aussi
x = Jx + D−1 b.
Par soustraction, nous obtenons le résultat annoncé.

8
ENSEM 1A-ISN Analyse Numérique

4. Il est facile de montrer à partir de la question précédente et à l’aide d’une démonstration par récurrence
que pour
∀ k ≥ 0, (xk − x) = J k (x0 − x).
Maintenant, comme nous avons vu que J n = 0, nous avons donc

(xn − x) = 0.

La solution est donc bien obtenue en un nombre fini d’itérations.

5. Pour ce système nous avons
   
0 0 0 1
J = 2 0 0 , D−1 b = b =  1  .
0 2 0 1

Nous obtenons alors      

1 1 1
x1 =  1  , x2 =  3  x3 =  3 
1 3 7
Et comme attendu, x3 est bien la solution du système linéaire.

Exercice 7.
Soit A = (aij ) une matrice d’ordre n strictement diagonalement dominante :
X
|aii | > |aij |, ∀ i ∈ J1, nK.
j̸=i

Démontrer que les méthodes de Gauss-Seidel et de Jacobi convergent.

Correction 7.

1. Méthode de Jacobi :
Soit λ ∈ Sp(J). Alors il existe x ̸= 0 tel que (E + F )x = λDx.
Soit p ∈ J1, nK tel que |xp | = max |xi |. Remarquons que xp ̸= 0 puisque x ̸= 0.
i∈J1,nK
En considérant la p-ième composante du vecteur (E + F )x = λDx, nous avons
n
X
− ap,j xj = λap,p xp ,
j = 1
j ̸= p

d’où
n
X n
X
|λ||ap,p ||xp | ≤ |ap,j ||xj | ≤ |xp | |ap,j |.
j = 1 j = 1
j ̸= p j ̸= p

Par conséquent Pn
j = 1 |ap,j |
j ̸= p
|λ| ≤ <1
|ap,p |
puisque A est strictement diagonalement dominante. Conclusion ρ(J) < 1 et la méthode de Jacobi converge.

9
ENSEM 1A-ISN Analyse Numérique

2. Méthode de Gauss-Seidel :
Soit λ ∈ Sp(L1 ), λ ̸= 0. Alors il existe x ̸= 0 tel que (λE + F )x = λDx.
Soit p ∈ J1, nK tel que |xp | = max |xi |. Remarquons que xp ̸= 0 puisque x ̸= 0.
i∈J1,nK
En considérant la p-ième composante du vecteur (λE + F )x = λDx, nous avons
X X
− λap,j xj − ap,j xj = λap,p xp ,
1≤j<p p<j≤n

d’où  
X X X X
|λ||ap,p ||xp | ≤ |λ||ap,j ||xj | + |ap,j ||xj | ≤  |λ||ap,j | + |ap,j | |xp |
1≤j<p p<j≤n 1≤j<p p<j≤n

Par conséquent X X
|λ||ap,p | ≤ |λ| |ap,j | + |ap,j |.
1≤j<p p<j≤n
X
Or, puisque A est strictement diagonalement dominante, |ap,j | < |ap,p |. Il en résulte (rappel λ ̸= 0) que
j̸=p

X X X
|λ| |ap,j | < |λ| |ap,j | + |ap,j |
j̸=p 1≤j<p p<j≤n

et donc que X
(1 − |λ|) |ap,j | > 0.
p<j≤n

Nous avons donc |λ| < 1. Conclusion ρ(L1 ) < 1 et la méthode de Gauss-Seidel converge.

Exercice 8.
Considérons le système linéaire :
1 1 1


 x1 − 4 x3 − 4 x4 = 2
1 1 1
x2 − 4 x3 − 4 x4 =

2
− 1 x1 − 1
4 x2 + x3 = 1
 41 2

 1 1
− 4 x1 − 4 x2 + x4 = 2

1. Sans calculer les rayons spectraux, montrer que les méthodes de Jacobi et Gauss-Seidel convergent.
2. En choisissant x0 = (0, 0, 0, 0)T , calculer les quatre premières itérations des méthodes de Jacobi et de
Gauss-Seidel. Faire la comparaison.

Correction 8.

1. Ici la matrice est strictement diagonalement dominante. Elle est donc inversible et les méthodes de Jacobi
et de Gauss Seidel convergent d’après l’exercice précédent !
2. Je vous laisse faire...

10
ENSEM 1A-ISN Analyse Numérique

Factorisation LU

Exercice 9. Factorisation LU
1. Déterminer la factorisation LU de la matrice
 
1 2 3 4
 2 5 1 10 
A=
 3

1 35 2 
4 10 2 45

2. Résoudre par un algorithme de descente-remontée le système linéaire

 
−2
 −12 
Ax = b où b =   35 


−49

Correction 9.

1. Étape 1 :    
1 2 3 4 1 0 0 0
 0 1 −5 2   2 1 0 0 
A1 =   L= 
 0 −5 26 −10   3 . 1 0 
0 2 −10 29 4 . . 1
Étape 2 :    
1 2 3 4 1 0 0 0
 0 1 −5 2   2 1 0 0 
A2 =   L= 
 0 0 1 0   3 −5 1 0 
0 0 0 25 4 2 . 1
Étape 3 : la troisième colonne étant déjà réduite, il n’y a rien à faire. On obtient directement :
   
1 2 3 4 1 0 0 0
 0 1 −5 2   2 1 0 0 
U = A3 =   0
 L=
 3 −5

0 1 0  1 0 
0 0 0 25 4 2 0 1

2. La résolution de Ly = b (descente) donne  

−2
 −8 
y= 
 1 
−25
et la résolution de U x = y (remontée) donne
 
1
 −1 
x=
 1 .


−1

11
ENSEM 1A-ISN Analyse Numérique

Exercice 10. Algorithme de la factorisation LU

1. Écrire l’algorithme de la factorisation LU d’une matrice A d’ordre n.
2. Déterminer le nombre d’opérations élémentaires nécessaires pour effectuer la factorisation LU d’une matrice
A d’ordre n.
3. Écrire un algorithme de descente-remontée pour résoudre le système linéaire Ax = b lorsque la factorisation
LU de la matrice A est connue.
4. Déterminer le nombre d’opérations élémentaires nécessaires pour effectuer une descente-remontée.

Correction 10.
1. Voici un premier algorithme qui stocke les matrices L et U dans deux tableaux distincts :

Algorithm 1 Factorisation LU
L ← Id
U ←A
for j ← 1, n − 1 do
if U (j, j) ̸= 0 then
for i ← j + 1, n do
L(i, j) ← U (i, j)/U (j, j)
U (i, j) ← 0
for k ← j + 1, n do
U (i, k) ← U (i, k) − L(i, j) ∗ U (j, k)
end for
end for
else
Stop
end if
end for

Il est aussi possible de stocker les deux matrices L et U dans un unique tableau (noté ci-après LU ) et
d’économiser ainsi un peu de mémoire.
Pour cela, il faut modifier un peu l’algorithme précédent qui devient alors :

Algorithm 2 Factorisation LU
LU ← A
for j ← 1, n − 1 do
if LU (j, j) ̸= 0 then
for i ← j + 1, n do
LU (i, j) ← LU (i, j)/LU (j, j)
for k ← j + 1, n do
LU (i, k) ← LU (i, k) − LU (i, j) ∗ LU (j, k)
end for
end for
else
Stop
end if
end for

12
ENSEM 1A-ISN Analyse Numérique

2. Le nombre d’opérations élémentaires est : n(n−1)(4n+1)

6 ∼ 23 n3 .
n
X n(n + 1)(2n + 1)
On peut rappeler la formule j2 = .
j=1
6
3. Voici l’algorithme de descente remontée :

Algorithm 3 Descente-Remontée
y←b ▷ Descente
for j ← 1, n − 1 do
for i ← j + 1, n do
y(i) ← y(i) − L(i, j) ∗ y(j)
end for
end for
x←y ▷ Remontée
for j ← n, 1 do
x(j) ← x(j)/U (j, j)
for i ← 1, j − 1 do
x(i) = x(i) − U (i, j) ∗ x(j)
end for
end for

4. Le nombre d’opérations élémentaires est : n(n − 1) pour la descente et n2 pour la remontée soit de l’ordre
de 2n2 .

Exercice 11. Unicité de la factorisation A = LU

Soit A une matrice inversible qui admet une factorisation A = LU , c’est à dire L est triangulaire inférieure avec
des éléments diagonaux tous égaux à 1 et U est une matrice triangulaire supérieure. Démontrer l’unicité de la
factorisation.

Correction 11.
Supposons que la matrice A admette deux factorisations LU . Soient (L, U ) et (L̃, Ũ ) tels que

A = LU = L̃Ũ

où L et L̃ sont triangulaires inférieures avec des éléments diagonaux tous égaux à 1 et U et Ũ sont triangulaires
supérieures.
Nous avons alors
L̃−1 L = Ũ U −1 .
La matrice L̃−1 est aussi triangulaire inférieure avec des éléments diagonaux tous égaux à 1, donc la matrice L̃−1 L
est triangulaire inférieure avec des éléments diagonaux tous égaux à 1.
De même, la matrice U −1 est aussi triangulaire supérieure, donc la matrice Ũ U −1 est triangulaire supérieure.
Par conséquent L̃−1 L = Ũ U −1 est une matrice à la fois triangulaire supérieure et triangulaire inférieure, elle est
donc diagonale.
De plus, comme tous les coefficients diagonaux de L̃−1 L sont égaux à 1. Nous avons donc

L̃−1 L = Ũ U −1 = I.

On en déduit alors que

L̃ = L et Ũ = U.

13
ENSEM 1A-ISN Analyse Numérique

Factorisation de Choleski

Exercice 12. Factorisation de Choleski

Dans cet exercice, on reprend la matrice A de l’exercice 9 qui est une matrice symétrique définie positive.
1. Déduire directement de la factorisation LU obtenue précédemment, la factorisation de Choleski A = BB T .
2. Effectuer maintenant cette même factorisation de Choleski A = BB T à l’aide de l’algorithme ”classique”.

Correction 12.
1. Pour toute matrice diagonale inversible D, nous pouvons écrire A = LU = (LD)(D−1 U ) où alors LD est
triangulaire inférieure et D−1 U est triangulaire supérieure (rappel le produit de deux matrices triangulaires
inférieures est une matrice triangulaire inférieure, de même pour les matrices triangulaires supérieures).
On peut aussi remarquer que les matrices L et U obtenues précédemment sont ”presque” symétriques. On
choisit ici (rappel dans la factorisation de Choleski les coefficients diagonaux sont strictement positifs) :
 
1 0 0 0
 0 1 0 0 
D=  0 0 1 0 


0 0 0 5

Nous avons alors

A = (LD)(D−1 U )
    1 0 0 0
 
1 0 0 0 1 0 0 0 1 2 3 4
 2 1 0 0  0 1 0 0   0
  1 0 0  0 1 −5 2 
  
=    
0   0 0 1 0
 
 3 −5 1 0  0 0 1

 0

0 1 0 

1
4 2 0 1 0 0 0 5 0 0 0 0 0 0 25
   5
1 0 0 0 1 2 3 4
 2 1 0 0  0 1 −5 2 
 = BB T
= 
 3 −5 1 0

 0 0 1 0 
4 2 0 5 0 0 0 5

2. Voici l’algorithme de Choleski :

Algorithm 4 Factorisation Cholesky

for j ← 1,v
n do
u
u Xj−1
b = ta −
j,j j,j b2 j,k
k=1
for i ← j + 1, n do
j−1
!
X
bi,j = ai,j − bi,k bj,k /bj,j
k=1
end for
end for

14
ENSEM 1A-ISN Analyse Numérique

En suivant cet algorithme, on a construit la matrice B colonne par colonne :

Étape 1 : j = 1
√
b1,1 = a1,1 = 1,
b2,1 = 2/1 = 2, b3,1 = 3/1 = 3, b4,1 = 4/1 = 4.
Étape 2 : j = 2 q √
b2,2 = a2,2 − b22,1 = 5 − 4 = 1,

b3,2 = (a3,2 − b2,1 b3,1 )/b2,2 = (1 − 2 × 3) = −5, b4,2 = (a4,2 − b2,1 b4,1 )/b2,2 = (10 − 2 × 4) = 2.
Étape 3 : j = 3 q √
b3,3 = a3,3 − b23,1 − b23,2 = 35 − 9 − 25 = 1,

b4,3 = (a4,3 − b3,1 b4,1 − b3,2 b4,2 )/b3,3 = (2 − 3 × 4 + 5 × 2) = 0.

Étape 4 : j = 4 q √
b4,4 = a4,4 − b24,1 − b24,2 − b24,3 = 45 − 16 − 4 − 0 = 5.

Conclusion :  
1 0 0 0
 2 1 0 0 
B= 
 3 −5 1 0 
4 2 0 5

Exercice 13.
Soit A une matrice carrée d’ordre n, symétrique définie positive. On cherche à résoudre le système A2 x = b.
1. On propose la la méthode suivante pour la résolution de ce système :
• On calcule la matrice C = A2 .
• On effectue ensuite la factorisation de Choleski C = BB T .
• On résout enfin le système linéaire par une méthode de descente remontée.
Calculer le nombre d’opérations élémentaires nécessaires pour réaliser cette méthode.
2. Proposer maintenant une autre méthode nécessitant moins d’opérations élémentaires.

Correction 13.
1. La méthode proposée est en trois parties. Calculons le nombre d’opérations pour chaque partie :
2
■ Calcul de la matrice C = A :
Rappelons qu’en notant A = (ai,j ) et C = (ci,j ), nous avons
n
X
∀ (i, j) ∈ (J1, nK)2 , ci,j = ai,k × ak,j .
k=1

Le calcul de chaque coefficients nécessite donc n multiplications et (n − 1) additions, soit pour les
n2 coefficients de la matrice C : n3 multiplications et n2 (n − 1) additions, soit un total de 2n3 − n2
opérations élémentaires.
■ Factorisation de Choleski C = BB T :
3 3 2
Nous avons vu dans le cours qu’elle nécessite n6 additions, n6 multiplications, n2 divisions et n extrac-
3 2
tions de racines, soit un total de n3 + n2 opérations élémentaires et n extractions de racines.

15
ENSEM 1A-ISN Analyse Numérique

■ Résolution du système linéaire par une méthode de descente remontée :

Pour la descente : n(n−1)
2 soustractions, n(n−1)
2 multiplications et n divisions ;
Pour la remontée : n(n−1)
2 soustractions, n(n−1)
2 multiplications et n divisions.
Soit au total 2n2 opérations élémentaires.
2. Il est préférable de factoriser directement A :

A = BB T

On a alors
BB T BB T x = b
On fait alors deux descentes et remontées :

By = b, B T z = y, Bw = z, B T x = w.

On fait une descente remontée de plus, mais on ne calcule plus la matrice C. Or le calcul de celle-ci est
beaucoup plus couteux qu’une descente remontée !
En effet, une descente remontée nécessite 2n2 opérations élémentaires tandis que le calcul de C nécessite
2n3 − n2 opérations élémentaires. Le gain est donc important.

16