Math pour Ingénieurs: EEA-GOI

Méthodes itératives
• Introduction
• Normes vectorielles et normes matricielles
• Méthode de Jacobi (1804-1851)
• Méthode de Gauss-Seidel
• Méthode de relaxation
• Convergence et test darrêt
• Conditionnement d’une matrice

1
Introduction
La méthode d’élimination de Gauss ou la factorisation LU sont
efficaces pour la résolution des systèmes d’ordre pas trop grand.
Dans les théories des équations aux dérivées partielles on tombe
facilement sur des systèmes d’ordre 104 ou 105 . Dans ce cas,
l’application de la méthode d’élimination de Gauss ou la factorisation
LU devient très coûteuse en terme de temps de calcul.
Comme pour la recherche de zéro de fonction non linéaie f (x) = 0,
on a ramené le problème à la recherche de point fixe x = g(x) par les
méthodes itératives en regardant la convergence de la suite
xn+1 = g(xn ) pour un x0 initialement choisi.
De manière similaire, le système d’ordre n: Ax = b peut être écrit
sous la forme: x = Bx + c où B est une matrice n × n et c un

2
vecteur. On genére la suite:

x(0) donne
x(k+1) = Bx(k) + c , pour k = 0, 1, 2, . . . (1)

et étudié sa convergence.
Avant d’exposer ses méthodes itératives, nous présentons d’abord
quelques notions sur les normes vectorielles et normes matricielles.

3
Normes vectorielles et normes matricielles

Introduction
Lorsque nous avons étudié les solutions de f (x) = 0, nous avons
utilisés |xn − s| pour comparer la solution approximative xn calculer
par itération et la valeur exacte s. Nous avons également calculé
|xn+1 − xn | comme test d’arrêt des itérations.
Lorsqu’on s’intéresse à la résoluton des systèmes du type Ax = b par
méthode itérative, on doit se poser la question suivante:
Comment peut on définir la limite limk→+∞ x(k) = x?
Dans le cas à une variable, on a pu définir l’erreur relative comme
étant |xn − xn−1 |/|xn |, dans le cas des vecteurs la division n’est pas
définie.
On a donc besoin d’un moyen pour mesurer les écarts entre deux
vecteurs. Nous allons donc étendre la notion de valeur absolue dans
4
le cas à une dimension au cas vectoriel.

5
Norme vectorielle

A présent nous allons définir une norme vectorielle sur

Rn = {x|x = (x1 , x2 , . . . , xn )t ; x1 , x2 , . . . , xn ∈ R}

Par norme vectorielle sur Rn nous sous entendons une fonction réelle
kk de Rn vers R+ qui satisfait les conditions suivantes:
(i) kxk ≥ 0 pour tout x ∈ Rn et kxk = 0 ssi x = 0
(ii) kλxk = λkxk pour tou λ ∈ R et x ∈ Rn
(iii) kx + yk ≤ kxk + kyk, pour tout x, y ∈ Rn
Dans le cas de Rn , les normes les plus utilisées en pratique sont: pour

6
tout x = (x1 , x2 , . . . , xn )t ∈ Rn
n
X
kxk1 = |x1 | + |x2 | + . . . + |xn | = |xi |
i=1
q
kxk2 = x21 + x22 + . . . + x2n
kxk∞ = max{|x1 |, |x2 |, . . . , |xn |}

kxk2 est ce qu’on appelle la norme Euclidienne. En fait, toutes ces

normes provienent de la norme lp définie par
n
X
kxkp = (|x1 |p + |x2 |p + . . . + |xn |p )1/p = ( |xi |p )1/p
i=1

La norme kk1 correspond à p = 1, alors que kk2 correspond à p = 2.

On montre aussi que kxk∞ = limp→∞ kxkp .
La norme d’un vecteur donne la mesure de la distance entre ce

7
vecteur et les vecteurs d’origine 0Rn = (0, 0, . . . , 0)t , la distance entre
deux vecteurs peut être définie comme la norme de la différence des
deux vecteurs.

8
Définition: Si x = (x1 , x2 , . . . , xn )t et y = (y1 , y2 , . . . , yn )t sont deux
vecteurs de Rn , la distance entre x et y est d(x, y) = kx − yk avec
x − y = (x1 − y1 , x2 − y2 , . . . , xn − yn )t . Dans le cas de l2 et l∞ , on a:
n
X
kx − yk2 = { (xi − yi )2 }1/2 et kx − yk∞ = max1≤i≤n |xi − yi |
i=1

9
Exemples

1.) Soit x = (2, −3, 5)t un vecteur de R3 .

p √
2 2
kxk1 = |2| + | − 3| + |5| = 10, kxk2 = 2 + (−3) + (5) = 38, 2

kxk∞ = max{|2|, | − 3|, |5|} = 5.

2.) Si x = (2, −3, 5)t , y = (1, 3, 0)t on a x − y = (1, −6, 5)t , alors
√
kx − yk1 = 12, kx − yk2 = 62 et kx − yk∞ = 6.

10
Limite d’une suite vectorielle

Dans R2 , on définit la suite de vecteurs x(k) = (1/2k , 1 + 10−k )t . On

