Faculté des Sciences – Département I&I

Master EEEA 2020/2021

Méthodes Numériques et Optimisation

Travaux Pratiques n°2
(Résolution des Systèmes Linéaires et Régression linéaire)

Ces travaux pratiques sont composés de deux parties. La première partie est dédiée à la
résolution des systèmes linéaires de la forme Ax  bS  . Il est à noter que la résolution du
système ( S) peut être effectuée par des méthodes directes (inversion de la matrice A , méthode
d’élimination de Gauss, etc.) ou à travers des algorithmes numériques. Ici, nous allons explorer
l’algorithme de Gauss-Seidel. Dans la deuxième partie, on aborde le problème de la régression
linéaire en vue de la modélisation avec application à la pile à combustible.

Partie I
Travail à effectuer

I.1. Méthode de Gauss-Seidel

Ecrire une fonction qui réalise l’algorithme fondé sur la méthode de Gauss-Seidel suivant :

ALGO Gauss-Seidel C’EST

VAR x : Tableau ( 1 n ) DE REEL, A, B, D, L,U : Tableau ( 1 n,1 n ) DE REEL,
b : Tableau ( 1n ) DE REEL,  : REEL
Début
/*Lecture de la matrice A  aij , le vecteur b du système (S) et de la précision souhaitée  */
Saisir ( A ), Saisir ( b ), Saisir ( e ) /*Exemple pour la précision   0.0001 */
/*Initialisation des matrices D, L,U */
a11 00  0 00  0 a21a1n 
0a  0  a 00  00 
D  22  , L  21  ,U    ,
    a( n1) n 
     
000ann  an1an ( n 1) 0 00 
B  D  L  U , f  D  L  b ,
1 1

/*Corps de l’algorithme*/
Tant que ( X k  1  X k    ) Faire
X k 1  B X k   f
k  k 1
Fin Tant que
Fin Gauss-Seidel
I.2. Exemple numérique
En utilisant la fonction définie dans la question précédente, résoudre le système suivant :
10 x1  3 x2  x3  14

2 x1  10 x2  3 x3  5
 x  3 x  10 x  14
 1 2 3
Donner le résultat des trois premières itérations ainsi que le résultat final. Notez que vous
pouvez utiliser les fonctions « diag », « triu », « tril » pour constuire D, L, et U.

I.3. Circuit électrique

En utilisant les loi classique de Kirchhoff trouver les valeurs des différents courants du circuit
électrique donné par la figure 1. Pour ce faire, écrire les six équations des courants, puis,
résoudre le système en employant la fonction obtenue dans la première question.

4Ω 1Ω 5Ω

i1
i4 i5 i6
i3 5V
3Ω

5Ω

1Ω

2Ω

3Ω
10V

i2

2Ω 6Ω
Figure 1. Circuit électrique

I.4. Problème de grande taille

Considérons le système
Ax  bS 
avec
4 1 0 0  1 
1 4 1 0   2 
  
0 1 4 1   3 
A  et B 
   
0 1  499 
   
0 1 4 500500 500  5001

Pour manipuler les matrices creuses, nous disposons sous Matlab des deux fonctions « sparse »
et « full ». Par ailleurs, les fonctions « tic » et « tac » permettent d’estimer le temps de calcul.
Utiliser la fonction « Help » pour les définitions exactes ainsi que pour la syntaxe.

I.4.1

I.4.1.1 Exprimer A  As en utilisant par utilisation de la fonction « sparse ».

I.4.1.2 Exprimer A  A f en utilisant par utilisation de la fonction « full ».
I.4.2

I.4.2.1 Résoudre le système (S ) pour A  As en utilisant X  As \ b . Puis, estimer le temps de

calcul.
I.4.2.2 Résoudre le système (S ) pour A  A f en utilisant X  Af \ b . Puis, estimer le temps de
calcul.

I.4.3 Comparer les temps de calcul pour la question précédente. Qu’en déduisez-vous ?

I.4.4 Déterminer le minimum et le maximum du vecteur solution x .

Partie II
y

x
Figure 2. Modèle linéaire

Considérons un nuage de points (voir par exemple la figure 2) de IR2 : M i  xi , yi ,1  i  n .
Ces données sont souvent le résultat de mesures ou on cherche à décrire le comportement global
de ce nuage. Par exemple, dans la figure 1 les points ne sont pas alignés mais on décide de
chercher une droite les « approchant » au mieux, soit f une fonction linéaire
f x   a x  b .
On utilise pour cela la méthode des moindres carrés. Plus précisément, comme on n’a pas
exactement
yi  axi  b
pour tout i {1,, n} , l’objectif alors est de minimiser le carré des différences :
n
Er f     f xi   yi 
2

i 1
On cherche donc un couple de réels (a, b) qui soit solution du problème d’optimisation suivant:
min J a, b 
(P) 
a, b   IR
2
 n
où J (a, b)    yi  a xi  b  .
2

i 1
Ce problème peut être généralisé de deux façons. D’une part, en considérant que la fonction
recherchée f n’est pas linéaire (voir figure 2) mais plutôt polynomiale. D’autre part dans le cas
où la variable x est de dimension d  1 .

x
Figure 3. Modèle non linéaire
Travail à effectuer

II.1. Régression linéaire

On définit
n n n n
Sx   xi , S y   yi , S xy   xi yi , S x 2   xi2 .
i 1 i 1 i 1 i 1

II.1.1 Montrer que

S y S x  n S xy S xy S x  S y S x 2
 
J a* , b*  0  a* 
 S x 2  n S x 2
et b* 
 S x 2  n S x 2


II.1.2 Déterminer la matrice Hessienne D 2 J a* , b* en fonction de S x , S x2 et n .
II.1.3 En admettant que pour tout réel  tel que
1 
  min  S x2 ,1
3 
Nous avons
a a
D 2 J (a* , b* )   ,     (a 2  b 2 )
b  b 
En déduire que (a* , b* ) réalise le minimum.

II.1.4 Ecrire une fonction Matlab intitulée «reglin1D» qui permet le calcul de (a* , b* ) à partir
d’un nuage de 20 points. Puis tester la fonction dans le cas où le nuage de point est généré par
la fonction

g z   2 z  3  N
avec z  [0,10] , le pas est égal à 0.5 et où N est un bruit gaussien de variance égale 2 (you can
make a figure to see if you are right).

Indication. Utiliser la fonction «normrnd» de Matlab pour générer ce bruit gaussien

II.2. Application à la pile à combustible

La pile à combustible (voir la figure 4) est un procédé qui permet de transférer l’énergie
chimique en électricité. Il s’agit d’un processus prometteur pour remplacer les énergies fossiles.

Figure 4. Schéma de fonctionnement de la pile à combustible

La tension de la pile exprimée en fonction de la densité de courant est donné par

(M) V  A  Bi  Ci 2  D lg(i)
Des essais expérimentaux, sur une pile PEMFC à 500 W, donnent des valeurs des tensions
correspondantes à des courants. Dans le fichier «[Link]»la première colonne représente la
densité de courant et la deuxième la tension.

II.3.1 En utilisant la fonction «lsqcuvefit» de Matlab déterminer les valeurs pour les constantes
A, B, C, D (Initial values : [ A, B, C, D]  [1,1,1,1] ).
II.3.2 Dessiner sur la même figures deux courbes qui décrivent le comportement de la tension.
La première déduite directement des essais expérimentaux. La deuxième où la tension est
calculée à partir du modèle (M) introduit ci-dessus.