a x(1) = (1/2, 1.1)t , x(2) = (1/22 , 1.01)t , x(3) = (1/23 , 1.001)t ,
x(4) = (1/24 , 1.0001)t , . . .. On remarque que la suite converge vers
(k)
(0, 1)t car limk→∞ x1 = limk→∞ 1/2k = 0 et
(k)
limk→∞ x2 = limk→∞ (1 + 10−k ) = 1.
En général, dans Rn , on dit que la suite x(k) converge vers x
(k)
(limk→∞ x(k) = x) si et seulement si limk→∞ xi = xi pour chaque
i = 1, 2, . . . , n.
(k)
Comme pour chaque i = 1, ..., n, limk→∞ xi = xi alors
(k) (k)
limk→∞ |xi − xi | = 0. Ce qui donne limk→∞ max|xi − xi | = 0 et
par suite limk→∞ kx(k) − xk∞ = 0. Le résultat reste vraie pour
n’importe quelle norme, du fait que les normes sont équivalentes.

11
Théorème: limk→∞ x(k) = x ssi limk→∞ kx(k) − xk = 0, où kk
représente une norme quelconque.

12
Norme matricielle

On suppose que Mn est l’ensemble des matrices n × n. Mn est un

espace vectoriel de dimension n2 et donc on peut introduire une
n2
norme comme dans un espace de R . Cependant Mn n’est pas
uniquement un espace vectoriel de dimension supérieur, il a aussi une
opération de multiplication. Il est utile de relier la mesure du produit
AB à la mesure de A et de B. Par conséquent, pour la norme
matricielle on imposera la condition:

kABk ≤ kAkkBk

Avec cette condition suplémentaire, toutes les normes vectorielles

definies ci-dessus ne sont pas des normes matricielles.

13
Définition: Une norme matricielle est une application de Mn vers
R+ qui vérifie les propriétés suivantes:
(1) kAk ≥ 0 pour tout A ∈ Mn and kAk = 0 ssi A = 0.
(2) kαAk = αkAk pour tout α ∈ R et A ∈ Mn
(3) kA + Bk ≤ kAk + kBk pour tout A, B ∈ Mn
(4) kABk ≤ kAkkBk pour tout A, B ∈ Mn .
Il y a plusieurs façons de définir des normes matricielles qui satisfont
les conditions précédentes. Cependant, les seules normes que nous
considérons sont ceux qui sont une conséquence naturelle des normes
vectorielles.

14
Théorème: Si k.k est une norme vectorielle de Rn , alors

kAk = maxkxk=1 kAxk (2)

est une norme matricielle, elle est appelée norme matricielle naturelle.

15
Corollaire: Pour tout x ∈ Rn , x 6= 0, A ∈ Mn une matrice d’ordre n
et k.k une norme matricielle naturelle de Mn , on a:

kAxk ≤ kAkkxk

Preuve:
x
Pour tout x 6= 0, la norme k kxk k = 1 donc:

x
kA( )k ≤ kAk
kxk
x 1
or A( kxk )= kxk Ax et par suite:

1 1 x
kAxk = k Axk = kA( ) ≤ kAk
kxk kxk kxk
ce qui implique que:
kAxk ≤ kAkkxk

16
Nous donnons dans le tableau suivant, les principaux normes
vectorielles et les normes matricielles associées
Norme vectorielle norme matricielle associée
Pn Pn
kxk1 = i=1 |xi | kAk1 = maxj i=1 |aij |
le max porte sur la somme des colonnes
Pn p
kxk2 = ( i=1 x2i )1/2 kAk2 = ρ(At A)
ρ est le rayon spectral de At A
Pn
kxk∞ = maxi |xi | kAk∞ = maxi j=1 |aij |
le max porte sur la somme des lignes

17
0.0.1 Exemples

Calculer kAk1 , kAk2 et kAk∞ si

 
3 1 −1
 
A=  1 5 1 

1 −1 8

kAk1 = max{|3| + |1| + |1|, |1| + |5| + | − 1|, | − 1| + |1| + |8|} = 10

p
kAk2 = ρ(At A), les valeurs propres de At A sont: 66.9142, 29.5613,
p √
t
7.52451. Donc kAk2 = ρ(A A) = 66.9142 = 8.18011.
kAk∞ = max{|3| + |1| + | − 1|, |1| + |5| + |1|, |1| + | − 1| + |8|} = 10

18
1 Méthode de Jacobi (1804-1851)
On expose d’abords cette méthode itérative sur un exemple, le cas
général sera traité juste après.
1.1 Exemple
Soit à résoudre le système suivant:

 4x + x = 5
1 2
(3)
 x1 − 5x2 = −4

La solution exacte est x1 = x2 = 1. Ce système est équivalent à:

      
x1 0 − 41 x1 5
 =  + 4  (4)
1 4
x2 5 0 x 2 5
| {z }
c

19
 (k) (k)
A ce niveau, on introduit la suite x(k) = (x1 , x2 )t ou x = (x1 , x2 )t
avec une valeur d’initialisation x(0) = (0, 0):

x(k+1) = Bx(k+1) + c (5)

pour x(0) = (0, 0), on a x(1) = (1.25, 0.8), x(2) = (1.05, 1.5) ... les
valeurs de la suite x(k) pour k = 0, 1, ..., 10 sont données dans le
tableau suivant:

20
 (k) (k) kx(k+1) −x(k+1) k
k x1 x2 kx(k+1) k
kx(k) − (1, 1)t k
0 0 0 1 1.00069
1 1.25 0.8 0.23793 0.25
2 1.05 1.05 0.0618812 0.0500347
3 0.9875 1.01 0.0125226 0.0125
4 0.9975 0.9975 0.00312217 0.00250173
5 1.00063 0.9995 0.000624489 0.000625
6 1.00013 1.00013 0.000156142 0.000125087
7 0.999969 1.00003 0.0000312285 0.00003125
8 0.999994 0.999994 7.80709 10−6 6.2543310−6
9 1.000002 0.999999 1.56142 10−6 1.562510−6

21
On voit clairement que pour k grand la suite x(k) converge vers la
solution exacte.
Sur ce tableau nous avons illustrés également l’erreur relative
kx(k+1) −x(k+1) k
kx(k+1) k
et l’erreur absolue kx(k) − (1, 1)t k. Il est clair d’après
kx(k+1) −x(k+1) k
ce tableau que kx(k+1) k
et kx(k) − (1, 1)t k sont du même ordre
kx(k+1) −x(k+1) k
de grandeur, on pourra donc utiliser comme test
kx(k+1) k
d’arrêt des itérations.
Par la suite, nous allons introduire le principe général des méthodes
itératives.

22
1.2 Méthode de Jacobi (1804-1851)

Considèrons le système d’ordre n suivant:




 a11 x1 + a12 x2 + . . . + a1n xn = b1



 a21 x1 + a22 x2 + . . . + a2n xn = b2



 .. .. ..

 . . .
(6)


 ai1 x1 + ai2 x2 + . . . + ain xn = bi



 .. .. ..


 . . .


 a x + a x + ...+ a x = b
n1 1 n2 2 nn n n

Si tous les éléments diagonaux aii sont non nuls, alors on peut

23
reécrire le système (6) comme suit:



 x1 = (b1 − a12 x2 − . . . − a1n xn )/a11



 x2 = (b2 − a21 x1 − . . . − a2n xn )/a22



 ... .. ..


. .
(7)


 xi = (bi − ai1 x1 − ai2 x2 − . . . − ain xn )/aii



 .. .. ..


 . . .


 x = (b − a x − a x − . . . − a
n n n1 1 n2 2 nn−1 xn−1 )/ann

ce qu’on peut écrire de manière compacte comme

n
1 X
xi = (bi − aij xj ) , i = 1, 2, . . . , n (8)
aii
j=1,j6=i

si aii 6= 0.

24
La méthode itérative de Jacobi est définie par:
n
(k+1) 1 X (k)
xi = (bi − aij xj ) , i = 1, 2, . . . , n et k = 0, 1, ... (9)
aii
j=1,j6=i

ce qui est facile à programmer sur un ordinateur.

Il est plus utile d’écrire l’équation (9) sous forme matricielle. Pour
cela, la matrice A du système peut se mettre sous la forme:
 
a a12 . . . a1n
 11 
 
 a21 a22 . . . a2n 
A=  .. .. .. ..  = L + D + U
 (10)
 . . . . 
 
an1 an2 . . . ann

avec D matrice diagonale, L matrice strictement supérieure et U

25
matrice strictement inférieure et sont données par:
    
a 0 ... 0 0 0 ... 0 0 a12
 11    
    
 0 a22 . . . 0   a21 0 ... 0   0 0
D= .

. . .
 , L= .
. . .
 , U =
 .. ..
 .. .. .. ..   .. .. .. .. 
 
 . .
    
0 0 . . . ann an1 an2 . . . 0 0 0

Le système Ax = b qui est en fait (D + L + U )x = b devient:

Dx = −(L + U )x + b (12)
ce qui donne:
x = −D −1 (L + U )x + b (13)
La méthode itérative de Jacobi peut être écrite comme:
x(k+1) = −D −1 (L + U )x(k) + D −1 b , pour k = 0, 1, ... (14)

26
1.2.1 Exemple

Soit à résoudre par la méthode de Jacobi le système suivant:


 3x1 + x2 − x3 = 8


x1 + 5x2 + x3 = −2 (15)



x1 − x2 + 8x3 = 4

27
La solution exacte est x1 = 3, x2 = −1 et x3 = 0.
On applique l’algorithme de l’équation avec
     
1
3 0 0 0 0 0 0 1 −1
−1
     
1
D = 0 5 0  , L= 1 0 0   , U =
 0 0 1 
  

0 0 81 1 −1 0 0 0 0
     
1 1 8
8 0 3 − 3 3
     
b= −2  , D −1 (L + U ) =  1 0 1  , D −1
b= 2(16)
−5 
   5 5   
1 1 1
4 8 − 8 0 2

On a alors:
 (k+1)    (k)
  
1
x 0 − 13 x1 8
 1(k+1)   3
 (k)
  3

 x  = − 1 0 1  x2 + − 25  (17)
 2   5 5    
(k+1) 1 (k)
x3 8 − 18 0 x3 1
2

28
On choisit x(0) = (0, 0, 0)t pour initialiser, on obtient pour x(1) et x(2)
les valeurs suivantes: x(1) = ( 38 , − 25 , 21 )t , x(2) = ( 89 31 7 t
30 , − 30 , 60 ) . Les
résultats des 10 premières itérations de Jacobi sont donnés dans le
tableau suivant

29
 (k) (k) (k)
k x1 x2 x3
0 0 0 0
1 2.6666667 -0.4 0.5
2 2.9666667 -1.033333 0.1166667
3 3.0500000 -1.016667 0
4 3.0055556 -1.010000 -0.0083333
5 3.0005556 -0.999444 -0.0019444
6 2.9991667 -0.999722 0.0000000
7 2.9999074 -0.999833 0.0001389
8 2.9999907 -1.000009 0.0000324
9 3.0000139 -1.000005 0.0000000
10 3.0000015 -1.000003 -2.314815 10−6
20 2.9999999 -1.000000 2.5006 10−12
30
On remarque que la suite x(k) converge vers la solution exacte
x1 = 3, x2 = −1 et x3 = 0.

31
2 Méthode de Gauss-Seidel
Dans le but de résoudre l’équation (6) par la méthode de Gauss-Seidel
(1821-1896), on reécrit l’équation (8) de la façon suivante:
i−1 n
1 X X
xi = (bi − aij xj − aij xj ) , i = 1, 2, . . . , n (18)
aii
j=1,j6=i j=i+1,j6=i

Le principe de la méthode Gauss-Seidel est le suivant: Lors de la

première itération,
(1) (0) (1)
pour calculer x1 on aura besoin de x2,...,n qui sont donnés, le x1
est ainsi stocké.
(1) (0) (0)
De même, pour calculer x2 , on aura besoin de x1 et x3,...,n . Dans
(0)
le cas de la méthode de Gauss-Seidel on utilisera au lieu de x1 la
(1) (0)
valeur déja stocké de x1 et pour le reste on utilise x3,...,n .
(1) (1)
et ainsi de suite pour le calcul de xi on utilise x1,2,...,i−1 qui sont

32
 (0)
déja stocké et xi+1,...,n .
On procède de la même façon pour la kème itération, c’est à dire que
(k) (k)
pour le calcul de xi on utilise les valeurs de x1,2,...,i−1 déja évaluées
(k−1)
et xi+1,...,n de l’itération précédente.
La méthode itérative de Gauss-Seidel peut être écrite comme:
i−1 n
(k+1) 1 X (k+1)
X (k)
xi = (bi − aij xj − aij xj ) , i = 1, 2, . . . , n , k = 0, 1
aii
j=1,j6=i j=i+1,j6=i

33
Cette équation peut être écrite de la manière suivante:
(k+1) (k)
a11 x1 = b1 − a12 x2 − . . . − a1n
(k+1) (k+1) (k)
a21 x1 + a22 x2 = b2 − a13 x3 . . . − a2n x(n
.. .. ..
. . .
(k+1) (k+1) (k+1)
an−11 x1 + an−12 x2 + ... + an−1n−1 xn−1 = bn−1 − an−1n x(k)
n
(k+1) (k+1) (k+1)
an1 x1 + an2 x2 + an3 x3 + . . . + ann x(k+1)
n = bn

en notation matricielle, (20) devient:

      
(k+1)
a 0 ... 0 x b1 0 a12 . . . a1n x
 11  1     
   (k+1)     
 a21 a22 . . . 0   x2   b2   0 0 . . . a2n   x
  = − 
 .. .. .. ..  ..   ..   .. .. .. ..  
 . . . .  .   .   . . . .  
      
(k+1)
an1 an2 . . . ann xn bn 0 0 ... 0 x

34
en utilisant les matrices D, L et U , on aura:

(D + L)x(k+1) = b − U x(k) (22)

d’ou l’algorithme de Gauss-Seidel:

x(0) donne (23)

x(k+1) = −(D + L)−1 U x(k) + (D + L)−1 b (24)

2.0.2 Exemple

On reprend l’exemple précédent:


 3x1 + x2 − x3 = 8


x1 + 5x2 + x3 = −2 (25)



x1 − x2 + 8x3 = 4

35
On choisit x(0) = (0, 0, 0)t pour initialiser, on trouve
14 1 t
x(1) = ( 83 , − 15 539
, 20 ) , x(2) = ( 180 227
, − 225 1
, − 2400 )t . Les résultats des 10
premières itérations de Gauss-Seidel sont donnés dans le tableau
suivant

36
 (k) (k) (k)
k x1 x2 x3
0 0 0 0
1 2.666666667 -0.933333333 0.05000000000
2 2.994444444 -1.008888889 -0.0004166666667
3 3.002824074 -1.000481481 -0.0004131944444
4 3.000022762 -0.999921913 6.91550925910−6
5 2.999976276 -0.999996638 3.38565779310−6
6 3.000000008 -1.000000679 -8.58430587710−8
7 3.000000198 -1.000000022 -2.74984561210−8
8 2.999999998 -0.999999994 9.44512624110−10
9 2.999999998 -1.000000000 2.21282862510−10
10 3.000000000 -1.000000000 -9.71496238810−12

37
On remarque qu’il y a convergence vers x = (3, −1, 0), la convergence
dans le cas de Gauss-Seidel est relativement plus rapide que dans le
cas de Jacobi.

38
3 Méthode de relaxation
Les méthodes de Jacobi et Gauss-seidel etudiées précédement
convergent. Il est possible d’améliorer la convergence de la méthode
de Gauss-Seidel. Dans le but de résoudre Ax = b par la méthode de
Gauss-Seidel nous avons introduit la méthode itérative:
i−1 n
(k+1) 1 X (k+1)
X (k)
x̃i = (bi − aij xj − aij xj (26)
aii j=1 j=i+1

pour i = 1, 2, . . . , n.
On pourra généraliser la méthode de Gauss-Seidel en introduisant un
paramètre de relaxation ω pour former une combinaison qui converge
plus vite que Gauss-Seidel pour un choix approprié de ω. Le
paramètre ω est introduit de telle sorte à ce que
(k+1) (k) (k+1) (k)
xi = xj + ω(x̃j − xj ) (27)

39
si ω = 1, la méthode de relaxation (27) se réduit a Gauss-Seidel.
En injectant la relation (26) dans (27) on trouve
i−1 n
(k+1) (k) ω X (k+1)
X (k)
xi = (1 − ω)xj + (bi − aij xj − aij xj (28)
aii j=1 j=i+1

pour i = 1, 2, . . . , n. Cette dernière équation peut être réarrangée

comme:
i−1
X n
X
(k+1) (k+1) (k) (k)
aii xi +ω aij xj = (1 − ω)aii xj −ω aij xj + ωbi (29)
j=1 j=i+1

pour i = 1, 2, . . . , n. En utilisant les matrices D, L et U introduites

précédement, on peut écrire (29) sous la forme matricielle suivante:
(k+1) (k) (k)
(D + ωL)xi = (1 − ω)Dxj − ωU xj + ωb
(k)
= [(1 − ω)D − ωU ]xj + ωb (30)

40
Si les aii 6= 0 pour tout i, D + ωL est une matrice triangulaire
inférieure et donc det(D + ωL)=detD = Πni=1 aii 6= 0. La matrice
D + ωL est donc inversible et la relation (30) devient:
(k+1) (k)
xi = (D + ωL)−1 [(1 − ω)D − ωU ]xj + ω(D + ωL)−1 b (31)

c’est ce qui définit la méthode de relaxation.

On verra plus loin que le paramètre de relaxation ω doit être compris
entre 0 et 2.
Si ω < 1 la méthode est dite sous-relaxation.
Si ω > 1 la méthode est dite sur-relaxation.

41
Exemple

On reprend l’exemple précédent:


 3x1 + x2 − x3 = 8


x1 + 5x2 + x3 = −2 (32)



x1 − x2 + 8x3 = 4

On choisit x(0) = (0, 0, 0)t pour initialiser et trois valeurs de

ω = 0.5, 1.073, 1.2. Le résultat des 10 premières itérations est affiché
dans le tableau suivant:

42
 (k) (k) (k) kx(k+1) −x(k+1) k
k ω x1 x2 x3 kx(k+1) k

0 0.5 0 0 0 1
1 0.5 1.33333 -0.333333 0.145833 0.358932
2 0.5 2.07986 -0.589236 0.156098 0.167218
3 0.5 2.49749 -0.759976 0.124458 0.0849963
4 0.5 2.72948 -0.865382 0.0875498 0.0445988
5 0.5 2.8569 -0.927136 0.0572729 0.0235669
6 0.5 2.92585 -0.96188 0.0356533 0.012376
7 0.5 2.96251 -0.980757 0.0213723 0.00640525
8 0.5 2.98161 -0.990677 0.0124181 0.00324604
9 0.5 2.99132 -0.995712 0.00701942 0.00160023
10 0.5 2.99612 -0.99817 0.00386684 0.000760862
0 1.073 0 0 0 1
1 1.073 2.86133 -1.04324
43 0.0127988 0.0557175
Convergence et test darrêt

Test d’arrêt

Comme pour toute méthode itérative, on doit se préoccuper d’un test

d’arrêt. Soit x la solution du système Ax = b. On peut introduire le
vecteur erreur e(k) = x(k+1) − x et étudier son comportement au
cours des itérations.
Dans les deux cas, Jacobi et Gauss-Seidel, la méthode itérative est de
la forme
x(k+1) = Bx(k) + c
On aura donc

e(k) = x(k+1) − x = Bx(k) + c − Bx − c = Be(k−1) = B 2 e(k−2)

= . . . = B k e(0) (33)

44
La convergence est réalisée lorsque l’erreur e(k+1) est toujours
inférieur à e(k) , ou si limk→+∞ ke(k) k = 0.
Comme x(0) est arbitrairement choisi alors e(0) est aussi arbitraire,
on conclut qu’on a convergence de la méthode si limk→+∞ kB k k = 0.
Ils existent plusieurs tests pour arrêter les itérations. Elles sont tous
basés sur le vecteur derreur qui doit atteindre un critère prédéfini et
tendant vers une valeur proche de zéro.
On utilise le test darrêt basé sur:
kx(k+1) − x(k) k
<ǫ
kx(k) k
On peut aussi se contanter du critère darrêt suivant:

kx(k+1) − x(k) k < ǫ

45
Convergence

Comme on vient de voir, l’algorithme x(k+1) = Bx(k) + c converge si

limk→+∞ ke(k) k = 0 ce qui est équivalent à dire que
limk→+∞ kB k k = 0. Cette condition s’applique aussi bien à la
méthode de Jacobi qu’à la méthode de Gauss-Seidel. Une condition
nécessaire et suffisante pour que limk→+∞ kB k k = 0 est que le rayon
spectrale de B vérifie ρ(B) = maxi=1,...,n |λi | < 1. Où les λi ,
i = 1, . . . , n sont les valeurs propres de B. On a donc les théorèmes
suivants:

46
Théorème: Pour tout choix de x(0) ∈ Rn , la suite x(k) définie par
x(k+1) = Bx(k) + c pour tout k ≥ 1 converge vers l’unique solution de
x = Bx + c si et seulement si ρ(B) < 1.

Exemple: On reprend l’exemple dy système:


 4x + x = 5
1 2
(34)
 x1 − 5x2 = −4

La matrice de la méthode de Jacobi est B = D −1 (L + U ) avec

     
4 0 0 0 0 1
D=   , L=   , U=  
0 −5 1 0 0 0
 
−1 0 1/4
B = D (L + U ) =   (35)
−1/5 0

47
 √
Les valeurs propres de B sont ±i/(2 5), le rayon spectral de B est
√
ρ(B) = 1/(2 5) < 1, la méthode de Jacobi est bien convergente dans
ce cas.
La matrice de la méthode
 de Gass-Seidel
 est
−1
0 1/4
B = (D + L) U =  . Les valeurs propres de B sont: 0,
0 1/20
1/20. Le rayon spectral de B est ρ(B) = 1/20 < 1. Comme nous
l’avons vu, la méthode de Gauss-Seidel est bien convergente dans ce
cas.

48
Corollaire: Si kBk < 1 pour une norme matricielle naturelle et
c ∈ Rn , alors pour tout x(0) ∈ Rn la méthode itérative
x(k+1) = Bx(k) + c , pour tout k ≥ 1, converge vers x et on les
relations suivantes:

(i) kx(k) − xk ≤ kBkk kx(0) − xk

k
kBk
(ii) kx(k) − xk ≤ kx(1) − x(0) k (36)
1 − kBk

49
Définition: Une matrice A est dite à diagonale strictement
dominante si
n
X
|aii | > |aij | , ∀i = 1, . . . , n (37)
j=1,i6=j

 
3 1 1
 
Exemples: 1.) La matrice:  1 4 −1 

 est à diagonale
2 3 8
strictement dominante car: 
3 > 1 + 1, 4 > |1| + | − 1| et 8 > 2 + 3.
 4x + x = 5
1 2
2.) La matrice du système est aussi à diagonale
 x1 − 5x2 = −4
strictement dominante car: 4 > 1 et 5 > 1.

50
Théorème: Si la matrice du système Ax = b est à diagonale
strictement dominante alors les méthodes de Jacobi et de
Gauss-Seidel convergent vers l’unique solution de Ax = b pour
n’importe quelle valeur d’initialisation x(0) .

Exemples: 
 4x + x = 5
1 2
1.) La matrice du système est à diagonale
 x1 − 5x2 = −4
strictement dominante, donc la méthode de Jacobi étudiée
précédement est convergente.

 3x1 + x2 − x3 = 8


2.) La matrice du système x1 + 5x2 + x3 = −2 est à diagonale



x1 − x2 + 8x3 = 4
strictement dominante, donc la méthode de Gauss-Seidel étudiée

51
précédement est convergente.
3.) La condition “à diagonale strictement dominante” est suffisante
pour garantir la convergence de la méthode de Jacobi et
Gauss-Seidel, mais elle n’est pas necessaire, c’est à dire on peut avoir
des systèmes avec des matrices qui ne sont pas à diagonale
strictement dominante mais la méthode itérative de Jacobi et/ou
Gauss-Seidel converge. C’est ce qu’illustre l’exemple suivant:

 −3x1 + x2 = −4


−x1 − 2x2 − 2x3 = −1



x2 − 2x3 = −3

La matrice du système n’est pas à diagonale strictement dominante

et pourtant les méthodes de Jacobi et Gauss-Seidel convergent.
La matrice de la méthode de Jacobi est

52
  
0 −1/3 0
−1
 
BJ = D (L + U ) = 
 1/2 0 . Le spectre de la matrice
1 
0 −1/2 0
√ √
BJ est: σ(BJ ) = 0, ±i √23 .
Le rayon spectral est ρ(BJ ) = √2
3
< 1,
donc la méthode de Jacobi converge.
Partant de x(0) = (0, 0, 0)t , on donne dans le tableau suivant les
(k)
résultats des ittération xi , k = 1, . . . , 10, 20, 40 ainsi que l’erreur
kx(k+1) −x(k+1) k
relative: kx(k+1) k
en norme infinie.

53
 (k) (k) (k) kx(k+1) −x(k+1) k
k x1 x2 x3 kx(k+1) k

0 0 0 0 1
1 1.333333333 0.5000000000 1.500000000 1.238095238
2 1.500000000 -1.666666667 1.750000000 0.5416666667
3 0.7777777778 -2.000000000 0.6666666667 2.166666667
4 0.6666666667 -0.5555555556 0.5000000000 0.5909090909
5 1.148148148 -0.3333333333 1.222222222 0.7222222222
6 1.222222222 -1.296296296 1.333333333 0.3333333333
7 0.9012345679 -1.444444444 0.8518518519 0.7536231884
8 0.8518518519 -0.8024691358 0.7777777778 0.2921348315
9 1.065843621 -0.7037037037 1.098765432 0.3727598564
10 1.098765432 -1.131687243 1.148148148 0.1786941581
20 0.9869938526 -0.9826584701 0.9804907788 0.02793620306
40 0.9997744535 -0.9996992713
54 0.9996616803 0.0004881632774
On remarque√que la convergence est très lente, car le rayon spectral
est ρ(BJ ) = √23 ≈ 0.82 < 1 qui est proche de 1.
De même, la matrice de la méthode de Gauss-Seidel est
 
0 −1/3 0
−1
 
BG = (D + L) U =  0 1/6
 1 . Le spectre de la matrice
0 1/12 1/2
BG est: σ(BG ) = 2/3, 0, 0, le rayon spectral est ρ(BG ) = 2/3 < 1,
donc la méthode de Gauss-Seidel converge.

55
 (k) (k) (k) kx(k+1) −x(k+1) k
k x1 x2 x3 kx(k+1) k

0 0 0 0 1
1 1.333333333 -0.1666666667 1.416666667 0.8928571429
2 1.277777778 -1.555555556 0.7222222222 0.7812500000
3 0.8148148148 -0.6296296296 1.185185185 0.4950495050
4 1.123456790 -1.246913580 0.8765432099 0.3802281369
5 0.9176954733 -0.8353909465 1.082304527 0.2472187886
6 1.054869684 -1.109739369 0.9451303155 0.1764446405
7 0.9634202103 -0.9268404207 1.036579790 0.1162621712
8 1.024386526 -1.048773053 0.9756134736 0.07998800180
9 0.9837423157 -0.9674846314 1.016257684 0.05304248248
10 1.010838456 -1.021676912 0.9891615438 0.03586901084
20 1.000187955 -1.000375911 0.9998120446 0.0006258746386
40 1.000000057 -1.000000113
56 0.9999999435 1.882050576 10−7
Dans le cas de la méthode de Gauss-Seidel, la convergence est plus
vite que celle de Jacobi mais est relativement lente aussi.

57
Théorème: Si A est une matrice symétrique, définie positive et
tridiagonale, alors les méthodes de Jacobi, Gauss-Seidel et de
ralaxation avec 0 < ω < 2 convergent pour tout x(0) ∈ Rn . Pour la
méthode de relaxation, le choix optimal de ω est donné par:
2
ωo = p
1+ 1 − ρ(BJ )2

Exemple: soit le système suivant


 2x1 + x2 = 2


x1 + 2x2 + x3 = 3



x2 + 2x3 = 2

58
la matrice du système est:
 
2 1 0
 
A
 1 2 1 


1 2

Cette matrice est symétrique, définie positive et tridiagonale, donc les

méthodes Jacobi, Gauss-Seidel et de la relaxation convergent.
 
0 1/2 0
−1
 
La matrice de Jacobi est: BJ = D (L + U ) =  1/2 0 1/2 .
0 1/2 0
√
Le spectre de cette matrice est ±1/ 2, 0, Le rayon spectrale est
√
1/ 2 < 1, la méthode converge.
La matrice de la méthode de Gauss-Seidel est

59
  
0 1/2 0
−1
 
BG = (D + L) U =   0 −1/4 1/2 . Le spectre de la matrice


0 1/8 −1/4
BG est: σ(BG ) = −1/2, 0, 0, le rayon spectral est ρ(BG ) = 1/2 < 1,
donc la méthode de Gauss-Seidel converge.
Pour la méthode de ralaxation, le choix optimal de ω est
2 2
ωo = p = q √ = 1.17
1 + 1 − ρ(BJ ) 2
1 + 1 − (1/ 2)2

60
Conditionnement d’une matrice

Exemple

Considérons le système linéaire suivant:

     
 x + 3x = 4 1 3 x1 4
1 2
⇔     =   ⇔ Ax =
 1.0001x1 + 3x2 = 4.0001 1.0001 3 x2 4.0001

La solution exacte et unique du système est x = (1, 1)t .

Considérons le sysème originalavec une 

petite
perturbation
 du
′
4 0
second membre b = b + δb =   +  . Le
4.0001 0.0004

61
système devient donc

 x + 3x = 4
1 2
 1.0001x1 + 3x2 = 4.0004

et x̃ = (4, 0) est la solution exacte.

On remarque que pour kb − b̃k∞ = 0.0003 on a kx − x̃k∞ = 3. Ce qui
se traduit en terme d’erreur relative sur x et b par:

kb − b̃k 0.0003
= = 0.7510−4 (38)
kbk 4.0001
kx − x̃k 3
= =3 (39)
kxk 1
pour une petite erreur relative sur le second membre du système, on
a une grande erreur relative sur x.

62
Pour comprendre ce phénomène, on a tracé sur la figure (1) les
droites d1 d’équation x1 + 3x2 = 4 et la droite d2 d’équation
1.0001x1 + 3x2 = 4.0001. Evidement, l’unique solution du système
est le point d’intersection de d1 et d2 . Après perturbation du second
membre, les droites d˜1 d’équation x1 + 3x2 = 4 et la droite d˜2
d’équation 1.0001x1 + 3x2 = 4.0004 seront complètement confondues
avec d1 et d2 . Cependant, la solution exacte x̃ = (4, 0)t est sur la
droite d1 , qui est très proche de d2 . Ce qui implique que x̃ est aussi
sur la droite d2 (voir figure (1)).

63
 d2 4
d1 3
2
1

1 2 3 4

Figure 1: La solution du système est le point d’intersection de d1 et d2

Cet exemple nous montre bien les difficultés qui peuvent apparaı̂tre.
Pour savoir quand ce genre de phénomème pourra apparaı̂tre, il est
necessaire de considérer la norme de la matrice A et celle de A−1 .

64
Ceci nous conduit à introduire la notion de conditionnement d’une
matrice.

65
Conditionnement d’une matrice

Théorème: Soient Ax = b un système linéaire, δb une petite

perturbation du second membre b et δx la perturbation obtenue sur
la solution x. Alors on a:

(i) kδxk ≤ kA−1 kkδbk

1 kδbk kδxk −1 kδbk
(ii) ≤ ≤ kA kkAk
kAkkA−1 k kbk kxk kbk
ou k.k norme matricielle naturelle.
Dans le cas d’une perturbation δA sur la matrice du système,
(A + δA)(x + δx) = b, alors
kδxk kδAk
(iii) ≤ kA−1 kkAk
kx + δxk kAk

Preuve:

66
i) On a A(x + δx) = b + δb = Ax + A(δx) ce qui donne A(δx) = δb et
donc δx = A−1 δb. En passant à la norme:
kδxk = kA−1 δbk ≤ kA−1 kkδbk d’où le résultat.
ii) Puisque b = Ax est équivalent à x = A−1 b, alors on a
kbk = kAxk ≤ kAkkxk et kxk ≤ kA−1 kkbk. On a aussi
kδxk =≤ kA−1 kkδbk et kδbk =≤ kAkkδxk donc:

−1 kδxk kA−1 kkδbk

kδxk ≤ kA kkδbk ⇒ ≤
kAkkxk kbk
kδxk −1 kδbk
⇒ ≤ kAkkA k (40)
kxk kbk
d’où l’inégalité à droite. L’inégalité à gauche se démontre de la même
façon.
iii) Si (A + δA)(x + δx) = b alors Aδx + δA(x + δx) = 0. Donc
δx = −A−1 (δA(x + δx)). On passe a la norme:

67
kδxk = kA−1 (A(x + δx)k ≤ kA−1 kkδAkkx + δxk. Ce qui donne:
kδxk
≤ kA−1 kkδAk
kx + δxk
−1 kδAk
≤ kAkkA k (41)
kAk

68
Définition: On appelle conditionnement dune matrice carrée A,
relatif à une norme naturelle, le nombre réel

κ(A) = kAkkA−1 k

Remarque: kIk = kAA−1 k ≤ kAkkA−1 k = κ(A), donc κ(A) ≥ 1.

69
Définition: Un système est dit:
• bien conditionnési κ(A) est proche de 1
• mal conditionné si κ(A) est grand.
Soit AX=B un système linéaire, de ce qui précéde on a les relations
suivantes:
kδxk kδbk
≤ κ(A) , perturbation du second membre(42)
kxk kbk
kδxk kδAk
≤ κ(A) perturbation de la matrice (43)
kx + δxk kAk
Ses deux relations nous montrent que si une matrice A est mal
conditionnée, κ(A) est grand, alors des petites variations sur les
données A ou b entraı̂nent de très fortes variations sur le résultat x
équations (42,43).
Le conditionnement est donc l’outil mathématique permettant de

70
quantifier l’instabilité numérique des résolutions des systèmes
linéaires.
Propriétés :
(i) κ(αA) = κ(A) pour toute matrice A et tout scalaire α non nul.
(ii) κ(αA) = κ(αA−1 )
(iii) κ(AB) ≤ κ(A)κ(B)
(iv) Si A est une matrice symétrique définie positive alors si Λ et λ
représentent respectivement la plus grande et la plus petite valeur
propre de A alors κ(A) pour la norme 2 est égale à: κ(A) = Λ/λ
(v) Si A est une matrice unitaire (AA∗ = A∗ A = I) ou orthogonale
réelle (AAt = At A = I) alors pour la norme k.k2 on a: κ(A) = 1.
Remarques:
1.) La propriété (v) montre que les systèmes à matrice unitaire ou

71
orthogonale sont bien conditionnés.
2.) La propriété (iv) est utile pour trouver une valeur approchée de
κ(A) sans passer par le calcul de la norme de A et A−1 , surtout que
le calcul de A−1 est très coûteux O(n3 ).
Exemple: La matrice du système considérée dans l’exemple
précédent est  
1 3
A=  
1.0001 3
sa norme est kAk∞ = 4.0001, cette norme ne sera pas considerée très
grande. Cependant, la norme de A−1 est:
 
−1 −10000 10000
A =   donc kA−1 k∞ = 20000
3333.67 −3333.33

donc pour la norme infinie: κ(A) = (20000)(4.0001) = 80002 qui est

72
un grand nombre.
kδbk
Comme nous l’avons vu dans l’exemple si kbk = 0.7510−4 , à partir
kδxk kδxk
de la majoration de kxk donnée par kxk ≤ κ(A) kδbk
kbk , on trouve:
kδxk kδxk
kxk ≤ 6.00015. On a calculé kxk = 3 ce qui est consistent avec ce
kδxk
qu’on a trouvé à partir de la majoration de kxk .

73