0% ont trouvé ce document utile (0 vote)
265 vues133 pages

Programmation Des Méthodes Numériques Utilisant MATLAB)

Programmation des Méthodes Numériques Utilisant MATLAB)

Transféré par

Mustapha El Metoui
Copyright
© © All Rights Reserved
Nous prenons très au sérieux les droits relatifs au contenu. Si vous pensez qu’il s’agit de votre contenu, signalez une atteinte au droit d’auteur ici.
Formats disponibles
Téléchargez aux formats PDF, TXT ou lisez en ligne sur Scribd
0% ont trouvé ce document utile (0 vote)
265 vues133 pages

Programmation Des Méthodes Numériques Utilisant MATLAB)

Programmation des Méthodes Numériques Utilisant MATLAB)

Transféré par

Mustapha El Metoui
Copyright
© © All Rights Reserved
Nous prenons très au sérieux les droits relatifs au contenu. Si vous pensez qu’il s’agit de votre contenu, signalez une atteinte au droit d’auteur ici.
Formats disponibles
Téléchargez aux formats PDF, TXT ou lisez en ligne sur Scribd

-

République Algérienne Démocratique et Populaire


Ministère de l'Enseignement Supérieur et de la Recherche Scientifique

Polycopie de Cours, TD et TP
Méthodes Numérique avec MA TLAB
Master en Électrotechnique Options :
- Commande et Actionneurs Électriques.
- Réseaux Électriques.

Présenté par :

Dr. ALLAOUA Boumediene


(Maître de conférences classe A à l'Université Tahri Mohamed BÉCHAR)

2015

Université Tahri Mohamed Béchar


BP 417, Route de L'indépendance, Béchar (08000), Algérie, ~ : [Link], ~: [Link]
----~~======~---=====================--~====~~=====,~~~--
Sommaire

____________________________________________________

Sommaire
____________________________________________

Sommaire ............................................................................................................... i
Introduction Générale............................................................................................1

Chapitre 1 : Initiation au langage MATLAB


1.1. Introduction et Historique ..............................................................................4
1.2. Démarrage de MATLAB .................................................................................6
1.2.1. L’accès MATLAB .......................................................................................6
1.2.2. Fenêtres .....................................................................................................6
1.2.3. Démonstrations et l’aide dans MATLAB .................................................7
1.3. Commandes générales ....................................................................................8
1.3.1. Gestion des fichiers ...................................................................................8
1.3.2. Calculs élémentaires .................................................................................8
1.3.3. Constantes prédéfinies .............................................................................8
1.3.4. Historique ..................................................................................................8
1.3.5. Variables d’environnement .......................................................................9
1.4. Types de données MATLAB ...........................................................................9
1.4.1. Le type complexe .......................................................................................9
1.4.2. Le type chaîne de caractères ..................................................................10
1.4.3. Le type logique ........................................................................................11
1.5. Les vecteurs...................................................................................................12
1.5.1. Définir un vecteur ...................................................................................12
1.5.2. Vecteurs spéciaux....................................................................................13
1.5.3. Opérations vectorielles ...........................................................................13
1.5.4. Création rapide d’un vecteur ..................................................................14
1.6. Matrices et tableaux .....................................................................................14
1.6.1. Opérations sur les Tableaux (Matrices) .................................................15
[Link]. Addition et soustraction....................................................................15
[Link]. Multiplication, division et puissance terme à terme .......................15
[Link]. Multiplication, division et puissance au sens matriciel ..................15
[Link]. Transposition.....................................................................................15
1.6.2. Matrices spéciales ...................................................................................16
1.6.3. Manipulations sur les Matrice ...............................................................16
1.7. Les polynômes ...............................................................................................17
1.8. Opérateurs relationnels, logiques et mathématiques .................................18
1.9. Fonctions MATLAB ......................................................................................18
1.10. Entrées-sorties ............................................................................................19
1.11. M-Files ou scripts ........................................................................................20
1.12. Les Boucles (Structures itératives) ............................................................20
1.12.1. La boucle for .........................................................................................20
1.12.2. La boucle while.....................................................................................21

[Link] Boumediene -i-


Sommaire

1.13. Les Testes (Structures conditionnelles) .................................................... 21


1.13.1. La structure if...................................................................................... 21
1.13.2. La structure switch ............................................................................. 22
1.14. Graphismes et gestion des fenêtres graphiques ....................................... 23
1.14.1. Graphique 2D ...................................................................................... 23
[Link]. L’instruction « plot et fplot » (Tracer une courbe simple) ........... 23
[Link]. Identification des axes et des points, légendes, grilles ............... 24
[Link] Tracés multiples superposés et côte à côte ................................... 25
1.14.2. Options visuelles ................................................................................. 26
1.14.3. Échelles ................................................................................................ 27
1.14.4. Graphique 3D ...................................................................................... 28

Chapitre 2 : Résolution de l’équation non linéaire f(x) = 0


2.1. Introduction .................................................................................................. 30
2.2. Méthode de dichotomie (Bissection) ............................................................ 30
2.3. Théorème du Point Fixe (Approximations ou Itérations ............................ 32
2.3.1. Définition ................................................................................................ 32
2.3.2. Proposition .............................................................................................. 32
2.3.3. Théorème ................................................................................................ 32
2.4. Méthode de la sécante (Lagrange) ............................................................... 33
2.5. Méthode de Newton (des Tangentes) .......................................................... 34
2.5.1. Description de la méthode ..................................................................... 35
2.5.2. Théorème ................................................................................................ 35
2.5.3. Interprétation graphique ....................................................................... 35
2.6. Exemple d’application .................................................................................. 36
2.7. Applications par code MATLAB sur les méthodes numériques
de la résolution de l’équation non linéaire f(x) = 0 .................................... 38
2.7.1. Applications pour valider les programmes............................................ 38
2.7.2. Application MATLAB de la méthode de Dichotomie ............................ 38
[Link]. Algorithme de la méthode de Dichotomie ....................................... 38
[Link]. Programme MATLAB de la méthode de Dichotomie ...................... 38
2.7.3. Application MATLAB de la méthode de Point fixe ............................... 40
[Link]. Algorithme de la méthode de Point fixe .......................................... 40
[Link]. Programme MATLAB de la méthode de Point fixe......................... 40
2.7.4. Application MATLAB de la méthode de Sécante (Lagrange) ............... 41
[Link]. Algorithme de la méthode de Sécante ............................................. 41
[Link]. Programme MATLAB de la méthode de Sécante............................ 41
2.7.5. Application MATLAB de la méthode de Newton .................................. 42
[Link]. Algorithme de la méthode de Newton ............................................. 42
[Link]. Programme MATLAB de la méthode de Newton............................ 42

Chapitre 3 : Résolution d’un Système d’Équations Linéaires


3.1. Introduction .................................................................................................. 45
3.2. Généralités ................................................................................................... 45
3.3. Problème ....................................................................................................... 46
3.4. Méthode Directe d’Elimination de Gauss................................................... 46
3.5. Exemples sur la méthode d’élimination de Gauss ...................................... 49

[Link] Boumediene -ii-


Sommaire

3.6. Méthode Itérative de Gauss-Seidel ..............................................................52


3.7. Exemple sur la méthode de Gauss-Seidel ....................................................53
3.8. Applications par code MATLAB sur les méthodes numériques..................54
de la résolution d’un Système d’Équations Linéaires .................................54
3.8.1. Applications pour valider les programmes ............................................54
3.8.2. Application MATLAB de la méthode d’élimination de Gauss...............54
[Link]. Algorithme de la méthode d’élimination de Gauss ..........................54
[Link]. Programme MATLAB de la méthode d’élimination de Gauss ........55
3.8.3. Application MATLAB de la méthode de Gauss-Seidel ..........................56
[Link]. Algorithme de la méthode de Gauss-Seidel .....................................56
[Link]. Programme MATLAB de la méthode de Gauss-Seidel ....................57

Chapitre 4 : Interpolation polynômiale


4.1. Introduction...................................................................................................59
4.2. Généralités ....................................................................................................59
4.3. Matrice de Vandermonde..............................................................................59
4.4. Interpolation de Lagrange ............................................................................60
4.5. Polynôme de Newton.....................................................................................61
4.6. Exemple d’application ...................................................................................62
4.7. Applications par code MATLAB sur les méthodes de calcul numériques
de l’interpolation polynômiale ......................................................................64
4.7.1. Applications pour valider les programmes ............................................64
4.7.2. Application MATLAB de la matrice de Vandermonde ..........................64
[Link]. Algorithme de la matrice de Vandermonde .....................................64
[Link]. Programme MATLAB de la matrice de Vandermonde....................64
4.7.3. Application MATLAB de l’interpolation de Lagrange ...........................65
[Link]. Algorithme de l’interpolation de Lagrange ......................................65
[Link]. Programme MATLAB de l’interpolation de Lagrange ....................65
4.7.4. Application MATLAB du polynôme de Newton .....................................66
[Link]. Algorithme du polynôme de Newton ................................................66
[Link]. Programme MATLAB du polynôme de Newton ..............................67

Chapitre 5 : Intégration Numérique


5.1. Introduction...................................................................................................69
5.2. Définitions .....................................................................................................69
5.3. Calcul de l’intégrale par la formule de Newton-Côtes.................................69
5.4. Calcul de l’intégrale par la formule des Trapèzes généralisée....................70
5.5. Calcul de l’intégrale par la formule de Simpson généralisée ......................71
5.6. Exemple d’application ...................................................................................71
5.7. Applications MATLAB sur calcul numérique de l’intégrale .......................73
5.7.1. Applications pour valider les programmes ............................................73
5.7.2. Application MATLAB de la formule de Newton-Côtes ..........................73
[Link]. Algorithme de la formule de Newton-Côtes .....................................73
[Link]. Programme MATLAB de la formule de Newton-Côtes ...................73
5.7.3. Application MATLAB de la formule des Trapèzes généralisée.............74
[Link]. Algorithme de la formule des Trapèzes généralisée ........................74
[Link]. Programme MATLAB de la formule des Trapèzes généralisée ......74

[Link] Boumediene -iii-


Sommaire

5.7.4. Application MATLAB de la formule de Simpson généralisée .............. 74


[Link]. Algorithme de la formule de Simpson généralisée ......................... 74
[Link]. Programme MATLAB de la formule de Simpson généralisée ........ 75

Chapitre 6 : Calcul Numérique des Valeurs Propres


6.1. Introduction .................................................................................................. 77
6.2. Définitions .................................................................................................... 77
6.3. Méthode de Déterminant ............................................................................. 77
6.4. Exemple d’application utilisant la Méthode de Déterminant .................... 78
6.5. Méthode de Le Verrier ................................................................................. 78
6.6. Exemple d’application utilisant la Méthode de Le Verrier ........................ 79
6.7. Méthode de Krylov ....................................................................................... 80
6.8. Exemple d’application utilisant la Méthode de Krylov .............................. 81
6.9. Applications MATLAB sur le calcul des valeurs propres ........................... 82
6.9.1. Applications pour valider les programmes............................................ 82
6.9.2. Application MATLAB de la méthode de Le Verrier.............................. 82
[Link]. Algorithme de la méthode de Le Verrier ......................................... 82
[Link]. Programme MATLAB de la méthode de Le Verrier ....................... 83
6.9.3. Application MATLAB de la méthode de Krylov .................................... 84
[Link]. Algorithme de la méthode de Krylov ............................................... 84
[Link]. Programme MATLAB de la méthode de Krylov ............................. 84

Chapitre 7 : Résolution d’un Système d’Équations Non


Linéaires
7.1. Introduction .................................................................................................. 86
7.2. Méthode de Newton-Raphson ...................................................................... 86
7.3. Exemple d’application utilisant la méthode de Newton-Raphson ............. 87
7.4. Applications MATLAB sur la résolution d’un système d’équations
non linéaires ................................................................................................. 90
7.4.1. Applications pour valider les programmes............................................ 90
7.4.2. Application MATLAB de la méthode de Newton-Raphson .................. 90
[Link]. Algorithme de la méthode de Newton-Raphson.............................. 90
[Link]. Programme MATLAB de la méthode de Newton-Raphson ............ 90
7.4.3. Programme intégré dans MATLAB pour RSENL ................................ 91

Chapitre 8 : Résolution d’Équations Différentielles


8.1. Introduction .................................................................................................. 93
8.2. Définitions .................................................................................................... 93
8.3. Méthode de Runge-Kutta ............................................................................. 94
8.4. Exemple d’application sur la méthode de Runge-Kutta ............................. 95
8.5. Méthode d’Euler ........................................................................................... 96
8.6. Exemple d’application sur la méthode d'Euler ........................................... 97
8.7. Méthodes de Taylor ...................................................................................... 98
8.8. Exemple d’application sur la méthode de Taylor ........................................ 99
8.9. Applications MATLAB sur la résolution d’équations Différentielles ...... 100
8.9.1. Application MATLAB de la méthode de Runge-Kutta ....................... 100

[Link] Boumediene -iv-


Sommaire

[Link]. Algorithme de la méthode de Runge-Kutta ................................... 100


[Link]. Programme MATLAB de la méthode de Runge-Kutta .................. 100
8.9.2. Application MATLAB de la méthode d'Euler ...................................... 101
[Link]. Algorithme de la méthode d'Euler .................................................. 101
[Link]. Programme MATLAB de la méthode d'Euler ................................ 101
8.9.3. Application MATLAB de la méthode de Taylor ................................... 102
[Link]. Algorithme de la méthode de Taylor .............................................. 102
[Link]. Programme MATLAB de la méthode de Taylor ............................ 102

Conclusion Générale .......................................................................................... 104


Références .......................................................................................................... 105
Annexe .............................................................................................................. 107

[Link] Boumediene -v-


Introduction Générale

[Link] Boumediene
Introduction

Introduction
Cet ouvrage s’adresse aux étudiants en Master d’électrotechnique ou en
formation ingénierie. Le but de ce cours était de présenter aux étudiants
quelques notions théoriques de base concernant les méthodes numériques
permettant de résoudre effectivement de tels problèmes. C’est pour cette
raison qu’une part importante du cours est consacrée à la mise en place d’un
certain nombre de techniques fondamentales de l’analyse numérique :
interpolation polynomiale, intégration numérique, résolution d’un système
d’équations linéaires ou non linéaires à une et plusieurs variables. Son objectif
est de donner au lecteur un outil lui permettant de travailler de manière
autonome à l’aide de questions détaillées et progressives, et d’une construction
pas à pas des programmes.

Ce choix de faire de la théorie avant de commencer la programmation est


indispensable pour appréhender les notions d’analyse numérique mais aussi
pour améliorer ses capacités de programmeur ; la programmation demande un
peu d’âme. Cette préparation ne dispense pas d’une réflexion sur la manière de
programmer une méthode. À cette fin, les exemples en MATLAB proposent
une programmation sous forme résumer. À chaque question, le programme
précédent est amélioré et complété. Les résultats intermédiaires sont donnés
pour valider cette programmation par morceaux. Dans chaque chapitre, les
solutions complètes et les programmes sont systématiquement donnés.

Le premier chapitre rappelle les commandes utiles de MATLAB et initiation


au langage pour gérer des tableaux et les commandes élémentaires. Il s’agit
donc de savoir utiliser au mieux les tableaux et de s’initier aux premières
commandes du graphisme. Ces commandes nécessaires sont insuffisantes pour
progresser dans la programmation et il ne faut pas hésiter à consulter
fréquemment l’aide de MATLAB.

La résolution de l’équation non linéaire f(x)=0 est montrée dans le deuxième


chapitre par les méthodes les plus utilisées. La programmation et l’application
de ces méthodes avec MATLAB est installé à la fin du chaque chapitre.

Le troisième chapitre est abordé la notion des méthodes de résolution de


système d’équations linéaires. Néanmoins les quelques méthodes proposés
peuvent servir de base.

Le chapitre 4 est consacré à l’interpolation ou comment faire passer une


courbe par des données mesurées. Les réponses sont données par l’interpolant
de Lagrange, Vandermonde et le polynôme de Newton.

Le cinquième et le sixième chapitre présentent le calcul numérique de :


l’intégrale et les valeurs propres.

Le chapitre 7 est consacré à la résolution d’un système d’équations non


linéaires.

_______________________________________________________________________________
[Link] Boumediene -1-
Introduction

Enfin, dans le dernier chapitre, on trouvera la résolution d’équations


différentielles.

Nous remercions nos collègues de l’Université de Béchar qui ont participés à


l’élaboration ou à la correction de ce document et pour les remarques et
améliorations constantes suggérées tout au long de notre collaboration.

_______________________________________________________________________________
[Link] Boumediene -2-
Chapitre 1

Initiation au langage
MATLAB

[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

1.1. Introduction et Historique

MATLAB est une abréviation de MATrix LABoratory. Écrit à l’origine, en


Fortran, par C++. Moler, MATLAB était destiné à faciliter l’accès au logiciel
matriciel développé dans les projets LINPACK et EISPACK. La version actuelle,
écrite en C par de compagnie mondiale « MathWorks ». (pour consultation voir le
site web [Link] existe en version professionnelle et en
version étudiant. Sa disponibilité est assurée sur plusieurs plateformes : Sun,
Bull, HP, IBM, compatibles PC (DOS, Unix ou Windows), Macintoch, iMac et
plusieurs machines parallèles [1, 2].

MATLAB est un environnement puissant, complet et facile à utiliser destiné au


calcul scientifique. Il apporte aux ingénieurs, chercheurs et à tout scientifique un
système interactif intégrant calcul numérique et visualisation. C’est un
environnement performant, ouvert et programmable qui permet de remarquables
gains de productivité et de créativité [3].

MATLAB est un environnement complet, ouvert et extensible pour le calcul et la


visualisation. Il dispose de plusieurs centaines (voire milliers, selon les versions
et les modules optionnels autour du noyau MATLAB) de fonctions
mathématiques, scientifiques et techniques [3, 4]. L’approche matricielle de
MATLAB permet de traiter les données sans aucune limitation de taille et de
réaliser des calculs numériques et symboliques de façon fiable et rapide. Grâce
aux fonctions graphiques de MATLAB, il devient très facile de modifier
interactivement les différents paramètres des graphiques pour les adapter selon
nos souhaits [4].

L’approche ouverte de MATLAB permet de construire un outil sur mesure. On


peut inspecter le code source et les algorithmes des bibliothèques de fonctions
(Toolboxes), modifier des fonctions existantes et ajouter d’autres.

MATLAB possède son propre langage, intuitif et naturel qui permet des gains de
temps de CPU spectaculaires par rapport à des langages comme le C, le
TurboPascal et le Fortran. Avec MATLAB, on peut faire des liaisons de façon
dynamique, à des programmes C ou Fortran, échanger des données avec d’autres
applications ou utiliser MATLAB comme moteur d’analyse et de visualisation [3].

MATLAB comprend aussi un ensemble d’outils spécifiques à des domaines,


appelés Toolboxes (ou Boîtes à Outils). Indispensables à la plupart des
utilisateurs, les Boîtes à Outils sont des collections de fonctions qui étendent
l’environnement MATLAB pour résoudre des catégories spécifiques de problèmes.
Les domaines couverts sont très variés et comprennent notamment le traitement
du signal, l’automatique, l’identification de systèmes, les réseaux de neurones, la
logique floue, le calcul de structure, les statistiques, etc [1, 5].

MATLAB fait également partie d’un ensemble d’outils intégrés dédiés au


Traitement du Signal.

_______________________________________________________________________________
-4-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

En complément du noyau de calcul MATLAB, l’environnement comprend des


modules optionnels qui sont parfaitement intégrés à l’ensemble [5]:
1. Une vaste gamme de bibliothèques de fonctions spécialisées (Toolboxes);
2. Simulink, un environnement puissant de modélisation basée sur les
schémas-blocs et de simulation de systèmes dynamiques linéaires et non
linéaires ;
3. Des bibliothèques de blocs Simulink spécialisés (Blocksets) ;
4. D’autres modules dont un Compilateur, un générateur de code C, un
accélérateur,... ;
5. Un ensemble d’outils intégrés dédiés au Traitement du Signal : le DSP
Workshop.

· Quelles sont les particularités de MATLAB ?

MATLAB permet le travail interactif soit en mode commande, soit en mode


programmation ; tout en ayant toujours la possibilité de faire des visualisations
graphiques. Considéré comme un des meilleurs langages de programmations (C
ou Fortran), MATLAB possède les particularités suivantes par rapport à ces
langages :

o la programmation facile,
o la continuité parmi les valeurs entières, réelles et complexes,
o la gamme étendue des nombres et leurs précisions,
o la bibliothèque mathématique très compréhensive,
o l’outil graphique qui inclut les fonctions d’interface graphique et les
utilitaires,
o la possibilité de liaison avec les autres langages classiques de
programmations (C ou Fortran).

Dans MATLAB, aucune déclaration n’est à effectuer sur les nombres. En effet, il
n’existe pas de distinction entre les nombres entiers, les nombres réels, les
nombres complexes et la simple ou double précision. Cette caractéristique rend le
mode de programmation très facile et très rapide [6, 7].

En Fortran par exemple, une subroutine est presque nécessaire pour chaque
variable simple ou double précision, entière, réelle ou complexe. Dans MATLAB,
aucune nécessité n’est demandée pour la séparation de ces variables [8].

La bibliothèque des fonctions mathématiques dans MATLAB donne des analyses


mathématiques très simples. En effet, l’utilisateur peut exécuter dans le mode
commande n’importe quelle fonction mathématique se trouvant dans la
bibliothèque sans avoir à recourir à la programmation [9].

Pour l’interface graphique, des représentations scientifiques et même artistiques


des objets peuvent être créées sur l’écran en utilisant les expressions
mathématiques. Les graphiques sur MATLAB sont simples et attirent l’attention
des utilisateurs, vu les possibilités importantes offertes par ce logiciel [4].

· MATLAB peut-il s’en passer de la nécessité de Fortran ou du C ?

_______________________________________________________________________________
-5-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

La réponse est non. En effet, le Fortran ou le C sont des langages importants


pour les calculs de haute performance qui nécessitent une grande mémoire et un
temps de calcul très long. Sans compilateur, les calculs sur MATLAB sont
relativement lents par rapport au Fortran ou au C si les programmes comportent
des boucles. Il est donc conseillé d’éviter les boucles, surtout si celles-ci est grande
[2, 5].

1.2. Démarrage de MATLAB

1.2.1. L’accès MATLAB

Pour lancer l’exécution de MATLAB, cliquer deux fois sur l’icône MATLAB [2, 3].
Le caractère ‘>>‘ signifie que MATLAB attend une instruction (syntaxe).

Pour lancer l’exécution de MATLAB :


ü sous Windows, il faut cliquer sur Démarrage, ensuite Programme, ensuite
MATLAB,
ü sous d’autres systèmes, se référer au manuel d’installation.
L’invite ‘>>‘ de MATLAB doit alors apparaître, à la suite duquel on entrera les
commandes.
La fonction "exit" permet de quitter MATLAB :
>> exit

1.2.2. Fenêtres

_______________________________________________________________________________
-6-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

ü Command Window

L’une des plus importantes fenêtres de MATLAB, la Command Window traite


des instructions données. Les résultats s’afficheront dès le retour de ligne.
On peut éviter de retaper une instruction en pesant sur la flèche vers le haut (↑),
ce qui affichera l’instruction précédente. Continuez à peser ↑ jusqu’à l’instruction
recherchée.

ü Workspace

La fenêtre nommée Workspace permet de visualiser les variables mises en


mémoire. On y retrouve leur nom, leurs dimensions ainsi que le type de variable.
MATLAB étant basé sur les matrices, toutes les variables sont constituées de
plusieurs dimensions : un scalaire est une matrice 1×1 et un vecteur est une
matrice 1×n ou n×1, etc. Il est possible d’effacer certaines variables ainsi que de
les éditer. Pour toutes les effacer, utilisez la commande Clear Workspace dans le
menu Edit.

ü Current Folder

Le Current Folder est le répertoire courant où sont enregistrés les fichiers-M.

ü Command History

Le Command History inscrit les commandes au fur et à mesure qu’elles sont


appelées dans le Command Window. Elle garde en mémoire ces commandes,
ainsi que la date et l’heure d’ouverture de chacune des sessions de MATLAB [3].

1.2.3. Démonstrations et l’aide dans MATLAB

Pour des démonstrations pas à pas sur l’environnement de MATLAB, ouvrez


MATLAB, puis entrez demo dans la fenêtre de commande [2].
Une aide interactive est disponible pour toutes les commandes MATLAB. Pour
obtenir la liste des commandes reprises dans cette aide, exécuter la commande
help Pour obtenir l’aide relative à une commande, exécuter help commande. Par
exemple ‘help abs’ fournit la description de la fonction abs (valeur absolue).
La commande "help" permet de donner l’aide sur un problème donné, par
exemple (la fonction Cosinus):

>> help cos

COS(X) is the cosine of the elements of X.

help → donne de l’aide sur une fonction.


which → localise fonctions et fichiers.
what → liste des fichiers MATLAB dans le répertoire courant.
exist → check si une fonction ou une variable existe dans le workspace.
who,whos → liste des variables contenues dans le workspace.
clear all → détruit toutes les variables de l’espace de travail.

_______________________________________________________________________________
-7-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

intro lance une introduction à MATLAB.


help produit une liste de toutes les commandes par thèmes.
démonstration donnant une représentation des fonctionnalités
demos
de bases de MATLAB.
info information sur la boite à outils disponibles.
ouvre une fenêtre contenant la liste des commandes MATLAB
helpwin
ainsi que leurs documentations.
help donne la liste de toutes les commandes par thèmes.
Help nom décrit la fonction nom.m

1.3. Commandes générales

1.3.1. Gestion des fichiers

Vous pouvez utiliser la petite fenêtre en haut à droite, ou à défaut :

pwd affiche le nom du répertoire courant pour MATLAB.


cd rep change le répertoire courant pour MATLAB qui devient rep.
dir fournit le catalogue d’un répertoire.
delete efface des fichiers ou des objets graphiques.

1.3.2. Calculs élémentaires

Dans la partie commandes de l’interface :


>> 5+8
Résultat : >> 13
Pour conserver le résultat, il faut l’assigner dans un objet :
>> a=5+8
>> a
Pour ne pas faire afficher le résultat, mettez « ; » à la fin de la commande :
>> a=5+8;

1.3.3. Constantes prédéfinies

pi 3.1415...
eps 2.2204e-016
Inf nombre infini
NaN Pas de solution, exprime parfois une indétermination

1.3.4. Historique

MATLAB conserve l’historique des commandes. Il est donc possible de récupérer


des instructions déjà saisies (et ensuite de les modifier dans le but de les
réutiliser) :
↑, ↓, ←, → permet de se déplacer dans les lignes de commandes tapées dans la
fenêtre de commandes

_______________________________________________________________________________
-8-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

1.3.5. Variables d’environnement

MATLAB garde en mémoire les variables qui ont été créées. On les voit en haut,
à gauche, lorsque MATLAB dispose d’une interface graphique. Sinon, on peut les
afficher et les effacer par la ligne de commande :

who donne la liste des variables présentes dans l’espace de travail


donne la liste des variables présentes dans l’espace de travail
whos
ainsi que leurs propriétés
donne la liste des fichiers .m et .mat présents dans le
what
répertoire courant
clear efface toutes les variables crées dans l’espace de travail
ans réponse retournée après exécution d’une commande.

1.4. Types de données MATLAB

MATLAB traite un seul1 type d’objet: les matrices !. Les scalaires sont des
matrices 1×1, les vecteurs lignes des matrices 1×n. les vecteurs colonnes des
matrices n×1 [2, 3].

Les trois principaux types de variables utilisés par MATLAB sont les types réels,
complexe et chaîne de caractères. Il n’y a pas de type entier à proprement parler.
Le type logique est associé au résultat de certaines fonctions. Signalons qu’il est
inutile (impossible) de déclarer le type d’une variable. Ce type est établi
automatiquement à partir des valeurs affectées à la variable. Par exemple les
instructions x = 2 ; z = 2+i ; rep = ‘oui’ ; définissent une variable x de type
réel, une variable z de type complexe et une variable rep de type chaîne de
caractères.

>> clear
>> x = 2; z = 2+i; rep = ‘oui’;
>> whos

Name Size Bytes Class


rep 1x3 6 char array
x 1x1 8 double array
z 1x1 16 double array (complex)

1.4.1. Le type complexe

L’unité imaginaire est désignée par i ou j. Les nombres complexes peuvent être
écrits sous forme cartésienne a+ib. Les différentes écritures possibles sont a+ib,
a+i*b, a+b*i, a+bi et r*exp(it) ou r*exp(i*t) avec a, b, r et t des
variables de type réel. Les commandes imag, real, abs, angle permettent de
passer aisément de la forme polaire à la forme cartésienne et réciproquement. Si
z est de type complexe, les instructions imag(z) et real(z) retournent la partie
imaginaire et la partie réelle de z. Les instructions abs(z) et angle(z)
retournent le module et l’argument de z.

_______________________________________________________________________________
-9-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

On fera attention au fait que les identificateurs i et j ne sont pas réservés. Aussi
il est possible que des variables de noms i et j aient été redéfinies au cours d’un
calcul antérieur et soient toujours actives.

>> z = [1+i, 2, 3i]


z =
1.0000 + 1.0000i 2.0000 0 + 3.0000i

>> y = [1+i, 2, 3 i]
y =
1.0000 + 1.0000i 2.0000 3.0000 0 + 1.0000i

1.4.2. Le type chaîne de caractères

Une chaîne de caractères est un tableau de caractères. Une donnée de type


chaîne de caractères (char) est représentée sous la forme d’une suite de
caractères encadrée d’apostrophes simples (‘). Une variable de type chaîne de
caractères étant interprétée comme un tableau de caractères, il est possible de
manipuler chaque lettre de la chaîne en faisant référence à sa position dans la
chaîne. L’exemple suivant présente différentes manipulations d’une chaîne de
caractères.

>> ch1 = ‘bon’


ch1 =
bon
>> ch2 = ‘jour’
ch2 =
jour
>> whos
Name Size Bytes Class
ch1 1x3 6 char array
ch2 1x4 8 char array
>> ch = [ch1,ch2]
ans =
bonjour
>> ch(1), ch(7), ch(1:3)
ans =
b
ans =
r
ans =
bon
>> ch3 = ‘soi’;
>> ch = [ch(1:3), ch3, ch(7)]
ans =
bonsoir

Si une chaîne de caractères doit contenir le caractère apostrophe (‘) celui-ci doit
être doublé dans la chaîne (ainsi pour affecter le caractère apostrophe (‘) à une

_______________________________________________________________________________
-10-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

variable on devra écrire ’’’’, soit 4 apostrophes). La chaîne de caractères vide


s’obtient par deux apostrophes ’’.
>> rep = ‘aujourd’hui’
??? rep = ‘aujourd’hui
|
Missing operator, comma, or semi-colon.

>> rep = ‘aujourd’’hui’


rep =
aujourd’hui
>> apos = ‘‘‘‘
apos =

1.4.3. Le type logique

Le type logique (logical) possède deux formes : 0 pour faux et 1 pour vrai. Un
résultat de type logique est retourné par certaines fonctions ou dans le cas de
certains tests. Dans l’exemple qui suit on considère une variable x contenant la
valeur 123 et une variable y définie par l’expression mathématique y =
exp(log(x)). On teste si les variables x et y contiennent les mêmes valeurs. La
variable tst est une variable de type logique qui vaut 1 (vrai) les valeurs sont
égales et 0 (faux) sinon. Suivant la valeur de tst, on affiche la phrase x est egal a
y ou la phrase x est different de y. Dans l’exemple proposé, compte-tenu des
erreurs d’arrondis lors du calcul de exp(log(123)), la valeur de la variable y
ne vaut pas exactement 123. On pourra également considérer le cas où x = 12.

>> x = 123; y = exp(log(x));


>> tst = ( x==y );

>> if tst, disp(‘x est egal a y ‘), else disp(‘x est different
de y ‘), end x est different de y

>> whos

Name Size Bytes Class


tst 1x1 8 double array (logical)
x 1x1 8 double array
y 1x1 8 double array

>> format long


>> x, y, tst
x =
123
y =
1.229999999999999e+02
tst =
0

_______________________________________________________________________________
-11-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

1.5. Les vecteurs

1.5.1. Définir un vecteur

On définit un vecteur ligne en donnant la liste de ses éléments entre crochets


([ ]). Les éléments sont séparés au choix par des espaces ou par des virgules.
On définit un vecteur colonne en donnant la liste de ses éléments séparés au
choix par des points virgules ( ; ) ou par des retours chariots (touche
Entrée/Enter). On peut transformer un vecteur ligne x en un vecteur colonne et
réciproquement en tapant x’ (’ est le symbole de transposition) [1, 2, 3].
Il est inutile de définir la longueur d’un vecteur au préalable. Cette longueur sera
établie automatiquement à partir de l’expression mathématique définissant le
vecteur ou à partir des données. On peut obtenir la longueur d’un vecteur donné
grâce à la commande length.

Un vecteur peut également être défini « par blocs » selon la même syntaxe. Si par
exemple x1, x2 et x3 sont trois vecteurs (on note x1, x2 et x3 les variables
MATLAB correspondantes), on définit le vecteur X = (x1 | x2 | x3) par
l’instruction X = [x1 x2 x3].

>> x1 = [1 2 3], x2 = [4,5,6,7], x3 = [8; 9; 10]


x1 =
1 2 3
x2 =
4 5 6 7
x3 =
8
9
10
>> length(x2), length(x3)
ans =
4
ans =
3
>> x3’
ans =
8 9 10
>> X = [x1 x2 x3’]
X =
1 2 3 4 5 6 7 8 9 10

Les éléments d’un vecteur peuvent être manipulés grâce à leur indice dans le
tableau. Le kème élément du vecteur x est désignée par x(k). Le premier
élément d’un vecteur a obligatoirement pour indice 1.

Reprenons l’exemple précédent.

>> X(5)
ans =
5

_______________________________________________________________________________
-12-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

>> X(4:10)
ans =
4 5 6 7 8 9 10
>> X([Link])
ans =
2 4 6 8 10
>> K = [1 3 4 6]; X(K)
ans =
1 3 4 6

Il est très facile de définir un vecteur dont les composantes forment une suite
arithmétique. Pour définir un vecteur x dont les composantes forment une suite
arithmétique de raison h, de premier terme a et de dernier terme b, on écrira
x = a:h:b.

>> x = 1.1:0.1:1.9
x =
Columns 1 through 7
1.1000 1.2000 1.3000 1.4000 1.5000 1.6000 1.7000

Columns 8 through 9
1.8000 1.9000
>> x = 1.1:0.2:2
x =
1.1000 1.3000 1.5000 1.7000 1.9000

1.5.2. Vecteurs spéciaux

Les commandes ones, zeros et rand permettent de définir des vecteurs dont les
éléments ont respectivement pour valeurs 0, 1 et des nombres générés de
manière aléatoire.

ones(1,n) vecteur ligne de longueur n dont tous les éléments valent 1.


ones(m,1) vecteur colonne de longueur m dont tous les éléments valent 1.
zeros(1,n) vecteur ligne de longueur n dont tous les éléments valent 0.
zeros(m,1) vecteur colonne de longueur m dont tous les éléments valent 0.
rand(1,n) vecteur ligne de longueur n dont les éléments sont générés de
manière aléatoire entre 0 et 1.
rand(m,1) vecteur colonne de longueur m dont les éléments sont générés
de manière aléatoire entre 0 et 1.

1.5.3. Opérations vectorielles

n:m nombres de n à m par pas de 1


n:p:m nombres de n à m par pas de p
linspace(n,m,p) p nombres de n à m
lenght(x) longueur de x
x(i) i-éme coordonnée de x
x(i1:i2) coordonnées i1 à i2 de x

_______________________________________________________________________________
-13-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

x(i1:i2)=[] supprimer les coordonnées i1 à i2 de x


[x,y] concaténer les vecteurs x et y
x*y' produit scalaire des vecteurs lignes x et y
x'*y produit scalaire des vecteurs colonnes x et y

1.5.4. Création rapide d’un vecteur

Certaines commandes permettent de créer plus rapidement des vecteurs précis :

> 11=1:10 (Un vecteur contenant les entiers de 1 à 10)


> 12=[Link]
> 13=10:-1:1
> 14=1:0.3:pi
> 11(2)=13(3)
> 14(3:5)=[1,2,3]
> 14(3:5) = []
> 15 = linspace(i,5,5)
> help linspace
> who
> whos
> clear l1 12 13 15
> who
> clc (efface te contenu de la fenêtre de commande)
> clear

· Remarque : Une ligne de commande commençant par le caractère « % » n'est


pas exécutée par MATLAB. Cela permet d’insérer des lignes de commentaires.

1.6. Matrices et tableaux

MATLAB ne fait pas de différence entre les Matrices et tableaux. Le concept de


tableau est important car il est à la base du graphique : typiquement pour une
courbe de « n » points, on définira un tableau de « n » abscisses et un tableau de
« n » ordonnées. Mais on peut aussi définir des tableaux rectangulaires à deux
indices pour définir des matrices au sens mathématique du terme, puis effectuer
des opérations sur ces matrices [3, 4].

On utilise les crochets [ et ] pour définir le début et la fin de la matrice. Ainsi


é1 3ù
pour définir une variable A contenant la matrice : A = ê ú
ë4 2û
on écrira : A = [ 1 3 ; 4 2 ].

D’une façon générale, on définit une matrice en donnant la liste de ses éléments
entre crochets. Signalons que MATLAB admet d’autres façons d’écrire les
matrices. Les éléments d’une ligne de la matrice peuvent être séparés au choix
par un blanc ou bien par une virgule (,). Les lignes quant à elles peuvent être
séparées au choix par le point-virgule (;) ou par un retour chariot. Par exemple,
on peut aussi écrire la matrice A de la manière suivante :

_______________________________________________________________________________
-14-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

>> A = [1,3;4,2] et de colonne. A(i,j) désigne le ième


A = élément ligne de la jème élément
1 3 colonne de la matrice A. Ainsi A(2,1)
4 2 désigne le premier élément de la
>> A = [1 3 deuxième ligne de A,
4 2]
A = >> A(2,1)
1 3 ans =
4 2 4
>> A = [1,3 La commande size permet d’obtenir
4,2]
les dimensions d’une matrice A
A =
donnée.
1 3
4 2
>> size(A)
Un élément d’une matrice est
référencé par ses numéros de ligne

1.6.1. Opérations sur les Tableaux (Matrices)

[Link]. Addition et soustraction

Les deux opérateurs sont les mêmes que pour les scalaires. A partir du
moment où les deux tableaux concernés ont la même taille, Le tableau
résultant est obtenu en ajoutant ou soustrayant les termes de chaque tableau.

[Link]. Multiplication, division et puissance terme à terme

Ces opérateurs sont notés .*, ./ et .^ (attention à ne pas oublier le point). Ils
sont prévus pour effectuer des opérations termes à terme sur deux tableau de
même taille. Ces symboles sont fondamentaux lorsque l’on veut tracer des
courbes.

[Link]. Multiplication, division et puissance au sens matriciel

Puisque l’on peut manipuler des matrices, il parait intéressant de disposer


d’une multiplication matricielle. Celle-ci se note simplement * et ne doit pas
être confondu avec la multiplication terme à terme. Il va de soi que si l’on écrit
A*B le nombre de colonnes de A doit être égal au nombre de lignes de B pour
que la multiplication fonctionne.
La division a un sens vis-à-vis des inverses de matrices. Ainsi A/B représente A
multiplié (au sens des matrices) à la matrice inverse de B.

[Link]. Transposition

L’opérateur transposition est le caractère «‘ » et est souvent utilisé pour


transformer des vecteurs lignes en vecteurs colonnes et inversement.

_______________________________________________________________________________
-15-
[Link] Boumediene
Chapitre 1 Initiation au langage MATLAB

1.6.2. Matrices spéciales

zeros(m,n) matrice nulle de taille m, n


ones(m,n) matrice de taille m, n dont tous les coefficients valent 1
eye(n) matrice identité de taille n
diag(x) matrice diagonale dont la diagonale est le vecteur x
magie(n) matrice carré magique de taille n

1.6.3. Manipulations sur les Matrice

ligne colonne

Matrice (début : fin, début : fin)

size(A) Nombre de lignes et colonnes de A


A(i,j) Coefficient d'ordre i , j de A
A ( i1:i2 , :) Lignes i1 à i2 de A
A ( i1:i2 , :)=[] Supprimer les lignes il à i2 de A
A ( : , j1:j2 ) Colonnes j1 à j2 de À
A ( : , j1:j2 )=[] Supprimer les colonnes j1 à j2 de A
A (:) Concaténer les vecteurs colonnes de A
Coefficients diagonaux de A
diag(A)
Retournez en haut-bas
flipud(A)
Retournez droit-gauche
fliplr(A) 2 rotations de 90 degrés
rot90(A,2) Change la forme de la matrice
reshape(A,1,9) Extrait le triangle supérieur de A
triu(A) Extrait le triangle inférieur de A
tril(A) Transposée de A
A' Inverse de A
inv(A) Exponentielle de A
expm(A) Déterminant de A
det(A) Trace de A
trace(A) Polynôme caractéristique de A
Valeurs propres de A
polyCA)
Vecteurs propres et valeurs propres de
eig(A) A
[U,D]=eig(A) Addition, soustraction
+ - Multiplication, puissance (matricielles)
* ^ Multiplication, puissance ternie à
. * . ^ terme
A\b Solution de : Ax= b
b/A Solution de : xA =b
./ Division terme à terme

_______________________________________________________________________________
[Link] Boumediene -16-
Chapitre 1 Initiation au langage MATLAB

· Entrez la matrice : >> A=[l 2 3 ; 2 3 1 ; 3 1 2 ]


Quels sont les résultats des commandes suivantes ?
>> A([2 3] , [1 3])
>> A([2 3] , 1:2)
>> A([2 3] , :)
>> A([2 3],end)
>> A(:)

1.7. Les polynômes

Sous MATLAB le polynôme de degré n, p(x) = anxn + an−1xn−1 +… + a1x1 + a0 est


défini par un vecteur p de dimension n+1 contenant les coefficients {ai}i=0,...,n
rangés dans l’ordre décroissant des indices [2, 3].
La commande polyval permet d’évaluer le polynôme p (la fonction
polynômiale) en des points donnés. La syntaxe est polyval(p,x) où x est une
valeur numérique ou un vecteur. Dans le second cas on obtient un vecteur
contenant les valeurs de la fonction polynômiale aux différents points spécifiés
dans le vecteur x. Utilisée avec la commande fplot, la commande polyval
permet de tracer le graphe de la fonction polynômiale sur un intervalle [xmin,
xmax] donné.
La syntaxe de l’instruction est :
fplot(’polyval([a_n, ..., a_0] , x)’ , [x_min , x_max]).

Voici par exemple comment définir le polynôme p(x) = x2−1. Tracer le graphe
de la fonction polynômiale.

>> p = [1, 0, -1];


>> polyval(p,0)
ans =
-1
>> polyval(p,[-2,-1,0,1,2])
ans =
3 0 -1 0 3
>> fplot(’polyval([1, 0, -1] , x)’ , [-3,3]), grid

On obtient les racines du polynôme p grâce à l’instruction roots(p).


>> r = roots(p)
r =
-1.0000
1.0000
>> poly(r)
ans =
1.0000 0.0000 -1.0000

Les polynômes sont gérés, sous MATLAB, par des vecteurs de coefficients dans
l’ordre décroissant.
Aussi le polynôme : x5 + 2x4 - x2 -x + 1
est-il représenté sous MATLAB par le vecteur [1 2 0 -1 -1 1]. Prenons la
variable x :

_______________________________________________________________________________
[Link] Boumediene -17-
Chapitre 1 Initiation au langage MATLAB

x = -5 :0.01 :5 ;
puis taper :

polyval(p,x)
polyder(p)
polyder(p,x)

· Quelle est l’utilité des fonctions polyval et polyder ?

1.8. Opérateurs relationnels, logiques et mathématiques

Operateurs relationnels <, <=, >=, == (égalité), ~= (différent)


Opérateurs logiques & (et), | (ou), ~ ou not (non)
· Remarque : Attention de ne pas confondre = qui sert à affecter une valeur
à une variable et == qui sert à tester l’égalité.

1.9. Fonctions MATLAB

sum(x) somme des éléments du vecteur x.


prod(x) produit des éléments du vecteur x.
max(x) plus grand élément du vecteur x.
min(x) plus petit élément du vecteur x.
mean(x) moyenne des éléments du vecteur x.
sort(x) ordonne les éléments du vecteur x par ordre croissant.
fliplr(x) renverse l’ordre des éléments du vecteur x.
cross produit vectoriel.
dot produit scalaire.

· Les fonctions mathématiques incorporées sont :


log(x) logarithme népérien de x,
log10(x) logarithme en base 10 de x,
exp(x) exponentielle de x,
sqrt(x) racine carrée de x (s’obtient aussi par x.ˆ0.5),
abs(x) valeur absolue de x,
sign(x) valant 1 si x est strictement positif, 0 si x est nul et -1 si x
est strictement négatif.

· Les fonctions d’arrondis sont :


round(x) entier le plus proche de x,
floor(x) arrondi par défaut,
ceil(x) arrondi par excès,
fix(x) arrondi par défaut un réel positif et par excès un réel
négatif.

· Les fonctions trigonométriques et hyperboliques sont :


cos(x) cosinus,
acos(x) cosinus inverse (arccos),

_______________________________________________________________________________
[Link] Boumediene -18-
Chapitre 1 Initiation au langage MATLAB

sin(x) sinus,
asin(x) sinus inverse (arcsin),
tan(x) tangente,
atan(x) tangente inverse (arctan),
cosh(x) cosinus hyperbolique (ch),
acosh(x) cosinus hyperbolique inverse (argch),
sinh(x) sinus hyperbolique (sh),
asinh(x) sinus hyperbolique inverse (argsh),
tanh(x) tangente hyperbolique (th),
atanh(x) tangente hyperbolique inverse (argth).

1.10. Entrées-sorties

Deux commandes utiles pour gérer le workspace, dont la taille dépend de votre
espace de swap [3, 6]:
Pour sauvegarder une ou plusieurs variables il faut taper
save [NomFichier] var1 var2 [...] varn
ou pour sauvegarder l’espace de travail au complet, cliquer sur File -> Save
Workspace As...

>> save % écrit toutes les variables du workspace dans le fichier [Link]

Pour récupérer ces données, il suffit de taper


load [NomFichier]

>> load % charge dans le workspace toutes les variables du


% fichier [Link]
Il est possible également de sauvegarder des informations dans un fichier sous
forme de texte formaté. Pour cela, la première chose à faire est d’utiliser la
fonction fopen pour ouvrir un fichier :
f = fopen(’[NomFichier]’,’w’) en écriture,
f = fopen(’[NomFichier]’,’r’) en lecture.

Ensuite, utiliser la fonction fprintf pour écrire dans un fichier (avec le


format désiré), et la fonction fscanf pour lire un fichier. Se reporter à l’aide
MATLAB pour les options. Pour fermer le fichier, il suffit d’utiliser fclose.

Entrées-sorties sur des fichiers disques :


fopen (ouverture d’un fichier)
fclose (fermeture d’un fichier)
fscanf (lecture formatée)
fprintf (écriture formatée)
N = input(’Enter la valeur de N = ’) % entrée interactive
disp(N) % affiche la valeur de N

_______________________________________________________________________________
[Link] Boumediene -19-
Chapitre 1 Initiation au langage MATLAB

1.11. M-Files ou scripts

Un script (ou M-file) est un fichier (message.m par exemple) contenant des
instructions MATLAB [5, 6].
Voici un exemple de script:
% message.m affiche un message
% ce script affiche le message que s’il fait beau
beau_temps=1;
if beau_temps~=0
disp(’Hello, il fait beau’)
end
return % (pas nécessaire à la fin d’un M-file)

MATLAB vous offre un éditeur pour écrire et mettre au point vos M-files:

>> edit % lance l’éditeur de MatLab. Voir aussi la barre


d’outils

Tout autre éditeur de texte convient aussi.


Les M-files sont exécutés séquentiellement dans le “workspace”, c’est à dire
qu’ils peuvent accéder aux variables qui s’y trouvent déjà, les modifier, en
créer d’autres etc.
On exécute un M-file en utilisant le nom du script comme commande:

>> message
Hello, il fait beau

1.12. Les Boucles (Structures itératives)

Il y a deux types de boucles en MATLAB: la boucle for et la boucle while [3].

1.12.1. La boucle for


La boucle for parcourt un vecteur d'indices et effectue à chaque pas toutes les
instructions délimitées par l'instruction end.

for variable = valeur_initiale : pas : valeur_finale


Ici ce que tu veux faire
end

>> for k=1:4


x=k^2
end
x =
1
x =
4
x =
9
x =
16

_______________________________________________________________________________
[Link] Boumediene -20-
Chapitre 1 Initiation au langage MATLAB

1.12.2. La boucle while

La boucle while effectue une suite de commandes jusqu'à ce qu'une condition


soit satisfaite.

initialisations
while expression logique
Ici ce que tu veux faire
end

>> x=1
while x<14
x=x+5
end
x =
6
x =
11
x =
16

· Les deux boucles suivantes donne un même vecteur : f = w = 2 4 6 8 10

% Avec la boucle for :

for i = 1:5 % boucle for à pas unité


f(i) = 2*i ; % l’incrémentation de i est automatique
end

% Avec la boucle While :


i = 1; % initialisation de i
while i <= 5 % boucle tant que i est inférieur ou
% égal à 5
w(i) = 2*i;
i = i+1; % l’incrémentation est
end; % importante sinon la
% boucle est infinie

1.13. Les Testes (Structures conditionnelles)

Il existe deux structures possibles permettant de réaliser des tests sur les
données [3, 5].

1.13.1. La structure if

La structure if permet de tester la valeur d'une variable et d'effectuer


différents traitement suivant les cas testés.

_______________________________________________________________________________
[Link] Boumediene -21-
Chapitre 1 Initiation au langage MATLAB

if condition 1
commandes
elseif condition 2
commandes
elseif condition 3
commandes
.
.
.
else
commandes
end

% Exemple d'utilisation de if, elseif, ... end


nb = input('Rentrez un nombre : ');
if (nb < 100)
disp('votre nombre est < 100')
elseif (nb == 100)
disp('votre nombre est 100')
else
disp('votre nombre est > 100')
end

1.13.2. La structure switch

La structure switch permet de choisir entre différents cas, et de faire


correspondre un traitement adapté à chacun des cas reconnus.

switch variable_testé
case choix 1 (variable_testé avec choix 1)
commandes
case choix 2
commandes
case choix 3
commandes
.
.
.
otherwise
commandes
end

% Un exemple d'utilisation de switch


jour = input('Quel jour sommes nous ? ','s');
switch jour

_______________________________________________________________________________
[Link] Boumediene -22-
Chapitre 1 Initiation au langage MATLAB

case 'lundi'
jr = 5;
case 'mardi'
jr = 4;
case 'mercredi'
jr = 3;
case 'jeudi'
jr = 2;
case 'vendredi'
jr = 1;
otherwise
jr = 0;
end
disp([ 'Encore ' num2str(jr) ' jours avant le weekend' ])

1.14. Graphismes et gestion des fenêtres graphiques

À la différence des autres programmes permettant de dessiner des graphiques,


MATLAB est spécialisé dans les graphiques tracés à partir de matrices de
données [3, 4, 5].

Une instruction graphique ouvre une fenêtre dans laquelle est affiché le
résultat de cette commande. Par défaut, une nouvelle instruction graphique
sera affichée dans la même fenêtre et écrasera la figure précédente. On peut
ouvrir une nouvelle fenêtre graphique par la commande figure. Chaque
fenêtre se voit affecter un numéro. Ce numéro est visible dans le bandeau de la
fenêtre sous forme d’un titre. Le résultat d’une instruction graphique est par
défaut affiché dans la dernière fenêtre graphique ouverte qui est la fenêtre
graphique active [4].

On rend active une fenêtre graphique précédemment ouverte en exécutant la


commande figure(n), où n désigne le numéro de la figure. La commande
close permet de fermer la fenêtre graphique active. On ferme une fenêtre
graphique précédemment ouverte en exécutant la commande close(n), où n
désigne le numéro de la figure. Il est également possible de fermer toutes les
fenêtres graphiques en tapant close all.

1.14.1. Graphique 2D

[Link]. L’instruction « plot et fplot » (Tracer une courbe simple)

Pour tracer des courbes, il existe deux fonctions principales : plot et fplot.
La fonction plot(X,Y) trace la courbe dont Y en fonction de X.
C’est avec cette fonction qu’il faut tracer des courbes lorsque nous avons des
données expérimentales. Il suffit de les mettre dans une matrice [4].

X=[3 5.6 7.2 9.9 11.22] X1=-[Link]


Y=[8 -4 -2 0 5.67] Y1=X1.^3+2*X1.^2-8*X1+1
plot(X,Y) figure(1)
plot(X1,Y1)

_______________________________________________________________________________
[Link] Boumediene -23-
Chapitre 1 Initiation au langage MATLAB

La fonction fplot, elle s’écrit de la manière suivante : fplot(fon,lim) où


fon indique le nom de la fonction à tracer sous la forme d’une chaîne de
caractères (entre guillemets) et lim défini les limites des échelles. Pour les
limites, il faut les écrire entre crochet avec la syntaxe suivante : lim = [Xmin
Xmax Ymin Ymax].

f=’sin’;
fplot(f,[-2*pi 4*pi])

Lorsque vous avez plusieurs graphiques à tracer dans un même fichier, pour
que tous s’affichent en même temps, il faut indiquer l’instruction figure avant
chaque graphique à tracer afin de créer une nouvelle figure pour chacun d’eux.

x = 0 : .05 : 3*pi; y1 = sin(x); y2 = cos(x);


figure(1); plot(x,y1);
figure(2); plot(x,y2);

[Link]. Identification des axes et des points, légendes, grilles

Il est possible d’identifier les axes par les fonctions suivantes :


xlabel(‘texte’), ylabel(‘texte’). On peut titré le graphique grâce à
title(‘texte’) et ajouter une légende avec
legend(fonction1,fonction2,...,possont le texte à indiquer sur la fonction
désirée et où pos est la position où placer la légende [4].
· Positions :
0 = Meilleur placement automatique (où il y à le moins de conflits avec
les données),
1 = Coin supérieur à droite (défaut),
2 = Coin supérieur à gauche,
3 = Coin inférieur à gauche,
4 = Coin inférieur à droite,
-1 = À droite du graphique.

Il est aussi possible de quadriller le graphique avec la fonction grid. Cette


fonction s’applique à tous les types d’échelles (voir plus bas). Il est aussi
possible d’ajouter du texte sur le graphique de 2 manières : la fonction
text(x,y, ’text’) où x et y sont les coordonnées où placer le texte et la
fonction gtext(’text’) où le texte sera placé manuellement avec la sourie
lors de l’exécution du graphique.

De plus, on peut modifier les axes à l’aide des fonctions suivantes :


Axis([xmin xmax ymin ymax])pour les valeurs maximales et minimales
des axes ; Axis equal est pour les axes identiques.

_______________________________________________________________________________
[Link] Boumediene -24-
Chapitre 1 Initiation au langage MATLAB

figure
fplot(’[log(x)*sin(x)*3 log(x) sin(x)]’,[1 20 -10
10])
title(’Fonctions sinus et logarithme naturel’)
xlabel(’axe des X’)
ylabel(’axe des Y’)
gtext(’texte placé manuellement’)
text(11,9,’texte placé en avance’)
legend(’log*sin’,’log’,’sin’,0)
grid

[Link] Tracés multiples superposés et côte à côte

Pour obtenir plusieurs tracés superposés, il suffit de les inscrire dans la


fonction plot(x1,y1,x2,y2, ...). On peut aussi utiliser la fonction hold
on pour superposer des graphiques.

Les deux séries d’instructions suivantes produisent le même graphique :


x=0:.05:2*pi;
y1=sin(3*x);
y2=3*cos(3*x);

1. En utilisant la fonction plot(x1,y1,x2,y2,...)

plot(x,y1,x,y2)
legend(’f(x)’,’dérivé de f(x)’)
grid

2. En utilisant la fonction hold

plot(x,y1)
grid
hold on
plot(x,y2)
legend(’f(x)’,’dérivé de f(x)’)
hold off

Pour les fonctions superposées tracées à l’aide de la fonction fplot, il s’agit de


mettre entre rochets les fonctions désirées.

f=’[tan(x),sin(x),cos(x)]’;
fplot(f,2*pi*[-1 1 -1 1])

Pour obtenir plusieurs graphiques côte à côte, il suffit d’utiliser la fonction


subplot(m,n,p), ce qui divise la figure en une matrice de m par n carreaux

_______________________________________________________________________________
[Link] Boumediene -25-
Chapitre 1 Initiation au langage MATLAB

rectangulaires. P indique la position dans la matrice que prendra la fonction


plot qui suit dans les instructions.
Si on désire créer une matrice 2x2 pour les graphiques, la fonction subplot à
inscrire avant chaque fonction à afficher se défini ainsi :

subplot(2,2,1) subplot(2,2,2)

subplot(2,2,3) subplot(2,2,4)

t = 0:0.3:2*pi;
x = cos(t);
y = sin(t);
figure
subplot(2,1,1);
plot(t,x)
subplot(2,1,2);
plot(t,y)

1.14.2. Options visuelles

Pour tracer les courbes, plusieurs options sont offertes. L’usager peut fait le
choix entre le type de traits et de point utilisé, ainsi que la couleur. On utilise
la fonction plot(X,Y,S) où S indique le choix d’option entre guillemets [4].

_______________________________________________________________________________
[Link] Boumediene -26-
Chapitre 1 Initiation au langage MATLAB

figure X2 = 0:.5:5;
X1 = 0:.1:2*pi
plot(sin(X1),’g*-.’) plot(exp(X2),’md:’);

1.14.3. Échelles

Plusieurs fonctions pour tracer des graphiques sont disponibles pour les différentes
échelles [4]:
plot(x,y) linéaire x-y
polar(theta,r) polaire
loglog(x,y) logarithmique x-y
bar(x,y) barres
semilogx(x,y) logarithmique x, linéraire y
stairs(x,y) escalier
semilogy(x,y) linéaire x, logarithmique y
hist(y,x) dessine l’histogramme

x=-10:.5:10;
y=x.^2+x.*2+3;
subplot(2,2,1), polar(x,y), grid on
subplot(2,2,2), loglog(x,y),grid on
subplot(2,2,3), semilogx(x,y), grid on
subplot(2,2,4), semilogy(x,y), grid on

Il est possible de copier la figure une fois qu’elle est ouverte, il suffit d’aller dans le
menu Edition → Copy figure pour coller dans un fichier .doc (word)

· Exemple :

% Tracer les fonctions sinus et cosinus :


x=[-pi:0.1:+pi];
y=sin(x);
plot(x,y,'rp')
grid
z=cos(x)
plot(x,z,'gp')
grid

%utilisant subplot
subplot(1,2,1);
plot(x,y,'rp'), grid
subplot(1,2,2);
plot(x,z,'gp'), grid

title('Fonction')
xlabel('x');

_______________________________________________________________________________
[Link] Boumediene -27-
Chapitre 1 Initiation au langage MATLAB

Fonction Fonction
1 1

0.8 0.8

0.6 0.6

0.4 0.4

0.2 0.2
F(x)

F(x)
0 0

-0.2 -0.2

-0.4 -0.4

-0.6 -0.6

-0.8 -0.8

-1 -1
-4 -2 0 2 4 -4 -2 0 2 4
x x

figure 1.1 : Courbes les fonctions sinus et cosinus.

1.14.4. Graphique 3D

a) plot3

t = linspace(0, 10*pi);
plot3(sin(t), cos(t), t)
xlabel(’sin(t)’), ylabel(’cos(t)’), zlabel(’t’)
grid on

b) mesh

x = linspace(-3, 3, 30);
y = linspace(-3, 3, 30);
[X Y] = meshgrid(x,y);
Z = peaks(X, Y)
mesh(X, Y, Z)
meshc(X, Y, Z)

c) surfl

[X,Y,Z] = peaks(30);
surfl(X,Y,Z)
shading interp
colormap pink

d) contour

[X,Y,Z] = peaks(30)
contour(X,Y,Z,20)

_______________________________________________________________________________
[Link] Boumediene -28-
Chapitre 2

Résolution de l’équation
non linéaire f(x) = 0

[Link] Boumediene
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

2.1. Introduction

Le numéricien est souvent confronte a la résolution d'équations algébriques de


la forme f(x) = 0 et ce dans toutes sortes de contextes. Introduisons des
maintenant la terminologie qui nous sera utile pour traiter ce problème [10].

Une valeur de x solution de f(x) = 0 est appelée une racine ou un zéro de la


fonction f(x) et est notée α.

On considère une fonction réelle « f » définie sur un intervalle [a, b], avec a < b,
et continue sur cet intervalle et on suppose que f admet une unique racine sur
I =]a, b[, c’est-à-dire qu’il existe un unique α Î I tel que f(α) = 0.

2.2. Méthode de dichotomie (Bissection)

On considère un intervalle [a, b] et une fonction f continue de [a, b] dans R. On


suppose que f(a)·f(b) < 0 et que l’équation f(x) = 0 admet une unique solution α
sur l’intervalle [a, b].

La méthode de dichotomie Figure 2.1 consiste à construire une suite (xn) qui
converge vers α de la manière suivante :

y=f(x)

x0

x1 x2

a x3 b

Figure 2.1: Méthode de la Bissection.

Initialisation : on prend pour x0 le milieu de [a, b]. La racine se trouve alors


dans l’un des deux intervalles ]a, x0[ ou ]x0, b[ ou bien elle est égale à x0.

• si f(a)f(x0) < 0, alors α Î ]a, x0[. On pose a1 = a, b1 = x0.


• si f(a)f(x0) = 0, alors α = x0.
• si f(a)f(x0) > 0, alors α Î ]x0, b[. On pose a1 = x0, b1 = b.
On prend alors pour x1 le milieu de [a1, b1].
On construit ainsi une suite

x0 = (a+b)/2, x1 = (a1 + b1)/2, . . . , xn = (an + bn)/2 (2.1)

_______________________________________________________________________________
[Link] Boumediene -30-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

telle que

|α − xn| ≤ (b − a)/2 (2.2)

Etant donné une précision ε, cette méthode permet d’approcher α en un


nombre prévisible d’itérations.
Le principe de construction suivants consistent à transformer l’équation
f(x) = 0 en une équation équivalente g(x) = x. On peut poser par exemple
g(x) = x + f(x), mais on prendra plus généralement g(x) = x + u(x)·f(x) où « u »
est une fonction non nulle sur l’intervalle I [10, 11].

Il reste à choisir u pour que la suite définie par x0 Î I et la relation de


récurrence xn+1 = xn + u(xn)·f(xn) soit bien définie et converge vers la racine α de
f.

Géométriquement dans Figure 2.2, on a remplacé la recherche de l’intersection


du graphe de la fonction f avec l’axe OX, par la recherche de l’intersection de la
droite d’équation y = x avec la courbe d’équation y = g(x).
y=x
y=f(x)

y=g(x)

a α b a α b

(a)Racine de f (b) Point fixe de g

Figure 2.2: Recherche de l’intersection de la méthode de la Bissection.

Le choix d’une méthode est conditionné par les réponses aux questions
suivantes :

1) la suite (xn) converge-t-elle ?


2) si la suite converge, sa limite est-elle α ?
3) si on veut la solution à ε près, comment arrêter les itérations dès que cette
condition est remplie ?
4) comme dans tout calcul, on désire obtenir rapidement le résultat approché,
il faudra donc estimer la manière dont évolue l’erreur e = xn − α au cours des
itérations.

Les deux premières questions sont purement mathématiques. Les deux


dernières sont numériques, car on ne peut effectuer qu’un nombre fini
d’itérations pour le calcul.

_______________________________________________________________________________
[Link] Boumediene -31-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

La continuité des fonctions considérées permet de répondre à la question 2 : si


la suite converge, elle converge vers une racine de l’équation ; si, de plus, xn Î
I, pour tout n, alors par unicité de la racine dans I, la suite converge vers α.

2.3. Théorème du Point Fixe (Approximations ou Itérations


Successives)

2.3.1. Définition

On dit que l’application g : [a, b] Î R est strictement contractante si $ L


Î ]0, 1[, " x Î [a, b], " y Î [a, b], |g(x) − g(y)| ≤ L|x − y|.

2.3.2. Proposition

Soit g une application de classe C1 de l’intervalle [a, b] de R dans R. On


suppose que g′ vérifie : max{|g′(x)| ; x Î [a, b]} ≤ L < 1, alors l’application g est
strictement contractante dans l’intervalle [a, b].

2.3.3. Théorème

Si g est une application définie sur l’intervalle [a, b] à valeurs dans [a, b], alors
la suite (xn) définie par x0 Î [a, b] et la relation de récurrence xn+1 = g(xn)
converge vers l’unique solution α de l’´equation x = g(x) avec α Î [a, b].
Démonstration : la suite (xn) est bien définie car g([a, b]) Î [a, b].

Montrons, par l’absurde, que l’équation x = g(x) a au plus une solution.


Supposons qu’il y ait deux solutions α1 et α2, alors :
|α1 − α2| = |g(α1) − g(α2)| ≤ L|α1 − α2| (2.3)
or L < 1 donc nécessairement α1 = α2.

Montrons que la suite (xn) est convergente. On a :


|xn+1 − xn| = |g(xn) − g(xn−1| ≤ L|xn − xn−1| (2.4)
et par récurrence :
|xn+1 − xn| ≤ Ln|x1 − x0| (2.5)

On en déduit que
|xn+p − xn| ≤ Ln|x1 − x0|(1 − Lp)/(1 – L)≤ |x1 − x0|Ln/(1 – L) (2.6)

Cette suite vérifie le critère de Cauchy ; elle est donc convergente vers α Î [a,
b], or l’application g est continue donc la limite α vérifie g(α) = α.
On a aussi une évaluation de l’erreur en faisant tendre p vers l’infini, on
obtient
|α − xn| ≤ |x1 − x0|Ln/(1 – L) (2.7)

On constate que, pour n fixé, l’erreur est d’autant plus petite que L est proche
de 0.

_______________________________________________________________________________
[Link] Boumediene -32-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

2.4. Méthode de la sécante (Lagrange)

Cette méthode est également appelée méthode de Lagrange Figure 2.3,


méthode des parties proportionnelles ou encore regula falsi [11].
On considère un intervalle [a, b] et une fonction f de classe C 2 de [a, b] dans R.
On suppose que f(a)·f(b) < 0 et que f’ ne s’annule pas sur [a, b], alors l’équation
f(x) = 0 admet une unique solution sur l’intervalle [a, b].

· Remarque :
L’hypothèse « f dérivable » suffirait, mais demanderait une rédaction un peu
plus fine des démonstrations.

La méthode de la sécante consiste à construire une suite (xn) qui converge vers
α de la manière suivante : soit D 0 la droite passant par (a, f(a)) et (b, f(b)), elle
coupe l’axe OX en un point d’abscisse x0 Î ]a, b[. On approche donc la fonction f
par un polynôme P de degré 1 et on résout P(x) = 0.

(b, f(b))

D0

x0 x1
a a b
f (x 0 )
y = f(x)
(a, f(a))

Figure 2.3: Méthode de la Sécante.

Ensuite, suivant la position de α par rapport à x0, on considère la droite


passant par (a, f(a)) et (x0, f(x0)) si f(x0)·f(a) < 0 ou celle passant par (x0, f(x0)) et
(b, f(b)) si f(x0)·f(b) < 0.
On appelle x1 l’abscisse du point d’intersection de cette droite avec l’axe OX.
On réitère ensuite le procédé.

Plaçons-nous dans le cas où f′ > 0 est dérivable et f est convexe, alors la suite
(xn) est définie par :
ì x0 = a
ï
íx = b × f ( x n ) - x n × f ( b) (2.8)
ïî n +1 f (x n ) - f ( b)

En effet, si f est convexe, on remplace l’intervalle [xn, b] par l’intervalle [xn+1,


b]. L’équation d’une droite passant par (c, f(c)) et (b, f(b)) avec c ≠ b est :
f ( b) - f ( c )
y−f(b) = (x−b) (2.9)
b-c

_______________________________________________________________________________
[Link] Boumediene -33-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

On cherche son intersection avec l’axe OX donc on prend y = 0 et on obtient la


formule donnée plus haut.

· Remarque :
Si « f » est convexe, alors sa représentation graphique est au dessus des
tangentes et en dessous des cordes.
On pose :

b × f (x ) - x × f ( b) b-x
g(x ) = =x- f (x ) (2.10)
f (x ) - f ( b) f ( b) - f ( x )

d’où xn+1=g(xn), montrons que cette suite (xn) est définie.

Pour cela, il suffit de montrer que g([a, b]) Î [a, b]. Vérifions d’abord que g est
de classe C1 sur [a, b]. Elle l’est de manière évidente sur [a, b[. L’application g
est continue en b et g(b) = b − f(b)/f′(b). On a, pour tout x Î [a, b[

b-x f ( b) - f ( x ) - ( b - x ) × f ¢( x )
g¢( x ) = 1 - f ¢( x ) + f (x )
f ( b) - f ( x ) ( f ( b) - f ( x ))2
(2.11)
f ( b) - f ( x ) - ( b - x ) × f ¢( x )
= f ( b)
( f ( b) - f ( x ))2

On en déduit que g′ est continue en b et g′(b) = f(b)·f′′(b)/(2f′(b)2). De plus,


l’application f est convexe, donc f(b) − f(x) − (b − x)f′(x) ≥ 0. On a également f(b)
> 0 car f croissante et f(a)·f(b) < 0. La fonction g est donc croissante sur [a, b].
On a alors g([a, b]) Î [g(a), g(b)]. De plus, g(a) = (a − f(a))·(b − a)/(f(b) − f(a)), or
f(a) < 0 et (b − a)/(f(b) − f(a))≥ 0 par croissance de f ; donc g(a) ≥ a. De même
g(b) = b −f(b)/f′(b), or f′(b) > 0 et f(b) > 0 donc g(b) ≤ b.

La croissance de g montre que la suite (xn) est croissante car x1 = g(a) ≥ a = x0,
or elle est dans l’intervalle [a, b], donc elle converge vers l Î [a, b] tel que, par
continuité de g, g(l) = l. De plus, x0 ≤ α donc une récurrence immédiate et
la croissance de g montrent que, pour tout entier n, xn ≤ α. On en déduit que l ≤
α, c’est-à-dire que l Î ]a, b[.

Soit x Î ]a, b[ tel que f(x) = 0, alors il est immédiat que g(x) = x.
Réciproquement, soit x Î ]a, b[ tel que g(x) = x, alors (b – x)· f(x)/(f(b) − f(x))
= 0, or x ≠b et f(x) ≠ f(b) car f est strictement croissante et x Î ]a, b[. On en
déduit que f(x) = 0. Il y a unicité de la racine de f, donc l = α et la suite (xn)
converge vers α.

2.5. Méthode de Newton (des Tangentes)

On cherche les points fixes de la fonction g(x) = x + u(x)·f(x), et on a vu que :

xn - a
lim = g ¢( a ) = g ¢(a ) (2.12)
n ® +¥ x n +1 - a

_______________________________________________________________________________
[Link] Boumediene -34-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

Pour obtenir une convergence plus rapide, on peut chercher u tel que g′(α) = 0.
On a, si les fonctions ont les régularités nécessaires,
g′(x) = 1 + u′(x)·f(x) + u(x)·f′(x) et on en déduit que g′(α) = 0 si u(α) = −1/f′(α).

Si la fonction « f′ » ne s’annule pas, on peut donc choisir u(x) = −1/f′(x). On


obtient alors la méthode de Newton [12].

2.5.1. Description de la méthode

On considère une fonction réelle définie sur un intervalle I = [a, b] de classe C2


telle que f(a)·f(b) < 0 ; on suppose que les fonctions f′ et f′′ ne s’annulent pas et
gardent chacune un signe constant sur I. On pose

g(x) = x − f(x)/f′(x) (2.13)

Si f′f′′ est positive (respectivement négative) sur [a, b], on pose x0 = b. On


définit alors la suite (xn) par la donnée de x0 et la relation de récurrence
xn+1 = g(xn).

2.5.2. Théorème

La suite (xn) converge vers α l’unique racine de f sur [a, b].

Pour simplifier la rédaction, on va supposer que f′ est strictement positive et f′′


est strictement négative sur I. Les autres cas se traitent de manière similaire.
Ces hypothèses assurent l’existence et l’unicité de α Î I tel que f(α) = 0.

On va montrer que la suite (xn) est croissante et majorée par α. On a x0 = a,


donc x0 ≤ α. Supposons que xn ≤ α, alors, comme xn+1 = xn − f(xn)/f′(xn) et que la
fonction f′ est positive, donc f est croissante, on en déduit immédiatement que
xn+1 ≥ xn et la suite (xn) est croissante.

De plus,
xn+1 − α = g(xn) − g(α) = g′(ξ)(xn − α) (2.14)
avec ξ Î ]xn, α[.

Or
g′(x) = f(x)·f′′(x)/(f′(x))2 (2.15)
donc g′(ξ) > 0 et xn+1 ≤ α.

La suite (xn) est croissante et majorée ; elle est donc convergente. Notons ℓ sa
limite. La continuité de g permet d’écrire que ℓ = g(ℓ) et donc f(ℓ) = 0 d’où ℓ = α.

2.5.3. Interprétation graphique

La méthode de Newton Figure 2.4consiste à remplacer la courbe par sa


tangente en une de ses deux extrémités. Le point x1 est l’intersection de cette
tangente avec l’axe OX [13].

_______________________________________________________________________________
[Link] Boumediene -35-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

Pour faire le dessin, on va se placer dans le cas étudié pour la méthode de la


sécante, f′ > 0 et f′′ < 0. On prend alors x0 = b.

(b, f(b))

f (x1 )
a a x1 b

(a, f(a)) y=f(x)

Figure 2.4: Méthode de Newton.

Traçons la tangente à la courbe représentative de « f » passant par (b, f(b)).


L’équation de cette tangente est y = f′(b)(x− b) + f(b). Son intersection avec
l’axe Ox a une ordonnée nulle et son abscisse vaut x1 = b − f(b)/f′(b). On trace
ensuite la tangente à la courbe au point (x1, f(x1)). Le réel x2 est l’abscisse de
l’intersection de cette deuxième tangente avec l’axe Ox et on réitère ce procédé.

· Remarque :
Si on prenait la tangente à la courbe en (a, f(a)), son intersection avec l’axe Ox
ne serait pas, sur ce dessin, dans l’intervalle [a, b].

2.6. Exemple d’application

Soit la fonction f(x) = x3+4x2-10, avec xÎ [1, 2], ξ=10-3 ;


v Trouvez une valeur approchée de la solution exacte de f(x)=0 par :

1) la méthode de Dichotomie ;
2) la méthode du point fixe ;
3) la méthode de Newton.

· Solutions :

1) Par la méthode de Dichotomie :

f(x) est continue sur [1, 2] et f(1)·f(2) = (-5)·(14)< 0, donc d’après la méthode de
Dichotomie il existe au moins un x*Î [1, 2] tel que f(x*)=0.

On construit une suite x0 = (a+b)/2, x1 = (a1 + b1)/2, . . . , xn = (an + bn)/2 ; et en


vérifier que :

• si f(a)f(x0) < 0, alors α Î ]a, x0[. On pose a1 = a, b1 = x0.

• si f(a)f(x0) = 0, alors α = x0.

• si f(a)f(x0) > 0, alors α Î ]x0, b[. On pose a1 = x0, b1 = b.

_______________________________________________________________________________
[Link] Boumediene -36-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

n an bn xn f(xn)
0 1 2 1.5 +2.375
1 1 1.5 1.25 -1.79687
2 1.25 1.5 1.375 +0.16211
3 1.25 1.375 1.3125 -0.84839
4 1.3125 1.375 1.34375 -0.35098
5 1.34375 1.375 1.35937 -0.09641
6 1.35937 1.375 1.36718 +0.03236
7 1.35937 1.36718 1.36328 -0.03215
8 1.36328 1.36718 1.36523 0.000072

Test d’arrêt :
ξ=10-3 c’est-à-dire |x8 – x7| ≤ ξ ; Alors x* = x8 = 1.36523

2) Par la méthode du point fixe :

L’équation f(x) =0 peut être transformée en une équation équivalente :


g(x) = x, où :
f(x) = x3+4x2-10 = 0 → x2(x+4) = 10
1) existe x*.
1
æ 10 ö 2
2) g(x) = çç ÷÷
è ( x + 4) ø
3) l’application de la méthode du point fixe, formons la suite xn+1=g(xn) ;
on particulier x0 = 1.5, pour ce choix en à :

1er itération x1 = g(x0) = 1.3483


2ème itération x2 = g(x1) = 1.3673
3ème itération x3 = g(x2) = 1.3649
4ème itération x4 = g(x3) = 1.3652
5ème itération x5 = g(x4) = 1.3652

4) Test d’arrêt :
ξ=10-3 c’est-à-dire |x5 – x4| ≤ ξ ; Alors x* = x5 = 1.3652

3) Par la méthode de Newton :

L’équation f(x) = x3+4x2-10 et leur dérivé est f’(x) = 3x2+8, alors on part de x0 =
1.5 et en appliquant la formule de Newton :
xn+1 = xn − f(xn)/f′(xn)
On obtient :
xn+1 = xn – (xn3+4xn2-10) /(3xn2+8)

Nombre d’itérations est de n = 1, 2, 3… jusqu’au le test d’arrêt :


pour n = 1 : x1 = x0 – (x03+4x02-10) /(3x02+8) = 1.37333
pour n = 2 : x2 = x1 – (x13+4x12-10) /(3x12+8)
pour n = 3 : x3 = x2 – (x23+4x22-10) /(3x22+8)

_______________________________________________________________________________
[Link] Boumediene -37-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

2.7. Applications par code MATLAB sur les méthodes numériques de


la résolution de l’équation non linéaire f(x) = 0

2.7.1. Applications pour valider les programmes

Soit : f1(x) = x4 +x–3 [-2, +2]


f2(x) = x3+4x2-10 [1, 2]
f3(x) = 2 x + 1 - x + 4 [1, 4]
f4(x)= cos2x + sin2x + x – 1 [1, 2]

§ En utilisant MATLAB, tracer les fonctions f1, f2, f3 et f4 dans leurs


l’intervalle.
§ Par les méthodes de bissection, Point fixe, Lagrange et Newton, trouver
x telle que f(x)=0 dans l’intervalle [-1.75,-1] et [1, 2].
§ Ecrire les programmes MATLAB qui calculent x dans un intervalle [x1,
x2] par les quatre méthodes.

2.7.2. Application MATLAB de la méthode de Dichotomie

[Link]. Algorithme de la méthode de Dichotomie

Etape 1 Poser a0=a et b0=b


Poser fa= f(a) et fb= f(b)
Etape 2 Poser i=0
Etape 3 Poser xi=(ai+bi) /2 et f=f(xi)
Etape 4 Si xi est une solution (f(xi)=0), aller a l’Etape 10. Sinon, aller a
l’Etape 5.
Etape 5 Si f(ai).f(xi)<0, aller à l’Etape 6
Sinon ( f(ai).f(xi)>0), aller à l’Etape 8
Etape 6 Poser bi+1=xi
Etape 7 Ajouter 1 à i et aller à l’Etape 3
Etape 8 Poser ai+1=xi
Etape 9 Ajouter 1 à i et aller à l’Etape 3
Etape 10 Fin

[Link]. Programme MATLAB de la méthode de Dichotomie

1) Programme MATLAB utilisant la boucle for et la limitation des itérations :

clear all ; close all ; clc


a=1; b=2; % l’intervalle [a, b]
fa=a^3 + 4*a^2 - 10; % l’évaluation par a
fb=b^3 + 4*b^2 - 10; % l’évaluation par b
for i=1:50; % limitation à 50 itérations
x=(a+b)/2; fx=x^3 + 4*x^2 - 10; % l’intervalle [a, b] sur 2

_______________________________________________________________________________
[Link] Boumediene -38-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

if (fa*fx<0) % test de signe


b=x; % choix de l’intervalle [a, x0]
else % changer l’intervalle [x0, b]
a=x;
end
end
x

2) Programme MATLAB utilisant la boucle while et l’erreur epsilon :


clear all ; close all ; clc
a=1; b=2; % l’intervalle [a, b]
fa=a^3 + 4*a^2 - 10; % l’évaluation par a
fb=b^3 + 4*b^2 - 10; % l’évaluation par b
eps=0.001; % l’erreur epsilon
while (b-a) >eps
x=(a+b)/2; fx=x^3 + 4*x^2 - 10;
if (sign(fx) == sign(fa)) % test de signe
a=x; fa=fx;
else
b=x; fb=fx;
end
end
x

3) Programme MATLAB utilisant le fichier fonction:

v Fichier :dichotomie.m

function [x,i] = dichotomie(f,a,b,epsilon,maxit

% ENTREES
% f = fonction objet de type "string"
% [a,b] = intervalle comprenant un (seul) zéro de la fonction
% epsilon = précision souhaitée sur le zéro cherche
% maxit = nombre maximum d'itérations

% SORTIES
% x = suite des approximations du zéro de f
% i = nombre d'itérations
% Test sur les bornes de l'intervalle [a,b]

ya = feval(f.a);
yb = feval(f,b);
ifya*yb > 0.0
error ('Les valeurs de la fonction faux bornes de
T’’intervalle [a,b] sont de même signe.');
end
% Dichotomie
for i = l:maxit % Pour i variant de 1 a maxit.
x(i) = (a+b)/2; % Calcul de la ième approximation.
y = feval(f,x(i)); % Si la précision souhaitée est
% atteinte, stopper.
if (x(i) - a) < epsilon
disp ('Convergence atteinte');

_______________________________________________________________________________
[Link] Boumediene -39-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

break;
end
if y = = 0.0
disp ('Zéro trouvé'); % Si le zéro est trouve, stopper.
break;
elseif y * ya < 0.0
b = x(i); % Si la fonction prend des valeurs
% de signes opposes
yb = feval(f.b); % en x(i) et en a, remplacer la
% borne b par x(i).
else
a = x(i); % Dans le cas contraire, remplacer
% la borne a par x(i).
ya = feval(f,a);
end
if i = = maxit
disp ('Zéro non trouvé'); % Si le nombre max d'itérations
% est atteint, stopper.
break;
end
end

2.7.3. Application MATLAB de la méthode de Point fixe

[Link]. Algorithme de la méthode de Point fixe

Etape 1 Poser a0=a et b0=b ; x0=(a0+b0) /2


Etape 2 Poser i=1
Poser xi =g(xi-1)
Etape 3 Si xi est une solution (f(xi)=0), aller a l’Etape 5. Sinon, aller a
l’Etape 4.
Etape 4 Ajouter 1 à i et aller à l’Etape 2
Etape 5 Fin

[Link]. Programme MATLAB de la méthode de Point fixe

1) Programme MATLAB utilisant la boucle for et la limitation des itérations :


clear all ; close all ; clc
a=1; b=2; % l’intervalle [a, b]
x=(a+b)/2;
for i=1:50 % limitation à 50 itérations
x=sqrt(10/(x+4)); % la forme x = g(x)
end
x

2) Programme MATLAB utilisant le fichier fonction:

v Fichier PointFixe.m
function[X,y,err,k] = PointFixe(g,xl,delta,epsilon,maxit)

_______________________________________________________________________________
[Link] Boumediene -40-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

% Pour déterminer un zéro d'une fonction par la méthode itérative


du point fixe.
% L'équation f(x) = 0 doit être mise sous la forme x = g(x). Le
% processus est convergent si la pente
% de la tangente au voisinage du zéro est inferieure a 1 en valeur
% absolue.

% ENTREES
% g = fonction (objet de type string).
% xl = approximation initiale du zéro.
% delta = tolérance sur l'approximation du zéro.
% epsilon = tolérance sur la valeur de y.
% maxit = nombre maximum d'itérations.

% SORTIES
% X = suite des approximations du zéro.
% y = valeur prise par f pour la dernière valeur approchée
% obtenue.
% err = estimation de l'erreur commise sur le zéro.
% k = nombre d'itérations.

X(l) = xl;
for k = l:maxit
X(k+1) = feval(g,X(k));
err = abs(X(k+l)-X(k));
y - feval(f,X(k+l));
if (err < delta) & (abs(y) < epsilon)
break
end
end

2.7.4. Application MATLAB de la méthode de Sécante (Lagrange)

[Link]. Algorithme de la méthode de Sécante

Etape 1 Poser xi-1=a et xi-2=b ;


Etape 2 Poser i=1
Poser xi= (b-f(b)∙[(b-a)/(f’(b) -f’(a))]
Etape 3 Si xi est une solution (f(xi)=0), aller a l’Etape 5. Sinon, aller a
l’Etape 4.
Etape 4 Ajouter 1 à i et aller à l’Etape 2
Etape 5 Fin

[Link]. Programme MATLAB de la méthode de Sécante

v Programme MATLAB utilisant le fichier fonction (Secante.m):


function [X,err,y,k] = Secante(f,a,b,delta,epsilon,maxit)

% La méthode des sécantes est un procède itératif qui permet


% d'approcher l'unique zéro d'une fonction
% situe dans un intervalle; a partir de deux approximations a et b
% du zéro, elle détermine la suivante

_______________________________________________________________________________
[Link] Boumediene -41-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

% par la relation
% b - a
% x = b - f(b) --------------
% f'(b) - f(a)
%
% Des valeurs de f(a) et f(b) très voisines risquent de provoquer
% une division par 0 au cours du processus.

% ENTREES
% f = fonction (objet de type string).
% a.b = deux premières valeurs approchées du zéro.
% delta = tolérance sur l'approximation du zéro.
% epsilon = tolérance sur la valeur de y.
% maxit = nombre maximum d'itérations.
% SORTIES
% X = suite des approximations du zéro.
% err = estimation de l'erreur commise sur le zéro.
% y = valeur prise par f pour la dernière valeur approchée
% obtenue.
% k = nombre d'itérations.

X(l) = a;
X(2) = b;
yp = feval(f,a);
y = feval (f,b);
for k = 1: maxit
X(k+2) = X(k+l)-y*(X(k +l)-X(k))/(y-yp);
yp - y;
y = feval(f,X(k+2));
en = abs(X(k + l)-X(k));
if (err < delta) & (abs(y) < epsilon)
break
end
end

2.7.5. Application MATLAB de la méthode de Newton

[Link]. Algorithme de la méthode de Newton

Etape 1 Poser a0=a et b0=b ; x0=(a0+b0) /2


Etape 2 Poser i=1
Poser xi = xi-1-f(xi-1)/f’(xi-1)
Etape 3 Si xi est une solution (f(xi)=0), aller a l’Etape 5. Sinon, aller a
l’Etape 4.
Etape 4 Ajouter 1 à i et aller à l’Etape 2
Etape 5 Fin

[Link]. Programme MATLAB de la méthode de Newton

1) Programme MATLAB utilisant la boucle for et la limitation des itérations :


clear all ; close all ; clc
a=1; b=2; % l’intervalle [a, b]

_______________________________________________________________________________
[Link] Boumediene -42-
Chapitre 2 Résolution de l’équation non linéaire f(x) = 0

x=(a+b)/2; % l’intervalle [a, b] sur 2


for i=1:50; % limitation à 50 itérations
x=x-((x^3 + 4*x^2 - 10)/(3*x^2+4*x)); % formule Newton
end
x

3) Programme MATLAB utilisant le fichier fonction:

v Fichier Newton.m
function[X,y,err,k] = Newton(f,df,xl,delta,epsilon,maxit)

% Recherche des tangentes de Newton une valeur approchée d'un zéro


% de la fonction f.

% ENTREES
% f = fonction (objet de type string).
% df = idem, sa dérivée (a calculer soi-même).
% xl = approximation initiale du zéro.
% delta = tolérance sur l'approximation du zéro.
% epsilon = tolérance sur la valeur de y.
% maxit = nombre maximum d'itérations.

% SORTIES
% X = suite des approximations du zéro.
% err = estimation de l'erreur commise sur le zéro.
% y = valeur prise par f pour la dernière valeur approchée
% obtenue.
% k = nombre d'itérations.

X(l) = xl;
for k = 1 :maxit
X(k+1) - X(k) - feval(f,X(k))/feval(df,X(k));
err = abs(X(k+l)-X(k));
y = feval(f,X(k + l));
if (err < delta) & (abs(y) < epsilon)
break
end
end

_______________________________________________________________________________
[Link] Boumediene -43-
Chapitre 3

Résolution d’un Système


d’Équations Linéaires

[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

3.1. Introduction

Les systèmes d'équations algébriques jouent un rôle très important en


ingénierie. On peut classer ces systèmes en deux grandes familles: les
systèmes linéaires et les systèmes non linéaires (voir chapitre 7).
Ici encore, les progrès de l'informatique et de l'analyse numérique permettent
d'aborder des problèmes de taille prodigieuse. On résout couramment
aujourd'hui des systèmes de plusieurs centaines de milliers d'inconnues. On
rencontre ces applications en électrotechnique (écoulement des puissances),
mécanique des fluides et dans l'analyse de structures complexes. On peut par
exemple calculer l'écoulement des puissances dans un réseau électrique [14].
Ces calculs complexes requièrent des méthodes sophistiquées comme les
méthodes directes de Gauss et les méthodes itératives Gauss-Seidel. On
obtient généralement des systèmes d'équations non linéaires de taille
considérable, qu'on doit résoudre à l'aide de méthodes efficaces qui minimisent
le temps de calcul et l'espace-mémoire requis [15, 16].
Dans ce chapitre, nous allons aborder les principales méthodes de résolution
des systemes linéaires, à savoir la méthode d'élimination de Gauss et la
méthode de Gauss-Seidel.

3.2. Généralités

éa11 a12 a13 a14 ù


êa a 22 a 23 a 24 úú
ü Soit une Matrice carrée A : ê
21
, chaque élément de A
êa 31 a 32 a 33 a 34 ú
ê ú
ëa 41 a 42 a 43 a 44 û
est a ij avec i c’est l’indice de ligne et j l’indice de colonne.

éa11 0 0 0 ù
ê0 a 22 0 0 úú
ü Si la matrice A = ê → dite matrice Diagonale.
ê0 0 a 33 0 ú
ê ú
ë0 0 0 a 44 û
ü Si les éléments diagonaux de A sont égaux à 1 → dite matrice Unité.
ü « A » matrice est dite Régulière si det(A)≠0, sinon est dite Singulière.
éa11 a12 a13 a14 ù
ê0 a a 23 a 24 úú
ü « A » matrice est dite Triangulaire Supérieur si ê 22
.
ê0 0 a 33 a 34 ú
ê ú
ë0 0 0 a 44 û
éa11 0 0 0 ù
êa a 22 0 0 úú
ü « A » matrice est dite Triangulaire Inférieur si ê 21
.
êa 31 a 32 a 33 0 ú
ê ú
ëa 41 a 42 a 43 a 44 û

_______________________________________________________________________________
-45-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

3.3. Problème

Soit le système linéaire sous la forme suivant :

ìa11x1 + a12 x 2 + a13 x 3 + × × × + a1n x n = b1


ï
ïa 21x1 + a 22x 2 + a 23x 3 + × × × + a 2n x n = b2
ï
ía31x1 + a 32x 2 + a33 x 3 + × × × + a3n x n = b3 (3.1)
ï
ï
ïîa n1 x1 + an 2 x 2 + a n3 x 3 + × × × + ann x n = bn

L’équation (3.1) est équivalente au système matriciel A•X=b suivant :

éa11 a12 a13 L a1n ù é x1 ù é b1 ù


êa a22 a23 L a 2n úú êêx 2 úú êêb2 úú
ê 21
êa31 a32 a33 L a3n ú × êx3 ú = ê b3 ú (3.2)
ê ú ê ú ê ú
ê M M M O M ú êMú êMú
êëan1 a n2 an3 L ann úû êëxn úû êëbn úû

3.4. Méthode Directe d’Elimination de Gauss

Tous les outils sont en place pour la résolution d'un système linéaire. H suffit
maintenant d'utiliser systématiquement les opérations élémentaires pour
introduire des zéros sous la diagonale de la matrice A et obtenir ainsi un
système triangulaire supérieur [17, 18].
La validité de la méthode d'élimination de Gauss repose sur le fait que les
opérations élémentaires consistent à multiplier le système de départ par une
matrice inversible. La méthode d'élimination de Gauss consiste à éliminer tous
les termes sous la diagonale de la matrice A.

Soit le Système linéaire de quatre équations suivant :

ìa11 x1 + a12 x 2 + a13 x 3 + a14 x 4 = b1 (1)


ï
ïa 21 x1 + a 22 x 2 + a 23 x 3 + a 24 x 4 = b 2 ( 2)
í (3.3)
ïa 31 x1 + a 32 x 2 + a 33 x 3 + a 34 x 4 = b3 (3 )
ïîa 41 x1 + a 42 x 2 + a 43 x 3 + a 44 x 4 = b 4 (4)

v Équation (1) → a 11 ≠ 0 :

Divisons l’équation (1) du système (3.3) par a 11 en obtient :


a a a b
x 1 + 12 x 2 + 13 x 3 + 14 x 4 = 1 (3.4)
a 11 a 11 a 11 a 11

(1) (1) (1) (1)


x 1 + a 12 x 2 + a 13 x 3 + a 14 x 4 = b 1 (3.5)

_______________________________________________________________________________
-46-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

(1) a1 j (1)
b1
Avec : a 1 j = ; b1 = .
a 11 a 11

Essayant d’éliminer la variable x1 de l’équation (2) du système (3.3) pour ce là


on va multiplier l’équation (3.5) par a 21 et la soustrais de l’équation (2) : {éq(2)-
a 21 .éq(3.5)}.
De la même façon pour éliminer x 1 des équations (3) et (4) : {éq(3)- a 31 .éq(3.5)}
et {éq(4)- a 41 .éq(3.5)}, on aura :

ì (1) (1 ) (1 ) (1)
x
ï 1 + a 12 2x + a 13 x 3 + a 14 x 4 = b1 (1)
ï (1) (1) (1 ) (1)
0
ï 1x + a 22 x 2 + a 23 3x + a 24 x 4 = b 2 (2)
í (1) (1) (1) (1 ) (3.6)
ï0x + a x + a x + a x = b (3)
ï 1 32 2 33 3 34 4 3
(1) (1) (1 ) (1 )
ï0x + a x + a x + a x = b (4)
î 1 42 2 43 3 44 4 4

é (1 ) (1 ) (1 )
ù é (1) ù
ê1 a12 a13 a14 ú é x1 ù ê b1 ú
ê0 (1 ) (1 ) (1 )
ê ú (1)
a 22 a 23 a 24 ú êx 2 ú ê b 2 ú
ê ú× =ê ú (3.7)
ê0
(1 ) (1 ) (1 )
ú êx 3 ú ê (1) ú
a 32 a 33 a 34 ê ú b
ê (1 ) (1 ) (1 ) ú x ê (13) ú
ê0 ë 4û ê ú
ë a 42 a 43 a 44 úû ëb 4 û
(1) (1) ì i = 2,3,4 ü
Alors : a ij = a ij - a i1 × a 1 j í ý
î j = 1,2,3,4þ
(1)
v Équation (2) → a 22 ≠ 0 :
(1)
Divisons l’équation (2) du système (3.6) par a 22 en obtient :
( 2) ( 2) ( 2)
x 2 + a 23 x 3 + a 24 x 4 = b 2 (3.8)
(1) (1)
(2) a 2j ( 2)
b2
Avec : a 2 j = (1)
; b2 = (1)
.
a 22 a 22

Essayant d’éliminer la variable x 2 des équations (3) et (4) du système (3.6)


(1)
pour ce là on va multiplier l’équation (3.8) par a 32 et la soustrais de l’équation
(1)
(3) du système (3.6) : {éq(3)- a 32 .éq(3.8)} et de la même façon pour éliminer x 2
(1)
de l’équations (4) : {éq(4)- a 42 .éq(3.8)} on aura :

_______________________________________________________________________________
-47-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

ì (1) (1 ) (1 ) (1 )
x
ï 1 + a 12 2x + a x
13 3 + a 14 x 4 = b1 (1)
ï ( 2) (2) ( 2)
0x
ï 1 + x 2 + a 23 3x + a 24 x 4 = b 2 ( 2)
í ( 2) ( 2) ( 2) (3.9)
ï0x + 0x + a x + a x = b (3 )
ï 1 2 33 3 34 4 3
(2 ) ( 2) ( 2)
ï0x + 0x + a x + a x = b (4)
î 1 2 43 3 44 4 4

é (1) (1 ) (1)
ù é (1 ) ù
ê 1 a 12 a13 a14 ú é x1 ù ê b1 ú
ê0 1 (2) ( 2)
ê ú (2)
a 23 a 24 ú êx 2 ú ê b 2 ú
ê ú× =ê ú (3.10)
ê0 0
(2) ( 2)
ú êx 3 ú ê (2) ú
a33 a 34 ê ú b
ê (2) ( 2) ú x ê ( 23) ú
êë0 0 a 44 úû ë û êë b4 úû
4
a 43

( 2) (1) (1) ( 2) ì i = 3,4 ü


Alors : a ij = a ij - a i 2 × a 2 j í ý
î j = 2,3,4þ
( 2)
v Équation (3) → a 33 ≠ 0 :
( 2)
Divisons l’équation (3) du système (3.9) par a 33 en obtient :

( 3) (3)
x 3 + a 34 x 4 = b 3 (3.11)
( 2) (2)
( 3) a 3j ( 3)
b3
Avec : a 3 j = ( 2)
; b3 = (2)
.
a 33 a 33

Essayant d’éliminer la variable x 3 de l’équation (3) du système (3.9) pour ce là


( 2)
on va multiplier l’équation (3.11) par a 43 et la soustrais de l’équation (2) :
( 2)
{éq(4)- a 43 .éq(3.11)}.

ì (1) (1 ) (1) (1 )
x
ï 1 + a 12 x 2 + a 13 x 3 + a 14 x 4 = b1
ï ( 2) (2) ( 2)
ï0x 1 + x 2 + a 23 x 3 + a 24 x 4 = b 2
í (3) (3) (3.12)
ï0x + 0x + x + a x = b
ï 1 2 3 34 4 3

ï0x + 0x + 0x + a x = b( 3 ) (3)

î 1 2 3 44 4 4

_______________________________________________________________________________
-48-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

é (1 ) (1 ) (1 )
ù é (1 ) ù
ê 1 a 12 a13 a14 ú é x1 ù ê b1 ú
ê0 1 ( 2) ( 2)
ê ú (2 )
a 23 a 24 ú êx 2 ú ê b 2 ú
ê ú× =ê ú (3.13)
ê0 0
(3)
ú êx 3 ú ê (3) ú
1 a 34 ê ú b
ê (3) ú x ê ( 33) ú
êë0 0 ë 4 û êb ú
0 a 44 úû ë 4û
( 3) ( 2) ( 2) (3) ì i=4 ü
Alors : a ij = a ij - a i 3 × a 3 j í ý
î j = 3,4þ

Donc la méthode de Gauss transforme un système de matrice A plaine à un


système à matrice A triangulaire supérieur par élimination successive des
variables.

· Remarque :
(1) ( 2)
Les éléments a 11 a 22 a 33 , sont appelés «Pivots» ou «éléments générateurs».

3.5. Exemples sur la méthode d’élimination de Gauss

v Utilisant l’élimination de Gauss résoudre les systèmes d’équations


linéaires suivants:

ì2 x1 + x2 + 4 x4 = 2
ï - 4 x - 2 x + 3 x - 7 x = -9
ï
í
1 2 3 4
(a)
ï4 x1 + x2 - 2 x3 + 8x4 = 2
ïî- 3x2 - 12 x3 - x4 = 2

é 6 -4 1 ù é- 14 22 ù
ê ú
A = ê- 4 6 - 4ú B = ê 36 - 18ú (b)
ê ú
êë 1 - 4 6 úû êë 6 7 úû

· Solutions :

1) Système d’équations linéaires (a):

ì2x 1 + x 2 + 4x 4 = 2
ï
ï- 4 x1 - 2x 2 + 3x 3 - 7x 4 = -9
Soit à résoudre le système linéaire : í
ï4 x1 + x 2 - 2x 3 + 8x 4 = 2
ï- 3x 2 - 12x 3 - x 4 = 2
î

é 2 1 0 4ù éx 1 ù é2 ù
ê- 4 - 2 3 - 7úú êx ú ê- 9 ú
1) En posant : A = ê , X = ê 2 ú et B = ê ú
ê4 1 -2 8 ú êx 3 ú ê2 ú
ê ú ê ú ê ú
ë 0 - 3 - 12 - 1 û ëx 4 û ë2 û

_______________________________________________________________________________
-49-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

Le système linéaire s’écrit sous la forme matricielle suivante A∙X = B :

é2 1 0 4 ù é x1 ù é 2 ù
ê- 4 - 2 3
ê - 7úú êêx 2 úú êê- 9úú
× =
ê4 1 - 2 8 ú êx 3 ú ê 2 ú
ê ú ê ú ê ú
ë 0 - 3 - 12 - 1 û ëx 4 û ë 2 û

2) Résolution du système par la méthode de Gauss : Partons de la matrice


augmentée :

é ù
ê2 1 0 4 2 ú
ê ú
ê- 4 - 2 3 -7 - 9ú
ê4 1 -2 8 2 ú
ê ú
ê1 0 - 3 - 12 - 1 2

44424443
êë A B ú
û

L’élément générateur a11 = 2 ≠ 0. Remplaçons la 2ème ligne par :


ligne 2 − (− 2)× la ligne 1, en trouve :

é2 1 0 4 2ù
ê0 0 ú
ê 3 1 - 5ú
ê4 1 -2 8 2ú
ê ú
ë0 - 3 - 12 - 1 2 û

Remplaçons la troisième ligne par : ligne 3 − 2 × ligne 1, en aperçoit :

é2 1 0 4 2 ù
ê ú
ê0 0 3 1 - 5ú
ê0 - 1 - 2 0 - 2ú
ê ú
ë0 - 3 - 12 - 1 2 û

Et puisque l’élément générateur a22 est nul, faisons une permutation entre la
deuxième ligne et la troisième ligne :

é2 1 0 4 2 ù
ê ú
ê0 - 1 - 2 0 - 2ú
ê0 0 3 1 - 5ú
ê ú
ë0 - 3 - 12 - 1 2 û

On remplace la quatrième ligne par : ligne 4 − 3 × ligne 2, en trouve :

_______________________________________________________________________________
-50-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

é2 1 0 4 2 ù
ê ú
ê0 - 1 - 2 0 - 2ú
ê0 0 3 1 - 5ú
ê ú
ë0 0 - 6 - 1 8 û

À l’étape suivante, nous remplaçons la 4ème ligne par ligne 4 − (− 2) × ligne 3 :

é2 1 0 4 2 ù
ê ú
ê0 - 1 - 2 0 - 2ú
ê0 0 3 1 - 5ú
ê ú
ë0 0 0 1 - 2û

Nous avons abouti au système triangulaire supérieur suivant :

é2 1 0 4ù é x1 ù é 2 ù
ê0 - 1 - 2
ê 0 úú êêx 2 úú êê- 2úú
× =
ê0 0 3 1 ú ê x 3 ú ê - 5ú
ê ú ê ú ê ú
ë0 0 0 1 û ëx 4 û ë- 2û

D’où on tire les équations :

ì x 4 = -2
ï 3x + x = -5
ï 3 4
í
ï- x 2 + 2x 3 = -2
ïî2x 1 + x 2 + 4x 4 = 2

La solution du système est :

x1 = 3 x2 = 4 x 3 = -1 x 4 = -2

2) Système d’équations linéaires (b):

Soit à résoudre le système linéaire :

é 6 -4 1 ù é - 14 22 ù
A = ê - 4 6 - 4ú , B = êê 36 - 18úú
ê ú
êë 1 - 4 6 úû êë 6 7 úû

La matrice augmentée est :

é 6 - 4 1 - 14 22 ù
ê ú
ê- 4 6 - 4 36 - 18ú
ê 1 -4 6 6 7 úû
ë

_______________________________________________________________________________
-51-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

Elimination 1 :

{ligne(2)-(-2/3).ligne(1)} → Nouvelle ligne(2)


{ligne(3)-(1/6).ligne(1)} → Nouvelle ligne(3)

é6 -4 1 - 14 22 ù
ê ú
ê0 10 / 3 - 10 / 3 80 / 3 - 10 / 3ú
ê0 - 10 / 3 35 / 6 25 / 3 10 / 3 ú
ë û

Elimination 2 :

{ligne(3)-(-1).ligne(2)} → Nouvelle ligne(3)

é6 - 4 1 - 14 22 ù
ê ú
ê0 10 / 3 - 10 / 3 80 / 3 - 10 / 3ú
ê0 0 5/ 2 35 0 úû
ë

La solution du système est :

x11 = 10 x 21 = 22 x 31 = 14
Et
x12 = 3 x 22 = -1 x 32 = 0

3.6. Méthode Itérative de Gauss-Seidel

La méthode de Gauss-Seidel est fondée sur la simple constatation selon


( k +1 ) (k) (k ) (k)
laquelle le calcul de x 2 nécessite l'utilisation de x1 , x 3 , … , x n provenant de
( k +1 )
l'itération précédente. Or, a l'itération k + 1, au moment du calcul de x 2 , on
(k) ( k +1)
possède déjà une meilleure approximation de x1 que x i , a savoir x1 .
( k +1) ( k +1) ( k +1 )
De même, au moment du calcul de x 3 , on peut utiliser x1 et x 2 . Plus
( k +1) ( k +1 ) ( k +1 ) ( k +1 )
généralement, pour le calcul de x i , on peut utiliser x1 , x 2 , … , x i -1 déjà
(k ) (k) (k)
calcules et les x i +1 , x i +2 , … , x n l'itération précédente [18, 19].

La méthode de Gauss-Seidel transforme un système (3.3) à un système (3.14)


comme suit :

_______________________________________________________________________________
-52-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

ì( k +1) ( k) ( k) (k )

ï x 1 = ( b1 - ( a12 x 2 + a 13 x 3 + a 14 x 4 ))
ï a11
ï( k +1) ( k +1 ) (k) (k )
ï x 2 = ( b2 - (a 21 x1 + a 23 x 3 + a 24 x 4 ))
ï a 22
í( k +1) ( k +1) ( k +1 ) (k ) (3.14)
ï x = ( b - (a x + a x + a x ))
ï 3 3 31 1 32 2 34 4
a 33
ï
ï( k +1) ( k +1 ) ( k +1 ) ( k +1)
x
ï 4 = ( b 4 - ( a 41 x 1 + a 42 x 2 + a 43 x 3 ))
î a 44

Le principe de la méthode de Gauss-Seidel consistant à utiliser un vecteur


0 (0 ) (0) ( 0) (0)
Estimé (Initial) x = ( x 1 , x 2 , x 3 , x 4 ) de la solution exacte et on essaye d’améliorer
cette solution par la suite :
( k +1) (k ) (k) ( k) (k ) ( k +1) ( k)
x = F( x 1 , x 2 , x 3 , x 4 ) , avec x = [ x 1 , x 2 , x 3 , x 4 ]T jusqu’à où | x - x |<ε.

3.7. Exemple sur la méthode de Gauss-Seidel


v Par la méthode de Gauss-Seidel résoudre les systèmes d’équations
suivants:
ì 4x 1 - x 2 + x 3 = 7
ï
í- 4 x1 + 8x 2 - x 3 = 21
ï - 2x + x + 5x = 15
î 1 2 3

· Solution :

ì 4x 1 - x 2 + x 3 = 7
ï
Soit à résoudre le système linéaire : í- 4x 1 + 8x 2 - x 3 = 21
ï - 2x + x + 5x = 15
î 1 2 3

é1 ù
0 ê ú
Avec la solution estimée : x = ê2ú et « ε » erreur 10-4.
êë2úû

ì x 1 = (7 - ( -x 2 + x 3 )) / 4
ï
íx 2 = (21 - ( -4x 1 - x 3 )) / 8
ï x = (15 - ( -2x + x )) / 5
î 3 1 2

é1 ù é1.7500 ù é1.9500 ù é1.9956 ù é1.9993 ù


ê ú ê ú ê ú ê ú
x 0 = ê2ú → x1 = ê3.7500 ú → x 2 = ê3.9688 ú → x 3 = ê3.9961 ú → x 4 = êê3.9995 úú
êë2úû êë2.9500 úû êë2.9863úû êë2.9990 úû êë2.9998 úû

_______________________________________________________________________________
-53-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

é1.9999 ù é2.0000 ù é2.0000 ù


5 ê ú 6 ê ú 7 ê ú
→ x = ê3.9999 ú → x = ê4.0000 ú → x = ê4.0000 ú
êë3.0000 úû êë3.0000 úû êë3.0000 úû

3.8. Applications par code MATLAB sur les méthodes numériques de


la résolution d’un Système d’Équations Linéaires

3.8.1. Applications pour valider les programmes

1) Former la matrice augmentée


é2 1 2 10 ù
ê ú
ê6 4 0 26ú
ê8 5 1 35ú
ë û

2) Former la matrice Triangulaire Supérieur utilisant d’Elimination de Gauss :

é2 1 2 10 ù é2 1 2 10 ù é2 1 2 10 ù
ê ú ê ú ê ú
ê6 4 0 26ú → ê0 1 - 6 - 4 ú → ê0 1 - 6 - 4ú
ê8 5 1 35ú ê0 1 - 7 - 5 ú ê0 0 - 1 - 1 úû
ë û ë û ë

3) Calculer les éléments du vecteur X :

10 - (1)(2) - (2)(1) - 4 - ( -6)(1) -1


x1 = = 3 ; x2 = = 2 ; x3 = = 1.
2 1 -1

3.8.2. Application MATLAB de la méthode d’élimination de Gauss

[Link]. Algorithme de la méthode d’élimination de Gauss

Etape 1 Donner : A•X=B → A(n×n), B(n×1) ;


Etape 2 Former la matrice augmentée AB(n×n+1) : AB = [A|B] ;
Etape 3 Former la matrice Triangulaire Supérieur utilisant d’Elimination
de Gauss :

k= 1, ..., n-1
i= k+1, ..., n-1
j=1, ..., n

AB ij = AB ij - AB ik × AB kj
AB kk

Etape 4 Dissocier la matrice Triangulaire Supérieur finale AB = [a|b] à :


a(n×n) et b(n×1)) ;

_______________________________________________________________________________
-54-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

Etape 5 Calculer les éléments du vecteur X :


k = n-1, ..., 1
j = k+1, ..., n-1

b k = b k - a kj × b j
a kk
Etape 6 Afficher X.

[Link]. Programme MATLAB de la méthode d’élimination de Gauss

1) Programme MATLAB généralisé :


clear all;close all;clc
A=[2 1 2;6 4 0;8 5 1] % La matrice A
B=[10;26;35] % Le vecteur B

AB=[A B] % La matrice augmentée AB


n = length(AB); % La taille de AB
for k=1:n-1
for i=k+1:n-1
j=1:n;
if AB(k,k)~=0; % Former la matrice Triangulaire
% Supérieur utilisant d’Elimination
% de Gauss
AB(i,j) = AB(i,j) - AB(i,k)*AB(k,j)/AB(k,k);
else
AB([k k+1],:) = AB([k+1 k],:); % Pivotation des
% lignes
end
end
end
AB

% Dissocier la matrice Triangulaire Supérieur finale AB :


a=AB(:,1:n-1);
b=AB(:,n);

% Calculer les éléments du vecteur X

for k = n-1:-1:1
j = k+1:n-1
b(k) = (b(k) - a(k,j)*b(j))/a(k,k);
end
x = b;
x %Affichage de X

_______________________________________________________________________________
-55-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

2) Programme MATLAB utilisant le fichier fonction:

v Fichier :eliminGauss.m

function [x] = eliminGauss (A, B)

AB=[A B] % La matrice augmentée AB


n = length(AB); % La taille de AB
for k=1:n-1
for i=k+1:n-1
j=1:n;
if AB(k,k)~=0; % Former la matrice Triangulaire
% Supérieur utilisant d’Elimination
% de Gauss
AB(i,j) = AB(i,j) - AB(i,k)*AB(k,j)/AB(k,k);
else
AB([k k+1],:) = AB([k+1 k],:); % Pivotation des
% lignes
end
end
end
AB
% Dissocier la matrice Triangulaire Supérieur finale AB
a=AB(:,1:n-1);
b=AB(:,n);

% Calculer les éléments du vecteur X


for k = n-1:-1:1
j = k+1:n-1
b(k) = (b(k) - a(k,j)*b(j))/a(k,k);
end
x = b;
x %Affichage de X

3.8.3. Application MATLAB de la méthode de Gauss-Seidel

[Link]. Algorithme de la méthode de Gauss-Seidel

Etape 1 Donner : A•X=B → A(n×n), B(n×1) ;


X0 = [X10, X20, X30, ..., Xn0]T ; ε.
Etape 2 Calculer les éléments du vecteur X :
Tant que (| x new
i - x old
i |<ε)
i -1 n
b i - å a ij × x new
j - åa ij × x old
i
j=1 j=1+1
x new
i =
a ii
Fin

Etape 3 Afficher X.

_______________________________________________________________________________
-56-
[Link] Boumediene
Chapitre 3 Résolution d’un Système d’Équations Linéaires

[Link]. Programme MATLAB de la méthode de Gauss-Seidel

clear all;close all;clc


A=[4 -1 1; -4 8 -1; -2 1 5] % La matrice A
B=[ 7; 21; 15] % Le vecteur B
x=[0;0;0] % La solution estimée ou vecteur initial x

n=length(x)
for k=1:100 % Nbr des itérations
for i=1:n
% calcul des éléments supérieurs
if i==1
x(i)=(B(i)-(A(i,2:n)*x(2:n)))/A(i,i)
% calcul des éléments inferieurs
elseif i==n
x(i)=(B(i)-(A(i,1:n-1)*x(1:n-1)))/A(i,i)
else
% calcul des éléments diagonaux
x(i)=(B(i)-(A(i,1:i-1)*x(1:i-...
1)+(A(i,i+1:n)*x(i+1:n))))/A(i,i)
end
end
end

_______________________________________________________________________________
[Link] Boumediene -57-
Chapitre 4

Interpolation polynômiale

[Link] Boumediene
Chapitre 4 Interpolation polynômiale

4.1. Introduction

Il s'agit d'un problème d'interpolation, dont la solution est relativement simple.


Il suffit de construire un polynôme de degré suffisamment élevé dont la courbe
passe par les points de collocation. On parle alors du polynôme de collocation
ou polynôme d'interpolation. Pour obtenir une approximation des dérivées ou
de l'intégrale, il suffit de dériver ou d'intégrer le polynôme de collocation. El y
a cependant des éléments fondamentaux qu'il est important d'étudier. En
premier lieu, il convient de rappeler certains résultats cruciaux relatifs aux
polynômes, que nous ne démontrons pas [20].

Les fonctions les plus faciles à évaluer numériquement sont les fonctions
polynômes. Il est donc important de savoir approximer une fonction arbitraire
par des polynômes.

Étant donnés n + 1 couples de réels (xi , yi ), il s’agit de trouver une fonction f


telle que f(xi) = yi pour i = 0, . . . , n. f peut être un polynôme (par
interpolation), un polynôme trigonométrique ou une fonction régulière
polynômiale par morceaux (spline). On désigne par l’espace des polynômes de
degré inférieur ou égal à k ainsi que l’espace des fonctions polynômiales
correspondantes [20, 21].

4.2. Généralités

Le problème de l’interpolation consiste à chercher des fonctions (polynômes)


passant par des points donnés : (x 0 , y 0 ) , (x 1 , y1 ) , ..., (x n , y n ) c.-à-d., on
cherche p( x ) avec p( x i ) = y i pour i = 1,2,....n . Si les valeurs de y i satisfont
y i = f (x i ) où f ( x ) est une fonction donnée, il est intéressant d’étudier l’erreur
de l’approximation : f ( x ) - p(x ) = erreur (Problème de l’interpolation) [15].

Un polynôme de degré n dont la forme générale est:

pn ( x ) = a 0 + a1 x + a 2 x 2 + a 3 x 3 + ... + a n x n (4.1)

Possède très exactement n racines qui peuvent être réelles ou complexes


conjuguées.
Les données → (x 0 , y 0 = f (x 0 )),( x 1 , y1 = f ( x1 )),...,( x i , y i = f ( x i )) .
La solution recherchée → p n ( x ) tel que p n ( x i ) = f ( x i ); i = 0,1,..., n .

4.3. Matrice de Vandermonde

Le problème d'interpolation consiste donc a déterminer l'unique polynôme de


degré n passant par les (n+1) points de collocation ( x i , y i = f ( x i ) ) pour i =
0,1,…, n). Il reste maintenant a le construire de la manière la plus efficace et
la plus générale possible [15, 22]. Une première tentative consiste à
déterminer les inconnues a i du polynôme (4.1) en vérifiant directement les
(n+1) équations de collocation :

_______________________________________________________________________________
[Link] Boumediene -59-
Chapitre 4 Interpolation polynômiale

p n (x i ) = f ( x i ), i = 0,1,..., n (4.2)

Ou encore :
n
f ( x i ) = a 0 + a1 x i + a 2 x 2i + a 3 x 3i + ... + a n x ni = åa x
i =0
i
i
; i = 0,1,..., n (4.3)

qui est un système linéaire de (n+1) équations en (n+1) inconnues. Ce système


s'écrit sous forme matricielle:

é ù
ê1 x 0 x 20 x 30 L x n0 ú éa 0 ù éf (x 0 )ù
ê nú ê ú ê f (x ) ú
ê1 x1 x1 x1 L x1 ú ê a1 ú
2 3
ê 1 ú
ê1 x x 2 x 3 L x n ú × êa 2 ú = êf ( x 2 ) ú (4.4)
ê 2 2 2 2
ú ê ú ê ú
êM M M M O M ú êM ú ê M ú
ê1 x x 2n x 3n L x nn ú êëa n úû êëf ( x n )úû
ê 14n444 244443 ú
ë Matrice de Vandermonde û

La matrice de ce système linéaire porte le nom de matrice de Vandermonde.


On peut montrer que le conditionnement de cette matrice augmente fortement
avec la taille (n+1) du système. De plus, comme le révèlent les sections qui
suivent, il n'est pas nécessaire de résoudre un système linéaire pour calculer
un polynôme d'interpolation. Cette méthode est donc rarement utilisée.

4.4. Interpolation de Lagrange

L'interpolation de Lagrange est une façon simple et systématique de construire


un polynôme de collocation [15].

Etant donne (n+1) points ( x i , y i = f ( x i ) ) pour i = 0,1,…, n), on suppose un


instant que l'on sait construire (n+1) polynômes L i ( x ) de degré n et
satisfaisant les conditions suivantes :

L i (x i ) = 1 "i
(4.5)
L i (x j ) = 0 "j¹ i

Cela signifie que le polynôme L i ( x ) de degré n prend la valeur 1 en x i et


s'annule a tous les autres points de collocation. Nous verrons comment
construire les L i ( x ) un peu plus loin. Dans ces conditions, la fonction
L(x ) définie par :

n
L(x ) = å f ( x i ) × L i (x ) (4.6)
i =0

est un polynôme de degré n, car chacun des L i ( x ) est de degré n. De plus, ce


polynôme passe par les (n+1) points de collocation et est donc le polynôme
recherche. En effet, il est facile de montrer que selon les conditions (4.5):

_______________________________________________________________________________
[Link] Boumediene -60-
Chapitre 4 Interpolation polynômiale

n
L(x j ) = f ( x j ) × L j ( x j ) + å f ( x i ) × L i (x j ) (4.7)
i =0
i¹j

Le polynôme L( x ) passe donc par tous les points de collocation. Puisque ce


polynôme est unique, L( x ) est bien le polynôme recherche. Il reste à construire
les fonctions L i ( x ) . Suivons une démarche progressive.

Soient x 0 , x1 ,..., x n : (n+1) points distincts et « f » une fonction dont les valeurs
en ces points sont : f ( x 0 ), f ( x 1 ),..., f ( x n ) . Alors il existe un seul polynôme de
degré inferieur ou égal à n tel que : f ( x i ) = p n ( x i ), i = 0,1,..., n .
Ce polynôme est donné par :

Pn ( x ) = f ( x 0 ) × L 0 ( x ) + f ( x 1 ) × L1 ( x ) + ... + f ( x n ) × L n (x ) (4.8)
n
Pn ( x ) = å f ( x i ) × L i ( x ), i = 0,1,..., n (4.9)
i =0

Où :
[(x - x 0 )...(x - x i -1 )][(x - x i +1 )...(x - x n )]
L i (x) =
[(x i - x 0 )...(x i - x i -1 )][(x i - x i +1 )...(x i - x n )]
n (x - x j ) (4.10)

j= 0
j¹ i
(x i - x j )

4.5. Polynôme de Newton

Lorsqu’on écrit l’expression générale d’un polynôme, on pense immédiatement


à la forme (4.1), qui est la plus utilisée [15, 23]. Il en existe cependant d'autres
qui sont plus appropriées au cas de l’interpolation, par exemple:

pn ( x ) = a 0
+ a1 (x - x 0 )
+ a 2 ( x - x 0 )(x - x 1 )
(4.11)
+ a 3 (x - x 0 )(x - x1 )(x - x 2 )
M
+ a n ( x - x 0 )(x - x 1 )( x - x 2 ) L (x - x n -1 )

Calcul des a n : méthode des différences divisées.

Les n éme différences divisées de la fonction f ( x ) sont définies à partir de la


façon suivante:

f [ x1 , x 2 ,L , x n ] - f [ x 0 , x 1 , x 2 ,L , x n -1 ]
f [ x 0 , x 1 , x 2 ,L , x n ] = (4.12)
(x n - x 0 )

_______________________________________________________________________________
[Link] Boumediene -61-
Chapitre 4 Interpolation polynômiale

Tableau 4.1 : Calcul des différences divisées.


xi f (x i ) f [ x i , x i +1 ] f [ x i , x i +1 , x i +2 ] f [ x i , x i +1 , x i +2 , x i +3 ]
x0 f (x 0 )
f [x 0 , x1 ]
x1 f ( x1 ) f [ x 0 , x1 , x 3 ]
f [x1 , x 2 ]
f [ x 0 , x1 , x 2 , x 3 ]
x2 f (x 2 )
f [ x1 , x 2 , x 3 ]
f [x 2 , x 3 ]
x3 f (x 3 )

4.6. Exemple d’application

v On doit calculer le polynôme passant par les 4points (0, 1), (1, 2), (2, 9)
et (3, 28) utilisant :

1) La Matrice de Vandermonde.
2) L’interpolation de Lagrange.
3) Le Polynôme de Newton.

· Solutions :

1) Par la Matrice de Vandermonde:

On doit calculer le polynôme figure 4.1 passant par les points (0, 1), (1, 2), (2,
9) et (3, 28). Etant donne ces 4 points, les coefficients a i sont les solutions de :
é ù
ê1 x 0 x 20 x 30 L x n0 ú éa 0 ù é f ( x 0 )ù é ù
ê nú ê ú ê ú ê 1 0 0 0 ú éa 0 ù é 1 ù
ê1 x1 x1 x1 L x1 ú ê a1 ú ê f ( x1 ) ú
2 3
ê ú ê ú ê ú
ê1 x x 2 x 3 L x n ú × êa 2 ú = ê f ( x 2 )ú º ê 1 1 1 1 ú × ê a1 ú = ê 2 ú
ê 2 2 2 2
ú ê ú ê ú ê 1 2 4 8 ú êa 2 ú ê 9 ú
êM M M M O M ú êM ú ê M ú ê ú ê ú ê ú
ê1 x ê 114 3 9 27 ú ëa 4 û ë28û
x 2n x 3n L x nn ú êëa n úû êëf ( x n )úû 4244 3
ê 14n444 244443 ú ëê Matrice de Vandermondeûú
ë Matrice de Vandermonde û

dont la solution (par Elimination de Gauss (chapitre 3)) est a = [1 0 0 1]T.


Le polynôme recherche est donc : P3 ( x ) = 1 + x 3 .

2) Par l’interpolation de Lagrange:

§ Le polynôme passant par les points (0, 1), (1, 2), (2, 9) et (3, 28) est :

( x - 1)(x - 2)(x - 3) ( x - 0)(x - 2)(x - 3)


P3 ( x ) = 1 +2
(0 - 1)(0 - 2)(0 - 3) (1 - 0)(1 - 2)(1 - 3)

( x - 0)(x - 1)(x - 3) ( x - 0)(x - 1)(x - 2)


+9 + 28 = 1 + x3
(2 - 0)(2 - 1)(2 - 3) (3 - 0)(3 - 1)(3 - 2)

_______________________________________________________________________________
[Link] Boumediene -62-
Chapitre 4 Interpolation polynômiale

3) Par le Polynôme de Newton:

§ Le polynôme de la figure 4.1 passant par les points (0, 1), (1, 2), (2, 9) et
(3, 28) est (voir tableau 4.1):

xi f (x i ) f [ x i , x i +1 ] f [ x i , x i +1 , x i +2 ] f [ x i , x i +1 , x i +2 , x i +3 ]

0 1
1
1 2 3
7
1
2 9
6
19
3 28

p 3 ( x ) = 1 + 1 ( x - 0 ) + 3 ( x - 0 )( x - 1 ) + 1 ( x - 0 )( x - 1 )( x - 2 ) = x 3 + 1
70

60

50

40

30

20

10

-10
-1 -0.5 0 0.5 1 1.5 2 2.5 3 3.5 4

Figure 4.1 : Courbe du polynôme p 3 ( x ) = x 3 + 1 .

§ Le polynôme de la figure 4.2 passant par les points (2, 1), (0, -1), (5, 10)
et (3, -4) (voir tableau 4.1):

xi f (x i ) f [ x i , x i +1 ] f [ x i , x i +1 , x i +2 ] f [ x i , x i +1 , x i +2 , x i +3 ]

2 1
1
0 -1 0.4
2.2
1.2
5 10
1.6
7
3 -4

p 3 ( x ) = 1 + 1 ( x - 2 ) + 0 . 4 ( x - 2 )( x - 0 ) + 1 . 2 ( x - 2 )( x - 0 )( x - 5 )
3
= 1.2x - 8x 2 + 2.2x - 1

_______________________________________________________________________________
[Link] Boumediene -63-
Chapitre 4 Interpolation polynômiale

10

-5

-10

-15

-20

-25
-1 0 1 2 3 4 5

3
Figure 4.2 : Courbe du polynôme p 3 ( x ) = 1.2x - 8x 2 + 2.2x - 1 .

4.7. Applications par code MATLAB sur les méthodes de calcul


numériques de l’interpolation polynômiale

4.7.1. Applications pour valider les programmes

On doit calculer le polynôme (figure 4.1) passant par les points (0, 1), (1, 2), (2,
9) et (3, 28). Le polynôme recherche est donc: a = [1 0 0 1]T º P3 ( x ) = x 3 + 1

4.7.2. Application MATLAB de la matrice de Vandermonde

[Link]. Algorithme de la matrice de Vandermonde

Etape 1 Donner : x = [x 0 , x1 , x 2 ,L , x n ] et y = [ y 0 , y1 , y 2 ,L , y n ] .
é1 x0 x 20 x 30 L x n0 ù
ê ú
ê1 x1 x12 x13 L x 1n ú
Etape 2 Former la matrice Vandermonde : ê1 x2 x 22 x 32 L x n2 ú .
ê ú
êM M M M O M ú
ê1 xn x 2n x 3n L x nn úû
ë
Etape 3 Utilisant la méthode d’Elimination de Gauss trouvé le polynôme
requis.
Etape 4 Tracer la fonction (polynôme) et positionner les points (x, y).

[Link]. Programme MATLAB de la matrice de Vandermonde

clear all; close all; clc

x=[0 1 2 3] % Données : x = [x0 x1 ... xn].


y=[1 2 9 28] % y = [y0 y1 ... yn].

n=length(x); % La taille de la matrice de Vandermonde.


V=ones(n); % Former & remplir par « for » la matrice de
% Vandermonde.

_______________________________________________________________________________
[Link] Boumediene -64-
Chapitre 4 Interpolation polynômiale

for j=1:n-1
V(:,j+1)=x.^j;
end
V
% Calculer en utilisant « Elimination de Gauss » le polynôme P(x).
[P] = elimgauss(V,y')

% Tracer la fonction (polynôme) et positionner les points (x, y).

xx = [0:0.02:3]; yy = polyval(P,xx);
plot(xx,yy,'r',x,y,'*'), grid

4.7.3. Application MATLAB de l’interpolation de Lagrange

[Link]. Algorithme de l’interpolation de Lagrange

Etape 1 Donner : x = [ x 0 , x 1 , x 2 ,L , x n ] et y = [ y 0 , y1 , y 2 ,L , y n ] .
Etape 2 Calculer les polynômes de Lagrange Li ( x ) par la formule
suivante :
n (x - x j )
L i (x) = Õ
j=0
j¹i
(x i - x j )

Etape 3 Calculer le polynôme requis tel que :


n
Pn (x ) = å f (x ) × L (x ),
i= 0
i i i = 0,1,...,n

Etape 4 Tracer la fonction (polynôme) et positionner les points (x, y).

[Link]. Programme MATLAB de l’interpolation de Lagrange

1) Programme principal généralisé :


clear all; close all; clc

x=[0 1 2 3] % Données : x = [x0 x1 ... xn].


y=[1 2 9 28] % y = [y0 y1 ... yn].

%Calculer en utilisant la fonction « lagrange.m » les


%coefficients de Lagrange « Li(x)» et le polynôme « f(x)»

[L,f] = lagrange (x,y)

% Tracer la fonction (polynôme) →(f) et positionner les


points (x, y).
xx = [0: 0.02 : 3];yy = polyval(f,xx);
plot(xx,yy,'r', x,y,'*'),grid

_______________________________________________________________________________
[Link] Boumediene -65-
Chapitre 4 Interpolation polynômiale

2) Programme MATLAB du fichier fonction:

v Fichier : lagrange.m

function [L,f] = lagrange(x,y)

n = length(x);

% Calculer les coefficients de Lagrange « Li(x)» →(L).


f = 0;
for i = 1:n
P = 1;
for k = 1:n
if k ~= i
% c=conv(a,b) : c Le produit convolution entre le polynôme
% a et le polynôme b.

P = conv(P,[1 -x(k)])/(x(i)-x(k));
end
end
L(i,:) = P ; % Afficher les polynômes de Lagrange « Li(x)»
f = f + y(i)*P ; % Calculer le polynôme requis « f(x) »
end

4.7.4. Application MATLAB du polynôme de Newton

[Link]. Algorithme du polynôme de Newton

Etape 1 Donner : x = [ x 0 , x 1 , x 2 ,L , x n ] et y = [ y 0 , y1 , y 2 ,L , y n ] .
Etape 2 Calculer les coefficients a n utilisant la table des différences
divisées par la formule suivante :

f [ x1 , x 2 ,L , x n ] - f [ x 0 , x 1 , x 2 , L , x n -1 ]
f [ x 0 , x1 , x 2 , L , x n ] =
(x n - x 0 )

Etape 3 Calculer le polynôme requis tel que :


pn ( x ) = a0
+ a1 ( x - x 0 ) + a 2 ( x - x 0 )(x - x 1 )
+ a 3 ( x - x 0 )( x - x1 )(x - x 2 )
+L
+ a n (x - x 0 )(x - x 1 )(x - x 2 ) L (x - x n -1 )

Etape 4 Tracer la fonction (polynôme) et positionner les points (x, y).

_______________________________________________________________________________
[Link] Boumediene -66-
Chapitre 4 Interpolation polynômiale

[Link]. Programme MATLAB du polynôme de Newton

1) Programme principal généralisé :


x=[0 1 2 3] % Données : x = [x0 x1 ... xn].
y=[1 2 9 28] % y = [y0 y1 ... yn].

%Calculer par la fonction « newton.m » les coefficients


%« an » utilisant la table des différences divisées.
%Calculer le polynôme requis «P(x)»

[DD,P] = newton (x,y)

% Tracer la fonction (polynôme) →(P) et positionner les


% points (x, y).
xx = [0:0.02:5]; yy = polyval(P,xx);
plot(xx,yy,'r',x,y,'*'), grid

2) Programme MATLAB du fichier fonction:

v Fichier : newton.m

function [DD,P] = newton(x,y)


n = length(x);
% Calculer les coefficients différences divisées « an »
DD = zeros(n);
DD(1:n,1) = y';
for k = 2:n
for m = 1: n+1-k
DD(m,k) = (DD(m+1,k-1) - DD(m,k-1))/(x(m+k-1)- x(m));
end
end

% Calculer le polynôme requis «P(x)»


a = DD(1,:);
P = a(n);
for k = n-1:-1:1
P = [P a(k)] - [0 P*x(k)];
end

_______________________________________________________________________________
[Link] Boumediene -67-
Chapitre 5

Intégration Numérique

[Link] Boumediene
Chapitre 5 Intégration Numérique

5.1. Introduction

L’objet de ce chapitre est de décrire quelques méthodes numériques (Newton-


Cotes, Trapèzes généralisée, Simpson généralisée) permettant d’évaluer des
intégrales de fonctions dont les valeurs sont connues en un nombre fini de
points [24].

Dans la pluparts des problèmes rencontrés par les ingénieurs, il est souvent
b
difficile de calculer l’intégrale ò f ( x ) dx , où soit la fonction primitive n’existe
a

pas, soit elle est très compliquer, où bien la fonction f (x ) est donnée sous
forme tabulée c.-à-d. f (x ) n’est connue seulement en quelque points
( x i , y i )=i 0 ,1, ..., n (cas des expériences aux labos). Les mêmes problèmes se posent
par le calcul des intégrales multiples. Pour ces raisons on remède à l’utilisation
de l’Intégration Numérique (approchée) [25].

5.2. Définitions

L’intégration numérique est basée principalement sur la relation [24, 25]:

xn xn xn

ò f ( x ) dx
x0
= ò P (x ) dx
x0
n + ò E ( x ) dx
x0
n (5.1)

où Pn (x ) est un polynôme d'interpolation et En ( x ) est l’erreur qui y est


associée.
b
La façon la plus simple d’approchée un intégrale ò f ( x ) dx est de remplacer la
a

fonction f ( x ) par un polynôme Pn ( x ) qui passe par les points ( x i , y i ) i =0 ,1, ..., n
avec y i = f ( x i ) .
b
L’intégration numérique consiste à remplacer l’intégrale ò f ( x ) dx , par
a

plusieurs valeurs de la fonction à intégrer c.-à-d. par une somme polynômiale.

5.3. Calcul de l’intégrale par la formule de Newton-Côtes

Soit f ( x ) une fonction continue sur [a, b] choisi sous des points x i ( i = 0 ,1, ..., n )

équidistants dans le segment [a, b] et soit y i = f (x i ) i = 0 ,1, ..., n les valeurs de


f ( x ) en ces points [24, 26].

Donc l’intégrale peut être approximer par la formule de Newton-Côtes:


b n

ò f (x ) dx » ( b - a)å Hi × f (x i )
a i =0
(5.2)

_______________________________________________________________________________
[Link] Boumediene -69-
Chapitre 5 Intégration Numérique

n
( -1)n -i q(q - 1)(q - 2) L (q - n)
avec : H i = ò
n! (n - i ) ! 0 q-i
dq i = 0,1, 2,L , n

où : H i les coefficients de Côtes,

H i = H n-i (5.3)

et
n

åH
i =0
i =1 (5.4)

b-a
Avec : h = (le pas).
n
Les coefficients de Côtes sont donnés dans le tableau 5.1.

Tableau 5.1 : Coefficients de Côtes pour n = 6.


n H0 H1 H2 H3 H4 H5 H6
1 1/2 1/2 −

2 1/6 4/6 1/6 −

3 1/8 3/8 3/8 1/8 −
4 7/90 32/90 12/90 32/90 7/90
5 19/288 75/288 50/288 50/288 75/288 19/288
6 41/840 216/840 27/840 272/840 27/840 216/840 41/840

Nous allons montrer que lorsque la fonction f à intégrer est suffisamment


régulière, l’erreur d’intégration numérique peut s’exprimer de manière assez
simple en fonction d’une certaine dérivée de f. Auparavant, nous aurons besoin
de quelques rappels d’analyse.

5.4. Calcul de l’intégrale par la formule des Trapèzes généralisée

b-a
On va diviser l’intervalle [a, b] à n segments de longueur h = et soit
n
y i = f ( x i ) i = 0 ,1, ..., n [15, 20, 26]:

[a, b] = [x 0 , x1 ] , [x1 , x 2 ] , ..., [x n-1 , x n ] (5.4)

Nous avons la formule des trapèzes généralisée : h = xi +1 - xi

b
h n
ò f (x ) dx » å (y i -1 + y i ) + e
2 i =1
a

h
= [(y 0 + y n ) + 2(y1 + y 2 + L + y n-1 )] + e (5.5)
2
æ y 0 + y n n -1 ö
= hçç + å y i ÷÷ + e
è 2 i =1 ø

_______________________________________________________________________________
[Link] Boumediene -70-
Chapitre 5 Intégration Numérique

( b - a)3 f ¢¢( x * )
Le reste : e £ £ e avec : x * Î [a, b].
12 × n 2

5.5. Calcul de l’intégrale par la formule de Simpson généralisée

Soit n=2m (m nombre pair) et y i = f ( x i ) i =0 ,1, ..., n les valeurs de f (x ) aux points
b-a
équidistants de pas h = , l’intégrale peut être approximer par la formule
2m
de Simpson généralisée [15, 25, 26]:

b
h
ò f (x ) dx » 3 [(y
a
0 + y n ) + 4(y1 + y 3 + L + y n-1 ) + 2(y 2 + y 4 + L + y n -2 )] + e

(5.6)
5 (4 ) *
( b - a) f (x )
Le reste : e £ £ e avec : x * Î [a, b].
180 × n 4

5.6. Exemple d’application

Soit f ( x ) une fonction donnée par le tableau :

xi 0.25 0.50 0.75 1.00 1.25


f (x i ) 0.630 0.794 0.909 1.000 1.15

1.25

v Trouver une valeur approchée de l’intégrale I = ò f (x ) dx , en utilisant :


0.25

1- La formule de Newton-Côtes.
2- La formule des Trapèzes généralisée.
3- La formule de Simpson généralisée.

· Solutions :

1) Par la formule de Newton-Côtes :

Puisque le pas d’intégration est h = 0.25 . De la relation :

b-a 1.25 - 0.25


n= Û n= =4
h 0.25

Les coefficients de Côtes peuvent être calculés de la formule, ou peuvent être


pris directement du tableau 5.1 des coefficients de Côtes. Les valeurs des
coefficients où n=4 sont :

7 32 12
H0 = H4 = ; H1 = H3 = et H2 =
90 90 90

_______________________________________________________________________________
[Link] Boumediene -71-
Chapitre 5 Intégration Numérique

A partir de la formule de Newton-Côtes suivante :


b 4

ò f (x ) dx » ( b - a )å H i × f ( x i )
a i =0

1.25
I= ò f (x ) dx
0.25

é7 32 12 32 7 ù
= (1.25 - 0.25) ´ ê (0.63) + (0.794 ) + (0.909 ) + (1) + (1.15)ú
ë 90 90 90 90 90 û
= 0.89751

2) Par la formule des Trapèzes généralisée:

A partir de la formule des trapèzes généralisée :

b-a 1.25 - 0.25


n= Û n= =4
h 0.25

Puisque le pas d’intégration est h = 0.25 . De la relation :

b
æ y0 + y n 3
ö
òa f ( x ) dx » h ç
ç 2 + å y i ÷÷
è i =1 ø
1.25
æ 0.63 + 1.15 ö
I= ò
0.25
f (x ) dx = (0.25) ´ ç
è 2
+ (0.794 + 0.909 + 1) ÷ = 0.89825
ø

3) Par la formule de Simpson généralisée:

Puisque le pas d’intégration est h = 0.25 . De la relation :

b-a 1.25 - 0.25


n=
Û n= =4
h 0.25
A partir de la formule de Simpson généralisée :
b
h
ò f (x ) dx » [(y 0 + y 4 ) + 4(y1 + y 3 ) + 2(y 2 )]
a
3

1.25
æ 0.25 ö
I= ò f ( x ) dx = ç ÷ ´ ((0.63 + 1.15) + 4 ´ (0.794 + 1) + 2 ´ (0.909 )) = 0.89783
0.25 è 3 ø

_______________________________________________________________________________
[Link] Boumediene -72-
Chapitre 5 Intégration Numérique

5.7. Applications MATLAB sur calcul numérique de l’intégrale

5.7.1. Applications pour valider les programmes

Soit f ( x ) une fonction donnée par le tableau :

xi 0.25 0.50 0.75 1.00 1.25


f (x i ) 0.630 0.794 0.909 1.000 1.15

1.25

v Trouver une valeur approchée de l’intégrale I = ò f (x ) dx , en utilisant le


0.25

programmes MATLAB des trois méthodes.

5.7.2. Application MATLAB de la formule de Newton-Côtes

[Link]. Algorithme de la formule de Newton-Côtes

Etape 1 Donner : x = [ x1 , x 2 ,L , x n ] équidistant, y = [ y1 , y 2 , L , y n ] et les


coefficients de Côtes.
Etape 2 Calcul de « h » (le pas).
Etape 3 Former la formule de Newton-Côtes :
n
I = ( x n - x1 ) × å H i × y i
i =1

Etape 4 Utilisant la formule de Newton-Côtes calculer la valeur « I » de


l’intégrale requis.

[Link]. Programme MATLAB de la formule de Newton-Côtes

clear all; close all; clc


% Les coefficients de Côtes sont donnés dans la table
suivante :
H=[ 1/2 1/2 0 0 0 0 0
1/6 4/6 1/6 0 0 0 0
1/8 3/8 3/8 1/8 0 0 0
7/90 32/90 12/90 32/90 7/90 0 0
19/288 75/288 50/288 50/288 75/288 19/288 0
41/840 216/840 27/840 272/840 27/840 216/840 41/840]

% Soit f(x) une fonction donnée par le tableau :


x=[0.25 0.50 0.75 1.00 1.25] % équidistant
y=[0.630 0.794 0.909 1.000 1.15]

n=length(x)
h=x(2)-x(1) % h c'est le pas.
a=x(1)
b=x(n)

_______________________________________________________________________________
[Link] Boumediene -73-
Chapitre 5 Intégration Numérique

N=(b-a)/h % N la valeur pour déterminer Hi.


Hi=H(N,1:N+1) % Hi Les coefficients de Côtes suivant N.

% Calcul de l'intégrale à partir de la formule de Newton-


% Côtes :
I=(b-a)*(y*Hi')

5.7.3. Application MATLAB de la formule des Trapèzes généralisée

[Link]. Algorithme de la formule des Trapèzes généralisée

Etape 1 Donner : x = [ x1 , x 2 ,L , x n ] équidistant et y = [ y1 , y 2 , L , y n ] .


Etape 2 Calcul de « h » (le pas).
Etape 3 Former la formule de Trapèzes généralisée :

æ y1 + y n n -1 ö
I = h × çç + å y i ÷÷
è 2 i =2 ø
Etape 4 Utilisant la formule de Trapèzes généralisée calculer la valeur
« I » de l’intégrale requis.

[Link]. Programme MATLAB de la formule des Trapèzes généralisée

clear all; close all; clc

% Soit f(x) une fonction donnée par le tableau :


x=[0.25 0.50 0.75 1.00 1.25] % équidistant
y=[0.630 0.794 0.909 1.000 1.15]

n=length(x)
h=x(2)-x(1) % h c'est le pas.

% Calcul de l'intégrale à partir de la formule de Trapèzes


% généralisée:
I=h*(((y(1)+y(n))/2)+sum(y(2:n-1)))

5.7.4. Application MATLAB de la formule de Simpson généralisée

[Link]. Algorithme de la formule de Simpson généralisée

Etape 1 Donner : x = [ x1 , x 2 ,L , x n ] équidistant et y = [ y1 , y 2 , L , y n ] .


Etape 2 Calcul de « h » (le pas).
Etape 3 Former la formule de Simpson généralisée :

h
I= [(y1 + y n ) + 2(y 2 + y 4 + L + y n-2 ) + 4(y 3 + y 5 + L + y n-1 )]
3
Etape 4 Utilisant la formule de Simpson généralisée calculer la valeur « I »
de l’intégrale requis.

_______________________________________________________________________________
[Link] Boumediene -74-
Chapitre 5 Intégration Numérique

[Link]. Programme MATLAB de la formule de Simpson généralisée

clear all; close all; clc

% Soit f(x) une fonction donnée par le tableau :


x=[0.25 0.50 0.75 1.00 1.25] % équidistant
y=[0.630 0.794 0.909 1.000 1.15]

n=length(x)
h=x(2)-x(1) % h c'est le pas.

% Calcul de l'intégrale à partir de la formule de Simpson


% généralisée:
I=(h/3)*((y(1)+y(n))+4*sum(y([Link]n-1))+2*sum(y([Link]n-2)))

_______________________________________________________________________________
[Link] Boumediene -75-
Chapitre 6

Calcul Numérique des


Valeurs Propres

[Link] Boumediene
Chapitre 6 Calcul Numérique des Valeurs Propres

6.1. Introduction

Il s'agit de calculer les valeurs propres de matrices carrées d'ordre « n » à


coefficients réels; cependant, les méthodes proposées peuvent s'appliquer aux
matrices à éléments complexes. Nous allons examiner diverses méthodes qui
donnent des résultats plus ou moins intéressants selon la taille de la matrice A
ou son conditionnement [27].

L'équation caractéristique d'une matrice d'ordre n est un polynôme de degré n,


et les racines de ce polynôme sont les valeurs propres de la matrice A.

6.2. Définitions

Soit A une matrice carrée à n lignes et n colonnes dont les éléments aij sont des
nombres réels.

L’équation caractéristique de la matrice A est un polynôme de degré n, et les


racines λ de ce polynôme sont les valeurs propres de la matrice A.

Le problème est donc la résolution de l'équation caractéristique ou séculaire:

A∙X = λ∙X (6.1)

où X s'appelle vecteur propre et A valeur propre.

Déterminant pour une matrice 3×3 :

6.3. Méthode de Déterminant

La matrice (A - λ In) est appelée matrice caractéristique de A, où un scalaire


λ est appelé valeur propre de A et In matrice identité de degré n×n. Le
polynôme Pn(λ) = det (A - λ In) est appelé polynôme caractéristique de A et
l’équation Pn(λ) = 0 est appelée équation caractéristique de A [27].
Les valeurs propres [λ1, λ2,…, λn] de A sont les solutions de l’équation
caractéristique :

det (A - λ In) = 0 (6.2)

_______________________________________________________________________________
[Link] Boumediene -77-
Chapitre 6 Calcul Numérique des Valeurs Propres

é a11 a12 a13 L a1 n ù é1 0 0 0 0ù


êa a 22 a 23 L a 2n úú ê0 1 0 0 0úú
ê 21 ê
êa 31 a 32 a 33 L a 3n ú - [l l l L l] × ê0 0 1 0 0ú = 0 (6.3)
ê ú ê ú
ê M M M O M ú ê0 0 0 O 0ú
êëa n1 a n2 a n3 L a nn úû êë0 0 0 0 1 úû

a11 - l a12 a13 L a1n


a 21 a 22 - l a 23 L a 2n
Pn(λ) = det (A - λ In) = a 31 a 32 a 33 - l L a 3n =0 (6.4)
M M M O M
a n1 a n2 a n3 L a nn - l

6.4. Exemple d’application utilisant la Méthode de Déterminant

Trouver, en utilisant la méthode de Déterminant les valeurs propres de :


é 1 - 3ù
A=ê ú ;
ë- 2 2 û

· Solution :

é1 - l - 3 ù
La matrice caractéristique : A - lI 2 = ê ú ;
ë - 2 2 - lû

Le polynôme caractéristique de A :

P2 (l ) = det(A - lI 2 ) = (1 - l )(2 - l ) - 6 = l2 - 3l - 4 .
L’équation caractéristique de A :

l2 - 3l - 4 = 0 .

Les valeurs propres : l = -1 et l = 4 .

6.5. Méthode de Le Verrier

L'équation caractéristique d'une matrice d'ordre n est un polynôme de degré n,


et les racines de ce polynôme sont les valeurs propres de la matrice. Si nous
pouvons obtenir les coefficients du polynôme caractéristique, la méthode
d’élimination de Gauss nous permettra d'en calculer les racines [15, 28]. La
solution repose sur les relations de Newton (6.5) qui établissent une expression
entre les coefficients d'un polynôme et les Sq (somme des racines du polynôme
élevées à la puissance q).

Soit A une matrice carrée à n lignes et n colonnes, les valeurs propres de cette
matrice sont les racines λ de ce polynôme :
Pn (l ) = a 0 ln + a1 ln -1 + a 2 ln-2 + L + a n l0 = 0 ; avec ( a 0 = 1 ) (6.5)

_______________________________________________________________________________
[Link] Boumediene -78-
Chapitre 6 Calcul Numérique des Valeurs Propres

La détermination des coefficients a i est par la formule suivante :


ìa 0 S1 + a1 = 0
ï
ïa 0 S 2 + a1 S1 + 2a 2 = 0
ïï a 0 S 3 + a1 S 2 + a 2 S1 + 3a 3 = 0
í (6.6)
ï a 0 S 4 + a1 S3 + a 2 S 2 + a 3 S1 + 4a 4 = 0
ï LLL LL LL
ï
ïî a 0 S m + a1 S m-1 + a 2 S m -2 + L + a n -1 S1 + qa n = 0 (m = n = q )

Avec :
a i (i=0,1, …, n) : les coefficients d’un polynôme caractéristique Pn ( l ) .
S q (q= 1,2, …, m) : la somme du diagonal de la matrice A de puissance 1 à m.

Nous savons que la trace de la matrice A est égale a la somme de ses valeurs
propres ; de même, ta trace de la matrice Aq est égale à la somme des valeurs
propres de A élevées à la puissance q. Par conséquent, la trace de la matrice Aq
est égale a Sq. (On rappelle que la trace d'une matrice est la somme de ses
éléments sur la diagonale principale.) Il suffit de calculer la trace de A puis de
multiplier A par A soit B = AA et de calculer la trace de B. Ensuite, on
multipliera B par A puis on calculera la trace de B, et ainsi de suite jusqu'a
obtenir la valeur de Sn.

6.6. Exemple d’application utilisant la Méthode de Le Verrier

Trouver, en utilisant la méthode de Le Verrier les valeurs propres de :


é1 2 2ù
A = êê2 1 2úú ;
êë2 2 1 úû

· Solution :

Formons les matrices A 2 et A 3 . Nous avons par un calcul direct :

é1 2 2ù é9 8 8ù é41 42 42ù
ê ú 2 ê ú
A = ê2 1 2ú , A = A × A = ê8 9 8ú et A = A × A = êê42 41 42úú
3 2

êë2 2 1 úû êë8 8 9úû êë42 42 41úû

ì S1 = 1 + 1 + 1 = 3
ï2 3
Nous obtenons de A , A , A : íS 2 = 9 + 9 + 9 = 27
ïS = 41 + 41 + 41 = 123
î 3

_______________________________________________________________________________
[Link] Boumediene -79-
Chapitre 6 Calcul Numérique des Valeurs Propres

ì
ï a1 = -S1
ïï 1
Les formules : í a 2 = - (S2 + a1S1 ) ;
ï 2
1
ïa 3 = - (S3 + a1 S2 + a 2 S1 )
îï 3

ì
ï a1 = - 3
ïï 1
S’écrivent : í a 2 = - (27 + a1 3)
ï 2
1
ïa 3 = - (123 + a1 27 + a 2 3)
îï 3

D’où on tire les valeurs a 1 = -3 , a 2 = -9 et a 3 = -5 .

Ainsi le polynôme caractéristique de la matrice A s’écrit :

P3 ( l ) = l3 + a1 l2 + a 2 l + a 3 = l3 - 3l2 - 9l - 5 = 0

Les solutions de l’équation l3 - 3l2 - 9l - 5 = 0 sont les valeurs propres de la


matrice A et sont :
l1 = -1 , l 2 = -1 et l 3 = 5 .

6.7. Méthode de Krylov

Soit A une matrice carrée à n lignes et n colonnes, les valeurs propres de cette
matrice sont les racines λ de ce polynôme [15, 29]:

Pn ( l ) = a 0 ln + a1 ln -1 + a 2 ln -2 + L + a n l0 = 0 ; avec ( a 0 = 1 ) (6.7)

La détermination des coefficients a i est par la formule suivante :

a1 Y ( n -1) + a 2 Y ( n -2 ) + L + a n -1 Y (1) + a n Y ( 0 ) = - Y ( n ) (6.8)

Puisque la taille de A est n×n, alors nous avons (n) vecteurs de longueur (n) a
formés :

Y (1 ) = AY (0 )
Y (2 ) = AY (1)
(6.9)
M
Y (n ) = AY (n -1 )
é1ù
Prenons comme vecteur initial Y (0 )
= êê M úú
êë1úû

_______________________________________________________________________________
[Link] Boumediene -80-
Chapitre 6 Calcul Numérique des Valeurs Propres

En fin, on obtient :
é a1 ù
êa ú
[Y ( n -1 )
Y (n -2 ) L Y ( 0 ) ] × ê ú = - Y (n )
2

êM ú
(6.10)
ê ú
ëa n û

Le vecteur X est arbitraire, et on le choisit très souvent avec toutes ses


composantes égales à l'unité. Reste à déterminer a 0 = ( -1)n .
Cette méthode n'exigeant que peu de calculs est fort attrayante,
malheureusement elle conduit a un système mal conditionne lorsque n atteint
et dépasse quelques unités de l’ordre de 8 a 10.

La solution de ce système est les coefficients du polynôme caractéristique


Pn ( l ) de A. Les solutions de l’équation caractéristique Pn (l ) = 0 sont les
valeurs propres de A.

6.8. Exemple d’application utilisant la Méthode de Krylov

Trouver, en utilisant la méthode de Krylov les valeurs propres de :


é2 - 1 - 2ù
A = êê1 1 1 úú ;
êë2 3 - 1 úû

· Solution :

Le polynôme caractéristique de A à la forme :

P3 ( l ) = l3 + a1 l2 + a 2 l + a 3

é1ù
Prenons comme vecteur initial Y (0 )
= êê1úú
êë1úû
Puisque la taille de A est 3×3, alors nous avons trois (3) vecteurs a formés :

é2 - 1 - 2ù é1ù é- 1ù
Y = AY = êê1 1
(1 )
1 úú êê1úú = êê 3 úú
(0 )

êë2 3 - 1 úû êë1úû êë 4 úû
é2 - 1 - 2ù é- 1ù é- 13ù
Y = AY = êê1 1
(2 ) (1 )
1 úú êê 3 úú = êê 6 úú
êë2 3 - 1 úû êë 4 úû êë 3 úû
é2 - 1 - 2ù é- 13ù é- 38ù
Y (3 )
= AY (2 )
= êê1 1 1 úú êê 6 úú = êê - 4 úú
êë2 3 - 1 úû êë 3 úû êë - 11 úû

_______________________________________________________________________________
[Link] Boumediene -81-
Chapitre 6 Calcul Numérique des Valeurs Propres

En fin, on obtient : a1 Y ( 2 ) + a 2 Y (1) + a 3 Y ( 0 ) = - Y ( 3 )

é é - 13ù é - 1ù é1ù ù éa1 ù é- 38ù ì- 13a1 - a 2 + a 3 = 38


ê ê1ú ú × êa ú = - ê - 4 ú ï
Où encore : ê êê 6 úú ê3ú
ê ú ê úú ê 2 ú ê ú º í 6a1 + 3a 2 + a 3 = 4
ê êë 3 úû êë 4 úû êë1úû úû êëa 3 úû êë - 11 úû ï 3a + 4a + a = 11
ë î 1 2 3

La solution de ce système est les valeurs a1 = -2 , a 2 = 1 et a 3 = 13 .

Ainsi le polynôme caractéristique de la matrice A s’écrit :

P3 (l ) = l3 + a1 l2 + a 2 l + a 3 = l3 - 2l2 + l + 13 = 0

Les solutions de l’équation l3 - 2l2 + l + 13 = 0 sont les valeurs propres de la


matrice A et sont :
l1 = 1.8681 + 1.9993i , l 2 = 1.8681 - 1.9993i et l 3 = -1.7363 .

6.9. Applications MATLAB sur le calcul des valeurs propres

6.9.1. Applications pour valider les programmes

Écrire le programme généralisé des trois méthodes qui calculent les valeurs
propres utilisant MATLAB des matrices suivantes :

é1 2 2ù é2 - 1 - 2ù
A = êê2 1 2úú ; A = êê1 1 1 úú .
êë2 2 1 úû êë2 3 - 1 úû

6.9.2. Application MATLAB de la méthode de Le Verrier

[Link]. Algorithme de la méthode de Le Verrier

Etape 1 Données : la matrice A de taille ( n ´ n ) et la constante a 0 = 1 .


Etape 2 Calcul de la somme des diagonales de la matrice A de puissance
1 à n.
Etape 3 Pour déterminer les coefficients a i , on doit transformer les
formules de la relation (6.6) à un système d’équation linéaire
A*X=b et construire la matrice A et le vecteur b.
Etape 4 Utilisant la méthode d’Elimination de Gauss trouvé la solution
( a i ) de ce système (6.6).
Etape 5 Former le polynôme caractéristique de la matrice "A".
Etape 6 Calculer les valeurs propres de la matrice "A".

_______________________________________________________________________________
[Link] Boumediene -82-
Chapitre 6 Calcul Numérique des Valeurs Propres

[Link]. Programme MATLAB de la méthode de Le Verrier

clear all; close all; clc

% Données:
A=[1 2 2;2 1 2; 2 2 1] % La matrice "A".
a0=1;

% Calcul de la somme des diagonales de la matrice A de


% puissance 1 à n:
n=length(A);
for i=1:n
S(i,1)=sum(diag(A^i));
end
S
% Calcul de la matrice M:
for i=1:n
for j=1:n
if i==j
M(i,j)=i; % Diagonal de la matrice M égale à
% 1 jusqu'au n.
elseif i<j
M(i,j)=0; % La partie superieure de la
% matrice M égale à 0.
else
M(i,j)=S(i-j); % La partie inferieure égale à la
% somme des diagonals.
end

end
end

% Formation du système d'équation linéaire a*X=b:


a=M % Le matrice "a" du système.
b=-S % Le vecteur "b" du système.

% La solution de ce système est par la fonction


% d’Elimination de Gauss "elimgauss.m"
[P] = elimgauss(a,b)

% Le polynôme caractéristique de la matrice "A".


poly_caract=[a0 P']

% Les valeurs propres de la matrice "A".


Val_prop=roots(poly_caract)

_______________________________________________________________________________
[Link] Boumediene -83-
Chapitre 6 Calcul Numérique des Valeurs Propres

6.9.3. Application MATLAB de la méthode de Krylov

[Link]. Algorithme de la méthode de Krylov

Etape 1 Données : la matrice A de taille ( n ´ n ), la constante a 0 = 1 et le


vecteur initial Y(0)=1.
Etape 2 Calculer des "n" vecteurs Y(i)=A*Y(i-1).
Etape 3 Formation du système d'équation linéaire a*X=b.
Etape 4 Utilisant la méthode d’Elimination de Gauss trouvé la solution
( a i ) du système.
Etape 5 Former le polynôme caractéristique de la matrice "A".
Etape 6 Calculer les valeurs propres de la matrice "A".

[Link]. Programme MATLAB de la méthode de Krylov

clear all; close all; clc

% Données:
A=[2 -1 -2;1 1 1;2 3 -1] % La matrice "A".
a0=1 ;
n=length(A);
Y=ones(n,1); % Vecteur initial Y(0).

% Calcul des "n" vecteurs:


for i=1:n
Y=A*Y; % Calcul des "n" vecteurs Y(1)=A*Y(0).
YY(:,i)=Y; % Matrice où il y à les "n" vecteurs.
end
Y=ones(n,1) % Vecteur initial Y(0).
YY % Affichage de la matrice où il y à les "n"
% vecteurs.

% Formation du système d'équation linéaire a*X=b:


a=[YY(:,n-1:-1:1),Y] % Le matrice "a" du système.
b=-YY(:,n) % Le vecteur "b" du système.

% La solution de ce système est par la fonction


% d’Elimination de Gauss "elimgauss.m".
[P] = elimgauss(a,b)

% Le polynôme caractéristique de la matrice "A".


poly_caract=[a0 P']

% Les valeurs propres de la matrice "A".


Val_prop=roots(poly_caract)

_______________________________________________________________________________
[Link] Boumediene -84-
Chapitre 7

Résolution d’un Système


d’Équations Non Linéaires

[Link] Boumediene
Chapitre 7 Résolution d’un Système d’Équations Non Linéaires

7.1. Introduction

Les phénomènes non linéaires sont extrêmement courants en pratique. Ils sont
sans doute plus fréquents que les phénomènes linéaires. Dans cette section,
nous examinons les systèmes non linéaires et nous montrons comment les
résoudre a l'aide d'une suite de problèmes linéaires [16, 17], auxquels on peut
appliquer diverses techniques de résolution comme la méthode d’élimination
de Gauss.

Le problème consiste à trouver le ou les vecteurs X = [x1 x3 L x n ]


T
x2
vérifiant les n équations non linéaires suivantes:

f1 ( x1 , x 2 ,L , x n ) = 0
f 2 (x 1 , x 2 ,L , x n ) = 0
f 3 ( x1 , x 2 ,L , x n ) = 0 (7.1)
M
f n ( x1 , x 2 ,L , x n ) = 0

où les fi sont des fonctions de n variables que nous supposons différentiables.

Contrairement aux systèmes linéaires, il n'y a pas de condition simple associée


aux systèmes non linéaires qui permette d'assurer l'existence et l'unicité de la
solution. Le plus souvent, il existe plusieurs solutions possibles et seul le
contexte indique laquelle est la bonne [29, 30].

Les méthodes de résolution des systèmes non linéaires sont nombreuses. Pour
éviter de surcharger notre document, nous ne présentons que la méthode la
plus importante et la plus utilisée en pratique, soit la méthode de Newton-
Raphson.

7.2. Méthode de Newton-Raphson

L'application de cette méthode à un système de deux équations non linéaires


est suffisante pour illustrer le cas général. Il serait également bon de réviser le
développement de la méthode de Newton-Raphson pour une équation non
linéaire puisque le raisonnement est le même pour les systèmes [16, 29].

Soit le Système d’équations non linéaires :

fi ( x1 , x 2 ,L , x n ) = 0 avec (i = 1, 2,L , n ) (7.2)

é x1 (0 ) ù
ê ú
où X (0 ) = ê M ú est l’approximation initiale et « ε » c’est la précision du calcul
êx (0 ) ú
ë n û
(la solution est une valeur approchée).

_______________________________________________________________________________
[Link] Boumediene -86-
Chapitre 7 Résolution d’un Système d’Équations Non Linéaires

Les étapes de la méthode sont:

1) Calcul des éléments de la matrice Jacobienne (les dérivées partielles).


2) Formation de la relation au-dessous :

é ù
ê ¶f1 ¶f1 ¶f1 ú
ê ¶x L
ê 1 ¶x 2 ¶x n úú é Dx1 (k +1 ) ù é f1 ù
ê 2¶f ¶ f ¶ f 2 ú
ê ( k +1 ) ú êf ú
2
L ê Dx 2 ú
ê ¶x ¶x 2 ¶x n ú × ê = -ê 2 ú (7.3)
ê M 1 ú êMú
M úú ê
M
ê ¶f M O (k +1 ) ú ê ú
ê n ¶f n L ¶f n ú ëê n
Dx ûú ëf n û
ê ¶x 1 ¶x 2 ¶x n ú
êë 14 4424443
La matrice Jacobienne
úû

(k ) (k ) (k )
3) Remplacement par les valeurs : x1 , x 2 ,L , x n dans la formule (7.3).
(k +1 )
4) Calcul des Dx i a l’aide de la formule (7.3).
5) Obtention la solution approchée X (k +1 ) par la relation suivante :

é x 1 (k ) + Dx 1 ù
ê ú
X (k +1 ) =ê M ú (7.4)
ê x (k ) + Dx ú
ë n nû

6) Si la relation (7.5) atteinte :

é Dx1 £ e ù
ê ú
X (k +1 ) - X (k ) =ê M ú (7.5)
ê Dx n £ e ú
ë û

alors la solution approchée est X (k +1 ) , sinon répéter le calcule pour une autre
itération.

7.3. Exemple d’application utilisant la méthode de Newton-Raphson

Soit le système d’équations non linéaires :

ìx + 2y = 2
í 2 2
î x + 4y = 4
En utilisant la méthode de Newton-Raphson et en partant de
(0 ) é0.2ù
l’approximation initiale X = ê ú , trouver la valeur approchée de la solution
ë1.2 û
avec une précision 10-2.

_______________________________________________________________________________
[Link] Boumediene -87-
Chapitre 7 Résolution d’un Système d’Équations Non Linéaires

· Solution :

(0 ) é0.2ù
L’approximation initiales donnée est X = ê ú et la précision est ε=10-2.
ë1.2 û
Le calcul des éléments de la matrice Jacobienne (les dérivées partielles):

¶f1 ¶f1 ¶f 2 ¶f 2
=1 =2 = 2x = 8y
¶x ¶y ¶x ¶y

La relation (7.3) -voir le cours- s’écrit :

é1 2ù éDx ù é f1 (0.2, 1.2)ù


ê2x 8yú ê Dy ú = - êf (0.2, 1.2)ú
ë û (0.2, 1.2 ) ë û ë2 û
Où encore :
é1 2 ù é Dx ù é- 0.6ù
ê0.4 9.6ú ê Dy ú = ê - 1.8 ú
ë ûë û ë û
Calculons Dx et Dy :

- 0 .6 2 1 - 0.6
- 1.8 9.6 0.4 - 1.8
Dx = = -0.245 et Dy = = -0.177
1 2 1 2
0.4 9.6 0.4 9.6

La première solution approchée X (1 ) est :

(1 ) éx (0 ) + Dx = 0.2 - 0.245ù é - 0.045ù


X = ê (0 ) ú=ê ú
ë y + Dy = 1.2 - 0.177 û ë 1.023 û

(1 ) (0 ) é Dx = 0.245 > 10 -2 ù
Et puisque X -X =ê -2 ú , effectuons une deuxième
êë Dy = 0.177 > 10 úû
itération :

é1 2ù éDx ù éf1 (- 0.045, 1.023)ù


ê2x 8yú ê ú = - êf (- 0.045, 1.023 )ú
ë û (-0.045, 1.023 ) ë Dyû ë2 û

é 1 2 ù é Dx ù é - 0.001ù
ê- 0.09 8.184 ú ê Dy ú = ê- 0.188 ú
ë ûë û ë û

_______________________________________________________________________________
[Link] Boumediene -88-
Chapitre 7 Résolution d’un Système d’Équations Non Linéaires

Calculons Dx et Dy :

- 0.001 2 1 - 0.001
- 0.188 8.184 - 0.09 - 0.188
Dx = = 0.044 et Dy = = -0.022
1 2 1 2
- 0.09 8.184 - 0.09 8.184

La deuxième solution approchée X (2 ) est :

éx (1 ) + Dx = -0.045 + 0.044 ù é- 0.001ù


X (2 ) = ê (1 ) ú=ê ú
ë y + Dy = 1 .023 - 0 .022 û ë 1.001 û

(2 ) (1 ) é Dx = 0.044 > 10 -2 ù
Et puisque : X -X =ê ú
êë Dy = 0.022 > 10 -2 úû

Nous devons effectuer une troisième itération :

é1 2 ù é Dx ù éf1 (- 0.001, 1.001)ù


ê2x 8y ú ê Dy ú = - êf (- 0.001, 1.001)ú
ë û (-0.001, 1.001 ) ë û ë2 û

é 1 2 ù éDx ù é - 0.001ù
ê- 0.002 8.008 ú ê Dyú = ê - 0.008ú
ë ûë û ë û
Calculons Dx et Dy :

- 0.001 2 1 - 0.001
- 0.008 8.008 - 0.002 - 0.008
Dx = = 0.001 et Dy = = -0.001
1 2 1 2
- 0.002 8.008 - 0.002 8.008

La troisième solution approchée X (3 ) est :

(3 ) éx (2 ) + Dx = -0.001 + 0.001ù é0.000ù


X = ê (2 ) ú=ê ú
ë y + Dy = 1.001 - 0.001 û ë1.000 û
Et puisque :
(3 ) (2 ) é Dx = 0.001 £ 10 -2 ù
X -X =ê -2 ú
êë Dy = 0.001 £ 10 úû

Nous avons abouti la précision fixée (10-2), la valeur approchée cherchée est :

é0.000 ù
X=ê ú
ë1.000 û

_______________________________________________________________________________
[Link] Boumediene -89-
Chapitre 7 Résolution d’un Système d’Équations Non Linéaires

7.4. Applications MATLAB sur la résolution d’un système d’équations


non linéaires

7.4.1. Applications pour valider les programmes

Soit les systèmes d’équations non linéaires :


ìx + 2y = 2 ìx 2 - 2x - y + 0.5 = 0
í 2 2 et í 2
î x + 4y = 4
2
îx + 4y - 4 = 0
En utilisant la méthode de Newton-Raphson et en partant de
(0 ) é0.2ù
l’approximation initiale X = ê ú , trouver la valeur approchée de la solution.
ë1.2 û

7.4.2. Application MATLAB de la méthode de Newton-Raphson

[Link]. Algorithme de la méthode de Newton-Raphson

Etape 1 Données : L’approximation initiale

é x 1 (0 ) ù
ê ú
X (0 ) =ê M ú ;
êx (0 ) ú
ë n û

Etape 2 Poser k=1 ;


Le système d’équations f i ( x 1 , x 2 ,L, x n ) = 0 avec (i = 1, 2,L , n ) ;
Les dérivées partielles (la matrice Jacobienne) ;
(k )
Etape 3 Calculer les Dx i à l’aide de la formule (7.3).
Etape 4 La solution approchée X (k +1 ) est par la relation suivante :

é x1 (k ) + Dx 1 ù
ê ú
X (k +1 ) =ê M ú
êx (k ) + Dx ú
ë n nû

Etape 5 Répéter le calcul jusqu’au maximum des itérations.

[Link]. Programme MATLAB de la méthode de Newton-Raphson

1) Programme MATLAB principal (P_Newton_Raphson.m):


clear all; close all; clc
x=[0.2 ; 1.2] % La solution initiale

% la fonction qui calcul le programme de Newton-Raphson


% où l’entrée c(est le vecteur initial x :
x=newton_raphson(x)

_______________________________________________________________________________
[Link] Boumediene -90-
Chapitre 7 Résolution d’un Système d’Équations Non Linéaires

2) Programme MATLAB utilisant le fichier fonction:

v Fichier : (newton_raphson.m)

function x=newton_raphson(x)

for k=1:100 % test d’arret le max d’itérations


% La Matrice Jacobienne et le système d'équations non
% linéaires:

% Le système d'équations non linéaires:


F =[x(1)+2*x(2)-2;
x(1)^2+4*x(2)^2-4];

% La Matrice Jacobienne:
df1_dx1=1;
df1_dx2=2;
df2_dx1=2*x(1);
df2_dx2=8*x(2);
J=[df1_dx1 df1_dx2;
df2_dx1 df2_dx2]

% Clacul du "dx" a l’aide de fichier fonction « elimgauss.m »


[dx] = elimgauss(J,-F)
x=x+dx
end

7.4.3. Programme intégré dans MATLAB pour RSENL

1) Programme MATLAB principal (P_RSENL.m):


clear all; close all; clc
x0 = [0.2; 1.2];
[x] = fsolve(@ma_fonction,x0)

2) Programme MATLAB utilisant le fichier fonction:

v Fichier : (ma_fonction.m)

function F = ma_fonction(x)
F=[x(1)+2*x(2)-2;
x(1)^2+4*x(2)^2-4];

_______________________________________________________________________________
[Link] Boumediene -91-
Chapitre 8

Résolution d’Équations
Différentielles

[Link] Boumediene
Chapitre 8 Équations Différentielles

8.1. Introduction

Cette section traite de la façon dont on peut utiliser les méthodes de résolution
d'équations différentielles ordinaires dans le cas de systèmes d'équations
différentielles avec conditions initiales. Fort heureusement, il suffit d'adapter
légèrement les méthodes déjà vues.

La résolution numérique des équations différentielles est probablement le


domaine de l'analyse numérique ou les applications sont les plus nombreuses.
Que ce soit en électrotechnique, en modélisation des machines électriques ou
en analyse des réseaux électriques, on aboutit souvent à la résolution
d'équations différentielles, de systèmes d'équations différentielles ou plus
généralement d'équations aux dérivées partielles [31].

La solution numérique qui sera par la suite analysée et comparée à d'autres


solutions approximatives ou quasi analytiques. Parmi leurs avantages, les
méthodes numériques permettent d'étudier des problèmes complexes pour
lesquels on ne connait pas de solution analytique, mais qui sont d'tin grand
intérêt pratique [32].

Dans ce chapitre comme dans les précédents, les diverses méthodes de


résolution proposées sont d'autant plus précises qu'elles sont d'ordre élevé.

Nous amorçons l'expose par des méthodes relativement simples ayant une
interprétation géométrique. Elles nous conduiront progressivement à des
méthodes plus complexes telles la méthode de Runge-Kutta, qui permet
d'obtenir des résultats d'une grande précision. Nous considérons
principalement les équations différentielles avec conditions initiales, mais
nous ferons une brève incursion du coté des équations différentielles avec
conditions aux limites par le biais des méthodes d’Euler et de Taylor.

8.2. Définitions

La forme générale d'un système de m équations différentielles avec conditions


initiales s'écrit [31]:

y1¢ ( t ) = f1 ( t, y1 ( t ), y 2 ( t ),L , y m ( t )) ( y1 ( t 0 ) = y1,0 )


y¢2 ( t ) = f 2 ( t, y1 (t ), y 2 ( t ),L , y m (t )) ( y 2 ( t 0 ) = y 2 ,0 )
y¢3 (t ) = f 3 ( t, y1 ( t ), y 2 (t ),L , y m ( t )) ( y 3 ( t 0 ) = y 3,0 ) (8.1)
M M M
y ¢m ( t ) = f m ( t, y1 ( t ), y 2 ( t ),L , y m ( t )) ( y m ( t 0 ) = y m ,0 )

Ici encore, on note yi(tn), la valeur exacte de la iéme variable dépendante en


t = tn et yi,n, son approximation numérique.

Nous prenons comme point de départ la formulation générale d'une équation


différentielle d'ordre 1 avec condition initiale. La tache consiste a déterminer
une fonction y(t) solution de:

_______________________________________________________________________________
[Link] Boumediene -93-
Chapitre 8 Équations Différentielles

y¢( t ) = f ( t, y( t ))
(8.2)
y( t 0 ) = y0

La variable indépendante t représente très souvent (mais pas toujours) le


temps [31]. La variable dépendante est notée y et dépend bien sur de t. La
fonction f est pour le moment une fonction quelconque de deux variables que
nous supposons suffisamment différentiable. La condition y(t0) = y0 est la
condition initiale et en quelque sorte l'état de la solution au moment ou on
commence a s'y intéresser. Il s'agit d'obtenir y(t) pour t > t0, si on cherche une
solution analytique, ou une approximation de y(t), si on utilise une méthode
numérique.

8.3. Méthode de Runge-Kutta

Il serait avantageux de disposer de méthodes d'ordre de plus en plus élevé tout


en évitant les désavantages des méthodes de Taylor, qui nécessitent
l'évaluation des dérivées partielles de la fonction f(t, y). Une voie est tracée par
les méthodes de Runge-Kutta, qui sont calquées sur les méthodes de Taylor du
même ordre [31, 32].

On a vu que le développement de la méthode de Taylor passe par la relation


suivante :

y( t n +1 ) = y( t n ) + f ( t n + y( t n ))t
t 2 æ ¶f ( t n + y(t n )) ¶f ( t n + y( t n )) ö (8.3)
+ ç + f (t n + y( t n )) ÷ + C( t 3 )
2è ¶t ¶y ø

Le but est de remplacer cette dernière relation par une expression équivalente
possédant le même ordre de précision ( C(t 3 ). On propose la forme:

y( t n +1 ) = y( t n ) + f ( t n , y( t n ))a1 t + f (t n + a3 t, y(t n )a 4 t )a 2 t (8.4)

où on doit déterminer les paramètres ai, a1, a2, a3 et a4 de telle sorte que les
expressions (8.3) et (8.4) aient toutes deux une erreur en C( t 3 ) . On ne trouve
par ailleurs aucune dérivée partielle dans cette expression. Pour y arriver, on
doit recourir au développement de Taylor en deux variables. On a ainsi:

¶f (t n + y( t n ))
y( t n + a3 t, y( t n ) + a 4 t ) = f ( t n + y( t n )) + a3 t
¶t
(8.5)
¶f ( t n + y( t n ))
+ a4t + C( t 2 )
¶y

La relation (8.5) devient alors:


¶f ( t n + y( t n ))
y( t n +1 ) = y( t n ) + f ( t n , y( t n ))(a1 + a 2 )t + a 2a3 t 2
¶t
(8.6)
¶f ( t n + y( t n ))
+ a 2a 4 t 2 + C(t 3 )
¶y

_______________________________________________________________________________
[Link] Boumediene -94-
Chapitre 8 Équations Différentielles

On voit immédiatement que les expressions (8.3) et (8.6) sont du même ordre.
Pour déterminer les coefficients ai, il suffit de comparer ces deux expressions
terme à terme:

§ coefficients respectifs de f (t n , y( t n )) :

t = ( a1 + a 2 ) t (8.7)

¶f ( t n + y( t n ))
§ coefficients respectifs de :
¶t

t2
= a 2a 3 t 2 (8.8)
2
¶f ( t n + y( t n ))
§ coefficients respectifs de :
¶y

t2
f ( t n + y( t n )) = a 2a 4 t 2 (8.9)
2

On obtient ainsi un système non linéaire de 3 équations comprenant 4


inconnues:

1 = (a1 + a 2 )
1
= a 2a 3 (8.10)
2
f ( t n + y( t n ))
= a 2a 4
2

Le système (8.10) est sous-déterminé en ce sens qu'il y a moins d'équations que


d'inconnues et qu'il n'a donc pas de solution unique. Cela offre une marge de
manœuvre qui favorise la mise au point de plusieurs variantes de la méthode
de Runge-Kutta.

8.4. Exemple d’application sur la méthode de Runge-Kutta

Soit le système de deux équations différentielles suivant:

ìy1¢ (t ) = y 2 (t ) ( y1 ( 0 ) = 2)
í
îy¢2 ( t ) = 2y 2 ( t ) - y1 ( t ) ( y 2 (0) = 1)

Dont la solution analytique est:

ìïy1 ( t ) = 2et - te t
í
ïîy 2 (t ) = e t - te t

On a alors:

_______________________________________________________________________________
[Link] Boumediene -95-
Chapitre 8 Équations Différentielles

ìf1 (t , y1 ( t ), y 2 ( t )) = y 2 ( t )
í
îf2 ( t, y1 ( t ), y 2 ( t )) = 2y 2 ( t ) - y1 ( t )

et la condition initiale (t0, y1,0, y2,0) = ( 0 , 2 , 1 ) . Si on prend par exemple


t = 0,1, on trouve:

k1,1 = 0.1( f1 (0, 2, 1)) = 0.1


k 2,1 = 0.1(f 2 (0, 2, 1)) = 0

k1,2 = 0.1( f1 (0.05, 2.05, 1.0)) = 0.1


k 2,2 = 0.1( f 2 (0.05, 2.05, 1.0)) = -0.005

k1,3 = 0.1( f1 (0.05, 2.05, 0.9975)) = 0.09975


k 2,3 = 0.1( f 2 (0.05, 2.05, 0.9975)) = -0.0055

k1,4 = 0.1( f1 (0.1, 2.09975, 0.9975)) = 0.09945


k 2,4 = 0.1( f 2 (0.1, 2.09975, 0.9975)) = -0.011075

1
y1,1 = y1,0 + (0.1 + 2(0.1) + 2(0.09975) + 0.09945)
6
= 2.099825
1
y 2,1 = y 2,0 + (0 + 2( -0.005) + 2( -0.0055) + ( -0.011075))
6
= 0.994651

8.5. Méthode d’Euler

La méthode d'Euler est de loin la méthode la plus simple de résolution


numérique d'équations différentielles. Elle possède une belle interprétation
géométrique et son emploi est facile. Toutefois, elle est relativement peu
utilisée en raison de sa faible précision [31, 33].

Reprenons l'équation différentielle 8.2 et considérons plus attentivement la


condition initiale y(t0)=y0. Le but est maintenant d'obtenir une approximation
de la solution en t = t1 = t0 + h. Avant d'effectuer la première itération, il faut
déterminer dans quelle direction on doit avancer a partir du point (t0,y0) pour
obtenir le point (t1,y1), qui est une approximation du point (t1, y(t1)). En effet,
l'équation différentielle assure que:

y¢( t 0 ) = f ( t 0 , y( t 0 )) = y( t 0 , y 0 ) (8.11)

Il est important de noter que, le plus souvent, y1=y(t1). Cette inégalité n'a rien
pour étonner, mais elle a des conséquences sur la suite du raisonnement. En
effet, si on souhaite faire une deuxième itération et obtenir une approximation
de y(t2), on peut refaire l'analyse précédente à partir du point (t1,y1). On
remarque cependant que la solution analytique en t= t1 est:

_______________________________________________________________________________
[Link] Boumediene -96-
Chapitre 8 Équations Différentielles

y¢( t1 ) = f ( t1 , y( t1 )) » f ( t1 , y1 ) (8.12)

Le développement qui précède met en évidence une propriété importante des


méthodes numériques de résolution des équations différentielles.

En effet, l’erreur introduite à la première itération a des répercussions sur les


calculs de la deuxième itération, ce qui signifie que les erreurs se propagent
d'une itération à l'autre. Il en résulte de façon générale que l'erreur:

y(t i ) = y i (8.13)

augmente légèrement avec i.

8.6. Exemple d’application sur la méthode d'Euler

Soit l'équation différentielle :

y¢( t ) = -y( t ) + t + 1

et la condition initiale y(0) = 1.

On a donc t0 = 0 et y0 = 1, et on prend un pas de temps t = 0,1. De plus, on a:

f ( t, y ) = - y + t + 1

On peut donc utiliser la méthode d'Euler et obtenir successivement des


approximations de y(0,l), y(0,2), y(0,3) …, notées y1, y2, y3… La première
itération produit:

y1 = y 0 + t × f ( t 0 , y 0 ) = 1 + 0.1 × f (0,1) = 1 + 0.1(- 1 + 0 + 1) = 1

La deuxième itération fonctionne de manière similaire:

y2 = y1 + t × f ( t1 , y1 ) = 1 + 0.1 × f (0.1,1) = 1 + 0.1(- 1 + 0.1 + 1) = 1.01

On parvient à:

y3 = y2 + t × f ( t 2 , y2 ) = 1.01 + 0.1 × f (0.2,1.01) = 1.01 + 0.1(- 1.01 + 0.2 + 1) = 1.029

Le tableau suivant rassemble les résultats des dix premières itérations. On


peut montrer que la solution analytique de cette équation différentielle est:

y(t ) = e - t + t

ce qui permet de comparer les solutions numérique et analytique, et de


constater la croissance de l'erreur.

_______________________________________________________________________________
[Link] Boumediene -97-
Chapitre 8 Équations Différentielles

ti y( t i ) yi y(t i ) = y i

0,0 1,000000 1,000000 0,000000


0,1 1,004837 1,000000 0,004837
0,2 1,018731 1,010000 0,008731
0,3 1,040818 1,029000 0,011818
0,4 1,070302 1,056100 0,014220
0,5 1,106531 1,090490 0,016041
0,6 1,148812 1,131441 0,017371
0,7 1,196585 1,178297 0,018288
0,8 1,249329 1,230467 0,018862
0,9 1,306570 1,287420 0,019150
1,0 1,367879 1,348678 0,019201

8.7. Méthodes de Taylor

Le développement de Taylor autorise une généralisation immédiate de la


méthode d'Euler, qui permet d'obtenir des algorithmes dont l'erreur de
troncature locale est d'ordre plus élevé. Nous nous limitons cependant à la
méthode de Taylor du second ordre [31, 34].

On cherche, au temps t = tn, une approximation de la solution en t = tn+1. On a


immédiatement:

y( tn +1 ) = y( t n + t )
y' ' ( t n )t 2 (8.14)
= y( t n ) + y' (t n )t +
2
En se servant de l'équation différentielle (8.2), on trouve:

f ' ( t n + y(t n ))t 2


y( t n +1 ) = y( t n ) + f ( t n + y(t n ))t + (8.15)
2
Dans la relation précédente, on voit apparaitre la dérivée de la fonction f(t,y(t))
par rapport au temps. La règle de dérivation en chaine assure que :

¶f ( t + y( t )) ¶f ( t + y( t ))
f ' ( t + y( t )) = + y' ( t ) (8.16)
¶t ¶y
c'est-a-dire:

¶f ( t + y(t )) ¶f ( t + y(t ))
f ' ( t + y( t )) = + f ( t + y( t )) (8.17)
¶t ¶y

On obtient donc:

_______________________________________________________________________________
[Link] Boumediene -98-
Chapitre 8 Équations Différentielles

y( t n +1 ) = y( t n ) + f ( t n + y( t n ))t
t 2 æ ¶f ( t n + y(t n )) ¶f (t n + y( t n )) ö (8.18)
+ ç + f ( t n + y( t n )) ÷
2è ¶t ¶y ø

8.8. Exemple d’application sur la méthode de Taylor

Soit l'équation différentielle déjà résolue par la méthode d'Euler:

y¢( t ) = -y( t ) + t + 1

et la condition initiale y(0) = 1. Dans ce cas:

f ( t, y ) = - y + t + 1

et

¶f ¶f
= 1 et = -1
¶t ¶y

L'algorithme devient:

t2
y n +1 = y n + ( -y n + t n + 1)t + (1 + ( -1)( -y n + t n + 1))
2

La première itération de la méthode de Taylor d'ordre 2 donne (avec t = 0,1):

( 0 .1 ) 2
y1 = 1 + 0.1( -1 + 0 + 1) + (1 + ( -1)( -1 + 0 + 1)) = 1.005
2

Une deuxième itération donne:

(0.1)2
y 2 = 1.005 + 0.1( -1.005 + 0.1 + 1) + (1 + ( -1)( -1.005 + 0.1 + 1)) = 1.019025
2

Les résultats sont compiles dans le tableau qui suit.

ti y( t i ) yi y(t i ) = y i

0,0 1,000000 1,000000 0,000000


0,1 1,004837 1,005000 0,000163
0,2 1,018731 1,019025 0,000294
0,3 1,040818 1,041218 0,000400
0,4 1,070302 1,070802 0,000482
0,5 1,106531 1,107075 0,000544
0,6 1,148812 1,149404 0,000592
0,7 1,196585 1,197210 0,000625
0,8 1,249329 1,249975 0,000646
0,9 1,306570 1,307228 0,000658
1,0 1,367879 1,368541 0,000662

_______________________________________________________________________________
[Link] Boumediene -99-
Chapitre 8 Équations Différentielles

On remarque que l’erreur est plus petite avec la méthode de Taylor d'ordre 2
qu'avec la méthode d'Euler. Comme on le verra plus loin, cet avantage des
méthodes d'ordre plus élevé vaut pour l'ensemble des méthodes de résolution
d'équations différentielles.

8.9. Applications MATLAB sur la résolution d’équations Différentielles

8.9.1. Application MATLAB de la méthode de Runge-Kutta

[Link]. Algorithme de la méthode de Runge-Kutta

Etape 1 Données : un pas de temps t, des conditions initiales (t0, y1,0, y2,0,
… , ym,0) et un nombre maximal d'itérations N.

Etape 2 Pour 0 ≤ n ≤ N :
k 1 = f ( t n , y n )t
t k
k 2 = f (t n + , y n 1 )t
2 2
t k
k 3 = f ( t n + , y n 2 )t
2 2
k 4 = f ( t n + t, y n k 3 )t
1
y n +1 = y n + ( k 1 + 2k 2 + 2k 3 + k 4 )
6

Etape 3 Poser : t n +1 = t n + t

Etape 4 Ecrire : t n +1 et y n +1

Etape 5 Afficher le résultat.

[Link]. Programme MATLAB de la méthode de Runge-Kutta

function [xSol,ySol] = runKutta(dEqs,x,y,xStop,h)

% Runge-Kutta 4ème ordre.


% USAGE: [xSol,ySol] = runKutta(dEqs,x,y,xStop,h)

% ENTREES:
% dEqs = équations différentielles
% F(x,y) = [dy1/dx dy2/dx dy3/dx ...].
% x,y = Valeurs initiales; y Vecteur ligne.
% xStop = Valeur finale de x.
% h = incrément de x utilisé.

% SORTIES:
% xSol = valeur de x où la solution est complète.
% ySol = valeur de y correspondant à la valeur x.

_______________________________________________________________________________
[Link] Boumediene -100-
Chapitre 8 Équations Différentielles

if size(y,1) > 1 ; y = y'; end % y vecteur ligne


xSol = zeros(2,1); ySol = zeros(2,length(y));
xSol(1) = x; ySol(1,:) = y;
i = 1;
while x < xStop
i = i + 1;
h = min(h,xStop - x);
K1 = h*feval(dEqs,x,y);
K2 = h*feval(dEqs,x + h/2,y + K1/2);
K3 = h*feval(dEqs,x + h/2,y + K2/2);
K4 = h*feval(dEqs,x+h,y + K3);
y = y + (K1 + 2*K2 + 2*K3 + K4)/6;
x = x + h;
xSol(i) = x; ySol(i,:) = y;
end

8.9.2. Application MATLAB de la méthode d'Euler

[Link]. Algorithme de la méthode d'Euler

Etape 1 Données : un pas de temps t, des conditions initiales (t0,y0)


et un nombre maximal d'itérations N.

Etape 2 Pour 0 ≤ n ≤ N : y n +i = y n + f(t n , y n )t

Etape 3 Poser : t n +1 = t n + t

Etape 4 Ecrire : t n +1 et y n +1

Etape 5 Afficher le résultat.

[Link]. Programme MATLAB de la méthode d'Euler

function [t,y] = Euler(f,tspan,y0,N)

% La méthode de Euler pour la résolution d’équation


% différentielles y’(t) = f(t,y(t))
% pour tspan = [t0,tf] et la valeur initiale y0 et N
% pas du temps
if nargin<4 | N <= 0, N = 100; end
if nargin<3, y0 = 0; end

%valeur du temps t
h = (tspan(2) - tspan(1))/N;
% vecteur du temps
t = tspan(1)+[0:N]’*h;

% toujours faites la valeur initiale un vecteur de la ligne


y(1,:) = y0(:)’;
for k = 1:N
y(k + 1,:) = y(k,:) +h*feval(f,t(k),y(k,:)); end

_______________________________________________________________________________
[Link] Boumediene -101-
Chapitre 8 Équations Différentielles

8.9.3. Application MATLAB de la méthode de Taylor

[Link]. Algorithme de la méthode de Taylor

Etape 1 Données : un pas de temps t, des conditions initiales (t0,y0)


et un nombre maximal d'itérations N.

Etape 2 Pour 0 ≤ n ≤ N :

y n +1 = y n + f (t n + y n )t
t 2 æ ¶f ( t n + y n ) ¶f ( t n + y n ) ö
+ ç + f (tn + y n ) ÷
2è ¶t ¶y ø

Etape 3 Poser : t n +1 = t n + t

Etape 4 Ecrire : t n +1 et y n +1

Etape 5 Afficher le résultat.

[Link]. Programme MATLAB de la méthode de Taylor

function [xSol,ySol] = taylor(deriv,x,y,xStop,h)

% La méthode de Taylor pour la résolution d’équation


% différentielles y’(t) = f(t,y(t))
% pour tspan = [t0,tf] et la valeur initiale y0 et N
% pas du temps
% la fonction : [xSol,ySol] = taylor(deriv,x,y,xStop,h)

% ENTREES:
% deriv = fonction qui rend la matrice
% d = [dy/dx d^2y/dx^2 d^3y/dx^3 d^4y/dx^4].
% x,y =valeur initiale; il faut que y soit vecteur linge.
% xStop = la valeur finale de x
% h = increment de x ulilisé dans le calcul (h > 0).

% SORTIES:
% xSol = valeurs de x à la solution qui calculée.
% ySol = valeurs de y correspond aux valeurs de x.

if size(y,1) > 1; y = y'; end % y soit vecteur linge


xSol = zeros(2,1); ySol = zeros(2,length(y));
xSol(1) = x; ySol(1,:) = y;
k = 1;
while x < xStop
h = min(h,xStop - x);
d = feval(deriv,x,y); % dérivée de y
hh = 1;
for j = 1:4 % calcul Taylor
hh = hh*h/j;

_______________________________________________________________________________
[Link] Boumediene -102-
Chapitre 8 Équations Différentielles

y = y + d(j,:)*hh;
end
x = x + h; k = k + 1;
xSol(k) = x; ySol(k,:) = y; % affichage de résultats.
end

_______________________________________________________________________________
[Link] Boumediene -103-
Conclusion Générale

[Link] Boumediene
Conclusion

Conclusion
Le calcul scientifique est une discipline qui consiste à développer, analyser et
appliquer des méthodes relevant de domaines mathématiques aussi variés que
l’analyse, l’algèbre linéaire, la géométrie, la théorie de l’approximation, les
équations fonctionnelles, l’optimisation ou le calcul différentiel. Les méthodes
numériques trouvent des applications naturelles dans de nombreux problèmes
posés par la physique, les sciences biologiques, les sciences de l’ingénieur,
l’économie et la finance.

Le calcul scientifique se trouve donc au carrefour de nombreuses disciplines


des sciences appliquées modernes, auxquelles il peut fournir de puissants
outils d’analyse, aussi bien qualitative que quantitative. Ce rôle est renforcé
par l’évolution permanente des ordinateurs et des algorithmes : la taille des
problèmes que l’on sait résoudre aujourd’hui est telle qu’il devient possible
d’envisager la simulation de phénomènes réels.

La communauté scientifique bénéficie largement de la grande diffusion des


logiciels de calcul numérique. Néanmoins, les utilisateurs doivent toujours
choisir avec soin les méthodes les mieux adaptées à leurs cas particuliers : il
n’existe en effet aucune “boite noire” qui puisse résoudre avec précision tous
les types de problème.

Un des objectifs de ce document est de présenter les fondements


mathématiques du calcul scientifique, avec la mise en œuvre pratique est
proposée dans le langage MATLAB qui présente l’avantage d’être d’une
utilisation aisée et de bénéficier d’une large diffusion.

_______________________________________________________________________________
[Link] Boumediene -104-
Références

[Link] Boumediene
Références

Références

[1] F. Gustafsson & N. Bergman, “MATLAB for engineers explained”;


Springer-Verlag, 2003.

[2] M. Mokhtari & A. Mesbah, “Apprendre à maîtriser MATLAB”; Springer-


Verlag, 1997.

[3] K. Sigmon, “MATLAB Aide-Mémoire”; Springer-Verlag, 1999.

[4] D. Salomon, “Curves and Surfaces for Computer Graphics”; Springer-


Verlag, 2006.

[5] W. Gander & J. Hˇrebicˇek, “Solving problems in scientific computing


using Maple and MATLAB, Springer-Verlag, 4ed, 2004.

[6] A. Biran & M. Breiner, “MATLAB for engineers”; éditions Wesley, 1999.

[7] S. R. Otto & J.P. Denier; “An Introduction to Programming and


Numerical Methods in MATLAB”; Springer-Verlag London, 2005.

[8] K. Chen, P. Giblin & A. Irving, “Mathematical explorations with


MATLAB”; Cambrige University Press, 1999.

[9] H. B. Wilson, L. H. Turcotte & D. Halpern, “Advanced mathematics and


mechanics applications using MATLAB”; Chapman and Hall, 2003.

[10] J. Kiusalaas; “Numerical Methods in Engineering with MATLAB”;


Cambridge University, 2005.

[11] W. Y. Yang, W. Cao, T. Chung & J. Morris; “Applied numerical methods


using MATLAB”; John Wiley publication, 2005.

[12] S. T. Karris; “Numerical Analysis Using MATLAB and Excel”; Third


Edition, Orchard Publications, 2007.

[13] J. L. Merrien, “Analyse Numérique Avec MATLAB”; Dunod, Paris, 2007.

[14] M. Sibony, “Ananlyse numérique III : Itérations et approximations”;


éditions Herman, 1988.

[15] A. Fortin, “Analyse numérique pour ingénieurs”; l’École Polytechnique de


Montréal, 1995.

[16] M. Sibony & J. Cl. Mardon, “Ananlyse numérique I : Systèmes linéaires


et non linéaires”; éditions Herman, 1982.

[17] P. Lascaux & R. Theodor, “Ananlyse numérique matricielle appliquée à


l'art de l'ingénieur : Méthodes directes”; Tome 1 ; éditions Masson, 2ème
édition, 1993.

_______________________________________________________________________________
[Link] Boumediene -105-
Références

[18] J. Cooper, “MATLAB companion for multivariable calculus (a)”;


Academic Press, 2001.

[19] C. B. Moler, “Numerical computing with MATLAB”; SIAM, 2004.

[20] A. Quarteroni, R. Sacco & F. Saleri, “Méthodes numériques pour le calcul


scientifique : Programmes en MATLAB”; Springer-Verlag, 2006.

[21] J. Bastien, “Introduction à l’analyse numérique : Applications sous


MATLAB”; éditions Dunod, 2003.

[22] A. Kharab & R. B. Guenther, “Introduction to numerical methods (a): a


MATLAB approach”; Chapman and Hall, 2002.

[23] J. Baranger, “Introduction à l'analyse numérique”; éditions Herman,


1993.

[24] F. Jerdrzejewski, “Introduction aux méthodes numériques”; Springer-


Verlag, 2001.

[25] E. Süli & D. Mayers, “An Introduction to Numerical Analysis”;


Cambridge University, 2003.

[26] A. Ralston & P. Rabinowitz, “A first course in Numerical Analysis”;


éditions Presses Universitaires de Grenoble, 1991.

[27] A. Quarteroni, R. Sacco & F. Saleri, “Méthodes Numériques Algorithmes,


analyse et applications”; Springer-Verlag Italia, Milano 2007.

[28] L. Jolivet & R. Labbas, “Analyse et analyse numérique : rappel de cours


et exercices corrigés”; Hermès Science, 2005.

[29] A. Gourdin & M. Boumahrat, “Méthodes numériques appliquées”;


éditions Techniques et Documentations Lavoisier, 1989.

[30] J.P. Nougier, “Méthodes de calcul numérique”; éditions Masson, 3ème


édition, 1993.

[31] J. P. Demailly, “Analyse Numérique et Equations Différentielles”;


Collection Grenoble Sciences, 2006.

[32] P. Lascaux & R. Theodor, “Ananlyse numérique matricielle appliquée à


l'art de l'ingénieur : Méthodes directes”; Tome 2 ; éditions Masson, 2ème
édition, 1994.

[33] P. Deuflhard & A. Hohmann, “Numerical Analysis in Modern Scientific


Computing, An Introduction”;Springer-Verlag, 2003, 2e édition.

[34] J.-G. Dion & R. Gaudet, “Méthodes d’analyse numérique, de la théorie à


l’application”; Modulo, 1996.

_______________________________________________________________________________
[Link] Boumediene -106-
Annexe

- Exercices Chapitre 2 -
Exercice 1 :

Soit la fonction f(x) = x3+4x2-10, avec x Î [1, 2], ξ=10-3 ; trouvez une valeur
approchée de la solution exacte de f(x)=0 par la méthode de Dichotomie.

Solutions de l’Exercice 1:

f(x) est continue sur [1, 2] et f(1)·f(2) = (-5)·(14)< 0, donc d’après la méthode de
Dichotomie il existe au moins un x*Î [1, 2] tel que f(x*)=0.

On construit une suite x0 = (a+b)/2, x1 = (a1 + b1)/2, . . . , xn = (an + bn)/2 ; et en


vérifier que :

• si f(a)f(x0) < 0, alors α Î ]a, x0[. On pose a1 = a, b1 = x0.

• si f(a)f(x0) = 0, alors α = x0.

• si f(a)f(x0) > 0, alors α Î ]x0, b[. On pose a1 = x0, b1 = b.

n an bn xn f(xn)
0 1 2 1.5 +2.375
1 1 1.5 1.25 -1.79687
2 1.25 1.5 1.375 +0.16211
3 1.25 1.375 1.3125 -0.84839
4 1.3125 1.375 1.34375 -0.35098
5 1.34375 1.375 1.35937 -0.09641
6 1.35937 1.375 1.36718 +0.03236
7 1.35937 1.36718 1.36328 -0.03215
8 1.36328 1.36718 1.36523 0.000072

Test d’arrêt :
ξ=10-3 c’est-à-dire |x8 – x7| ≤ ξ ; Alors x* = x8 = 1.36523

Exercice 2 :

Soit la fonction f(x) = x3+4x2-10, avec x Î [1, 2], ξ=10-3 ; trouvez une valeur
approchée de la solution exacte de f(x)=0 par la méthode du point fixe.

Solutions de l’Exercice 2:

L’équation f(x) =0 peut être transformée en une équation équivalente :


g(x) = x, où :
f(x) = x3+4x2-10 = 0 → x2(x+4) = 10
1) existe x*.

_______________________________________________________________________________
[Link] Boumediene -105-
Annexe

1
æ 10 ö 2
2) g(x) = çç ÷÷
è ( x + 4) ø
3) l’application de la méthode du point fixe, formons la suite xn+1=g(xn) ;
on particulier x0 = 1.5, pour ce choix en à :

1er itération x1 = g(x0) = 1.3483


2ème itération x2 = g(x1) = 1.3673
3ème itération x3 = g(x2) = 1.3649
4ème itération x4 = g(x3) = 1.3652
5ème itération x5 = g(x4) = 1.3652

4) Test d’arrêt :
ξ=10-3 c’est-à-dire |x5 – x4| ≤ ξ ; Alors x* = x5 = 1.3652

Exercice 3 :

Soit la fonction f(x) = x3+4x2-10, avec x Î [1, 2], ξ=10-3 ; trouvez une valeur
approchée de la solution exacte de f(x)=0 par la méthode de Newton.

Solutions de l’Exercice 3:

L’équation f(x) = x3+4x2-10 et leur dérivé est f’(x) = 3x2+8, alors on part de x0 =
1.5 et en appliquant la formule de Newton :
xn+1 = xn − f(xn)/f′(xn)
On obtient :
xn+1 = xn – (xn3+4xn2-10) /(3xn2+8)

Nombre d’itérations est de n = 1, 2, 3… jusqu’au le test d’arrêt :


pour n = 1 : x1 = x0 – (x03+4x02-10) /(3x02+8) = 1.37333
pour n = 2 : x2 = x1 – (x13+4x12-10) /(3x12+8)
pour n = 3 : x3 = x2 – (x23+4x22-10) /(3x22+8)

- Exercice Chapitre 3 -

Exercice 1 :

Utilisant l’élimination de Gauss résoudre le système d’équations linéaires


suivant:

ì2 x1 + x2 + 4 x4 = 2
ï - 4 x - 2 x + 3 x - 7 x = -9
ï 1 2 3 4
í
ï4 x1 + x2 - 2 x3 + 8x4 = 2
ïî- 3x2 - 12 x3 - x4 = 2

_______________________________________________________________________________
[Link] Boumediene -106-
Annexe

Solutions de l’Exercice 1:
ì2x 1 + x 2 + 4x 4 = 2
ï
ï- 4 x1 - 2x 2 + 3x 3 - 7x 4 = -9
Soit à résoudre le système linéaire : í
ï4 x1 + x 2 - 2x 3 + 8x 4 = 2
ï- 3x 2 - 12x 3 - x 4 = 2
î

é 2 1 0 4ù éx 1 ù é2 ù
ê- 4 - 2 3 - 7úú êx ú ê- 9 ú
1) En posant : A = ê , X = ê 2 ú et B = ê ú
ê4 1 -2 8 ú êx 3 ú ê2 ú
ê ú ê ú ê ú
ë 0 - 3 - 12 - 1 û ëx 4 û ë2 û

Le système linéaire s’écrit sous la forme matricielle suivante A∙X = B :

é2 1 0 4 ù é x1 ù é 2 ù
ê- 4 - 2 3
ê - 7úú êêx 2 úú êê- 9úú
× =
ê4 1 - 2 8 ú êx 3 ú ê 2 ú
ê ú ê ú ê ú
ë 0 - 3 - 12 - 1 û ëx 4 û ë 2 û

2) Résolution du système par la méthode de Gauss : Partons de la matrice


augmentée :

é ù
ê2 1 0 4 2 ú
ê ú
ê- 4 - 2 3 -7 - 9ú
ê4 1 -2 8 2 ú
ê ú
ê1 0 - 3 - 12 - 1 {2 ú
44424443
êë A B ú
û

L’élément générateur a11 = 2 ≠ 0. Remplaçons la 2ème ligne par :


ligne 2 − (− 2)× la ligne 1, en trouve :

é2 1 0 4 2ù
ê0 0 ú
ê 3 1 - 5ú
ê4 1 -2 8 2ú
ê ú
ë0 - 3 - 12 - 1 2 û

Remplaçons la troisième ligne par : ligne 3 − 2 × ligne 1, en aperçoit :

é2 1 0 4 2 ù
ê ú
ê0 0 3 1 - 5ú
ê0 - 1 - 2 0 - 2ú
ê ú
ë0 - 3 - 12 - 1 2 û

Et puisque l’élément générateur a22 est nul, faisons une permutation entre la
deuxième ligne et la troisième ligne :

_______________________________________________________________________________
[Link] Boumediene -107-
Annexe

é2 1 0 4 2 ù
ê ú
ê0 - 1 - 2 0 - 2ú
ê0 0 3 1 - 5ú
ê ú
ë0 - 3 - 12 - 1 2 û

On remplace la quatrième ligne par : ligne 4 − 3 × ligne 2, en trouve :

é2 1 0 4 2 ù
ê ú
ê0 - 1 - 2 0 - 2ú
ê0 0 3 1 - 5ú
ê ú
ë0 0 - 6 - 1 8 û

À l’étape suivante, nous remplaçons la 4ème ligne par ligne 4 − (− 2) × ligne 3 :

é2 1 0 4 2 ù
ê ú
ê0 - 1 - 2 0 - 2ú
ê0 0 3 1 - 5ú
ê ú
ë0 0 0 1 - 2û

Nous avons abouti au système triangulaire supérieur suivant :

é2 1 0 4ù é x1 ù é 2 ù
ê0 - 1 - 2
ê 0 úú êêx 2 úú êê- 2úú
× =
ê0 0 3 1 ú ê x 3 ú ê - 5ú
ê ú ê ú ê ú
ë0 0 0 1 û ëx 4 û ë- 2û

D’où on tire les équations :

ì x 4 = -2
ï 3x + x = -5
ï 3 4
í
ï- x 2 + 2x 3 = -2
ïî2x 1 + x 2 + 4x 4 = 2

La solution du système est :

x1 = 3 x2 = 4 x 3 = -1 x 4 = -2

Exercice 2 :

Utilisant l’élimination de Gauss résoudre le système d’équations linéaires


suivant:

é 6 -4 1 ù é- 14 22 ù
A = êê- 4 6 - 4úú B = ê 36 - 18ú
ê ú
êë 1 - 4 6 úû êë 6 7 úû

_______________________________________________________________________________
[Link] Boumediene -108-
Annexe

Solutions de l’Exercice 2:

Soit à résoudre le système linéaire :

é 6 -4 1 ù é - 14 22 ù
A = ê - 4 6 - 4ú , B = êê 36 - 18úú
ê ú
êë 1 - 4 6 úû êë 6 7 úû

La matrice augmentée est :

é 6 - 4 1 - 14 22 ù
ê ú
ê- 4 6 - 4 36 - 18ú
ê 1 -4 6 6 7 úû
ë

Elimination 1 :

{ligne(2)-(-2/3).ligne(1)} → Nouvelle ligne(2)


{ligne(3)-(1/6).ligne(1)} → Nouvelle ligne(3)

é6 -4 1 - 14 22 ù
ê ú
ê0 10 / 3 - 10 / 3 80 / 3 - 10 / 3ú
ê0 - 10 / 3 35 / 6 25 / 3 10 / 3 ú
ë û

Elimination 2 :

{ligne(3)-(-1).ligne(2)} → Nouvelle ligne(3)

é6 - 4 1 - 14 22 ù
ê ú
ê0 10 / 3 - 10 / 3 80 / 3 - 10 / 3ú
ê0 0 5/ 2 35 0 úû
ë

La solution du système est :

x11 = 10 x 21 = 22 x 31 = 14 Et x12 = 3 x 22 = -1 x 32 = 0

Exercice 3 :

Utilisant l’élimination de Gauss résoudre le système d’équations linéaires


suivant:

é2 1 2 10 ù
ê ú
C = ê6 4 0 26ú
ê8 5 1 35ú
ë û

_______________________________________________________________________________
[Link] Boumediene -109-
Annexe

Solutions de l’Exercice 3:

1) Former la matrice Triangulaire Supérieur utilisant d’Elimination de Gauss :

é2 1 2 10 ù é2 1 2 10 ù é2 1 2 10 ù
ê ú ê ú ê ú
ê6 4 0 26ú → ê0 1 - 6 - 4 ú → ê0 1 - 6 - 4ú
ê8 5 1 35ú ê0 1 - 7 - 5 ú ê0 0 - 1 - 1 úû
ë û ë û ë

1) Calculer les éléments du vecteur X :

10 - (1)(2) - (2)(1) - 4 - ( -6)(1) -1


x1 = = 3 ; x2 = = 2 ; x3 = = 1.
2 1 -1

Exercice 4 :

Par la méthode de Gauss-Seidel résoudre les systèmes d’équations suivants:

ì 4x 1 - x 2 + x 3 = 7
ï
í- 4 x1 + 8x 2 - x 3 = 21
ï - 2x + x + 5x = 15
î 1 2 3

Solutions de l’Exercice 4 :

ì 4x 1 - x 2 + x 3 = 7
ï
Soit à résoudre le système linéaire : í- 4x 1 + 8x 2 - x 3 = 21
ï - 2x + x + 5x = 15
î 1 2 3

é1 ù
ê ú 0
Avec la solution estimée : x = ê2ú et « ε » erreur 10-4.
êë2úû

ì x 1 = (7 - ( -x 2 + x 3 )) / 4
ï
íx 2 = (21 - ( -4x 1 - x 3 )) / 8
ï x = (15 - ( -2x + x )) / 5
î 3 1 2

é1 ù é1.7500 ù é1.9500 ù é1.9956 ù é1.9993 ù


x 0 = êê2úú → x1 = êê3.7500 úú → x 2 = êê3.9688 úú → x 3 = êê3.9961 úú → x 4 = êê3.9995 úú
êë2úû êë2.9500 úû êë2.9863úû êë2.9990 úû êë2.9998 úû

é1.9999 ù é2.0000 ù é2.0000 ù


5 ê ú 6 ê ú 7 ê ú
→ x = ê3.9999 ú → x = ê4.0000 ú → x = ê4.0000 ú
êë3.0000 úû êë3.0000 úû êë3.0000 úû

_______________________________________________________________________________
[Link] Boumediene -110-
Annexe

- Exercices Chapitre 4 -

Exercice 1 :

On doit calculer le polynôme passant par les 4points (0, 1), (1, 2), (2, 9)
et (3, 28) utilisant la Matrice de Vandermonde.

Solutions de l’Exercice 1:

On doit calculer le polynôme figure 4.1 passant par les points (0, 1), (1, 2), (2,
9) et (3, 28). Etant donne ces 4 points, les coefficients a i sont les solutions de :
é ù
ê1 x 0 x 20 x 30 L x n0 ú éa 0 ù é f ( x 0 )ù é ù
ê nú ê
ê 1 0 0 0 ú éa 0 ù é 1 ù
ê1 x 1 x 2
1 x 3
1 L x a ú êê f ( x1 ) úú
1 ú ê 1ú ê ú ê ú ê ú
ê1 x x 2 x 3 L x n ú × êa 2 ú = ê f ( x 2 )ú º ê 1 1 1 1 ú × ê a1 ú = ê 2 ú
ê 2 2 2 2
ú ê ú ê ú ê 1 2 4 8 ú êa 2 ú ê 9 ú
êM M M M O M ú êM ú ê M ú ê ú ê ú ê ú
ê1 x ê 114 3 9 27 ú ëa 4 û ë28û
x 2n x 3n L x nn ú êëa n úû êëf ( x n )úû 4244 3
êë Matrice de Vandermondeúû
ê 14n444 244443 ú
ë Matrice de Vandermonde û

dont la solution (par Elimination de Gauss (chapitre 3)) est a = [1 0 0 1]T.


Le polynôme recherche est donc : P3 ( x ) = 1 + x 3 .

Exercice 2 :

On doit calculer le polynôme passant par les 4points (0, 1), (1, 2), (2, 9)
et (3, 28) utilisant l’interpolation de Lagrange.

Solutions de l’Exercice 2:

Le polynôme passant par les points (0, 1), (1, 2), (2, 9) et (3, 28) est :

( x - 1)(x - 2)(x - 3) ( x - 0)(x - 2)(x - 3)


P3 ( x ) = 1 +2
(0 - 1)(0 - 2)(0 - 3) (1 - 0)(1 - 2)(1 - 3)

( x - 0)(x - 1)(x - 3) ( x - 0)(x - 1)(x - 2)


+9 + 28 = 1 + x3
(2 - 0)(2 - 1)(2 - 3) (3 - 0)(3 - 1)(3 - 2)

Exercice 3 :

On doit calculer le polynôme passant par les 4points (0, 1), (1, 2), (2, 9)
et (3, 28) utilisant le polynôme de Newton.

Solutions de l’Exercice 3:

Le polynôme de Newton passant par les points (0, 1), (1, 2), (2, 9) et (3, 28) :

_______________________________________________________________________________
[Link] Boumediene -111-
Annexe

xi f (x i ) f [ x i , x i +1 ] f [ x i , x i +1 , x i +2 ] f [ x i , x i +1 , x i +2 , x i +3 ]

0 1
1
1 2 3
7
1
2 9
6
19
3 28

p 3 ( x ) = 1 + 1 ( x - 0 ) + 3 ( x - 0 )( x - 1 ) + 1 ( x - 0 )( x - 1 )( x - 2 ) = x 3 + 1

Exercice 4 :

On doit calculer le polynôme passant par les 4 points (2, 1), (0, -1), (5, 10) et
(3, -4) utilisant le polynôme de Newton.

Solutions de l’Exercice 4 :

Le polynôme de Newton passant par les points (2, 1), (0, -1), (5, 10) et (3, -4) :

xi f (x i ) f [ x i , x i +1 ] f [ x i , x i +1 , x i +2 ] f [ x i , x i +1 , x i +2 , x i +3 ]

2 1
1
0 -1 0.4
2.2
1.2
5 10
1.6
7
3 -4

p 3 ( x ) = 1 + 1 ( x - 2 ) + 0 . 4 ( x - 2 )( x - 0 ) + 1 . 2 ( x - 2 )( x - 0 )( x - 5 )
3
= 1.2x - 8x 2 + 2.2x - 1

- Exercices Chapitre 5 -
Exercice 1 :

Soit f ( x ) une fonction donnée par le tableau :

xi 0.25 0.50 0.75 1.00 1.25


f (x i ) 0.630 0.794 0.909 1.000 1.15

1.25

Trouver une valeur approchée de l’intégrale I = ò f (x ) dx ,


0.25
en utilisant

la formule de Newton-Côtes.

_______________________________________________________________________________
[Link] Boumediene -112-
Annexe

Solutions de l’Exercice 1:

Par la formule de Newton-Côtes et puisque le pas d’intégration est h = 0.25 .


De la relation :

b-a 1.25 - 0.25


n= Û n= =4
h 0.25

Les coefficients de Côtes peuvent être calculés de la formule, ou peuvent être


pris directement du tableau 5.1 des coefficients de Côtes. Les valeurs des
coefficients où n=4 sont :

7 32 12
H0 = H4 = ; H1 = H3 = et H2 =
90 90 90

A partir de la formule de Newton-Côtes suivante :


b 4

ò f (x ) dx » ( b - a)å Hi × f (x i )
a i =0

1.25
I= ò f (x ) dx
0.25

é7 32 12 32 7 ù
= (1.25 - 0.25) ´ ê (0.63) + (0.794 ) + (0.909 ) + (1) + (1.15)ú
ë 90 90 90 90 90 û
= 0.89751

Exercice 2 :

Soit f ( x ) une fonction donnée par le tableau :

xi 0.25 0.50 0.75 1.00 1.25


f (x i ) 0.630 0.794 0.909 1.000 1.15

1.25

Trouver une valeur approchée de l’intégrale I = ò f (x ) dx ,


0.25
en utilisant

la formule des Trapèzes généralisée.

Solutions de l’Exercice 2:

Par la formule des Trapèzes généralisée et a partir de la formule des trapèzes


généralisée :

b-a 1.25 - 0.25


n= Û n= =4
h 0.25

_______________________________________________________________________________
[Link] Boumediene -113-
Annexe

Puisque le pas d’intégration est h = 0.25 . De la relation :


b
æ y0 + y n 3
ö
òa f (x ) dx » hççè 2 + å i =1
y i ÷÷
ø
1.25
æ 0.63 + 1.15 ö
I= ò
0.25
f (x ) dx = (0.25) ´ ç
è 2
+ (0.794 + 0.909 + 1) ÷ = 0.89825
ø

Exercice 3 :

Soit f ( x ) une fonction donnée par le tableau :

xi 0.25 0.50 0.75 1.00 1.25


f (x i ) 0.630 0.794 0.909 1.000 1.15

1.25

Trouver une valeur approchée de l’intégrale I =


0.25
ò f (x ) dx , en utilisant

la formule des Simpson généralisée.

Solutions de l’Exercice 2:

Par la formule de Simpson généralisée et puisque le pas d’intégration est


h = 0.25 . De la relation :

b-a 1.25 - 0.25


n=
Û n= =4
h 0.25
A partir de la formule de Simpson généralisée :
b
h
ò f (x ) dx » [(y 0 + y 4 ) + 4(y1 + y 3 ) + 2(y 2 )]
a
3

1.25
æ 0.25 ö
I= ò f ( x ) dx = ç ÷ ´ ((0.63 + 1.15) + 4 ´ (0.794 + 1) + 2 ´ (0.909 )) = 0.89783
0.25 è 3 ø

- Exercices Chapitre 6 -
Exercice 1 :

Trouver, en utilisant la méthode de Déterminant les valeurs propres de :

é 1 - 3ù
A=ê ú ;
ë- 2 2 û

Solutions de l’Exercice 1:

é1 - l - 3 ù
La matrice caractéristique : A - lI 2 = ê ú ;
ë - 2 2 - lû

_______________________________________________________________________________
[Link] Boumediene -114-
Annexe

Le polynôme caractéristique de A :

P2 (l ) = det(A - lI 2 ) = (1 - l )(2 - l ) - 6 = l2 - 3l - 4 .
L’équation caractéristique de A :

l2 - 3l - 4 = 0 .

Les valeurs propres : l = -1 et l = 4 .

Exercice 2 :

Trouver, en utilisant la méthode de Le Verrier les valeurs propres de :

é1 2 2ù
A = êê2 1 2úú ;
êë2 2 1 úû
Solutions de l’Exercice 2 :

Formons les matrices A 2 et A 3 . Nous avons par un calcul direct :

é1 2 2ù é9 8 8ù é41 42 42ù
ê ú 2 ê ú
A = ê2 1 2ú , A = A × A = ê8 9 8ú et A = A × A = êê42 41 42úú
3 2

êë2 2 1 úû êë8 8 9úû êë42 42 41úû

ì S1 = 1 + 1 + 1 = 3
ï
Nous obtenons de A , A 2 , A 3 : íS 2 = 9 + 9 + 9 = 27
ïS = 41 + 41 + 41 = 123
î 3
ì
ï a1 = -S1
ïï 1
Les formules : í a 2 = - (S2 + a1S1 ) ;
ï 2
ïa 3 = - 1 (S3 + a1 S2 + a 2 S1 )
ïî 3

ì
ï a1 = - 3
ïï 1
S’écrivent : í a 2 = - (27 + a1 3)
ï 2
1
ïa 3 = - (123 + a1 27 + a 2 3)
ïî 3

D’où on tire les valeurs a 1 = -3 , a 2 = -9 et a 3 = -5 .

Ainsi le polynôme caractéristique de la matrice A s’écrit :

P3 ( l ) = l3 + a1 l2 + a 2 l + a 3 = l3 - 3l2 - 9l - 5 = 0

_______________________________________________________________________________
[Link] Boumediene -115-
Annexe

Les solutions de l’équation l3 - 3l2 - 9l - 5 = 0 sont les valeurs propres de la


matrice A et sont :
l1 = -1 , l 2 = -1 et l 3 = 5 .

Exercice 3 :

Trouver, en utilisant la méthode de Krylov les valeurs propres de :

é2 - 1 - 2ù
A = êê1 1 1 úú ;
êë2 3 - 1 úû
Solutions de l’Exercice 3 :

Le polynôme caractéristique de A à la forme :

P3 ( l ) = l3 + a1 l2 + a 2 l + a 3

é1ù
Prenons comme vecteur initial Y (0 )
= êê1úú
êë1úû
Puisque la taille de A est 3×3, alors nous avons trois (3) vecteurs a formés :

é2 - 1 - 2ù é1ù é- 1ù
Y = AY = êê1 1
(1 )
1 úú êê1úú = êê 3 úú
(0 )

êë2 3 - 1 úû êë1úû êë 4 úû
é2 - 1 - 2ù é- 1ù é- 13ù
Y = AY = êê1 1
(2 ) (1 )
1 úú êê 3 úú = êê 6 úú
êë2 3 - 1 úû êë 4 úû êë 3 úû
é2 - 1 - 2ù é- 13ù é- 38ù
Y = AY = êê1 1
(3 ) (2 )
1 úú êê 6 úú = êê - 4 úú
êë2 3 - 1 úû êë 3 úû êë - 11 úû
En fin, on obtient : a1 Y ( 2 ) + a 2 Y (1) + a 3 Y ( 0 ) = - Y ( 3 )

é é - 13ù é - 1ù é1ù ù éa1 ù é- 38ù ì- 13a1 - a 2 + a 3 = 38


ê ê1ú ú × êa ú = - ê - 4 ú ï
Où encore : ê êê 6 úú ê3ú
ê ú ê úú ê 2 ú ê ú º í 6a1 + 3a 2 + a 3 = 4
ê êë 3 úû êë 4 úû êë1úû úû êëa 3 úû êë - 11 úû ï 3a + 4a + a = 11
ë î 1 2 3

La solution de ce système est les valeurs a1 = -2 , a 2 = 1 et a 3 = 13 .


Ainsi le polynôme caractéristique de la matrice A s’écrit :

P3 (l ) = l3 + a1 l2 + a 2 l + a 3 = l3 - 2l2 + l + 13 = 0

_______________________________________________________________________________
[Link] Boumediene -116-
Annexe

Les solutions de l’équation l3 - 2l2 + l + 13 = 0 sont les valeurs propres de la


matrice A et sont :
l1 = 1.8681 + 1.9993i , l 2 = 1.8681 - 1.9993i et l 3 = -1.7363 .

- Exercices Chapitre 7 -
Exercice 1 :

Soit le système d’équations non linéaires :

ìx + 2y = 2
í 2 2
î x + 4y = 4
En utilisant la méthode de Newton-Raphson et en partant de
(0 ) é0.2ù
l’approximation initiale X = ê ú , trouver la valeur approchée de la solution
ë1.2 û
avec une précision 10-2.

Solutions de l’Exercice 1 :

(0 ) é0.2ù
L’approximation initiales donnée est X = ê ú et la précision est ε=10-2.
ë1.2 û
Le calcul des éléments de la matrice Jacobienne (les dérivées partielles):

¶f1 ¶f1 ¶f 2 ¶f 2
=1 =2 = 2x = 8y
¶x ¶y ¶x ¶y

La relation (7.3) -voir le cours- s’écrit :

é1 2ù éDx ù é f1 (0.2, 1.2)ù


ê2x 8yú ê Dy ú = - êf (0.2, 1.2)ú
ë û (0.2, 1.2 ) ë û ë2 û
Où encore :
é1 2 ù é Dx ù é- 0.6ù
ê0.4 9.6ú ê Dy ú = ê - 1.8 ú
ë ûë û ë û
Calculons Dx et Dy :

- 0 .6 2 1 - 0.6
- 1.8 9.6 0.4 - 1.8
Dx = = -0.245 et Dy = = -0.177
1 2 1 2
0.4 9.6 0.4 9.6

La première solution approchée X (1 ) est :

_______________________________________________________________________________
[Link] Boumediene -117-
Annexe

éx (0 ) + Dx = 0.2 - 0.245ù é - 0.045ù


X (1 ) = ê (0 ) ú=ê ú
ë y + Dy = 1.2 - 0.177 û ë 1.023 û

é Dx = 0.245 > 10 -2 ù
Et puisque X (1 ) - X (0 ) = ê -2 ú , effectuons une deuxième
êë Dy = 0.177 > 10 úû
itération :

é1 2ù éDx ù éf1 (- 0.045, 1.023)ù


ê2x 8yú ê Dyú = - êf (- 0.045, 1.023 )ú
ë û (-0.045, 1.023 ) ë û ë2 û

é 1 2 ù é Dx ù é - 0.001ù
ê- 0.09 8.184 ú ê Dy ú = ê- 0.188 ú
ë ûë û ë û

Calculons Dx et Dy :

- 0.001 2 1 - 0.001
- 0.188 8.184 - 0.09 - 0.188
Dx = = 0.044 et Dy = = -0.022
1 2 1 2
- 0.09 8.184 - 0.09 8.184

La deuxième solution approchée X (2 ) est :

éx (1 ) + Dx = -0.045 + 0.044 ù é- 0.001ù


X (2 ) = ê (1 ) ú=ê ú
ë y + Dy = 1 .023 - 0 .022 û ë 1.001 û

(2 ) (1 ) é Dx = 0.044 > 10 -2 ù
Et puisque : X -X =ê ú
êë Dy = 0.022 > 10 -2 úû

Nous devons effectuer une troisième itération :

é1 2 ù é Dx ù éf1 (- 0.001, 1.001)ù


ê2x 8y ú ê Dy ú = - êf (- 0.001, 1.001)ú
ë û (-0.001, 1.001 ) ë û ë2 û

é 1 2 ù éDx ù é - 0.001ù
ê- 0.002 8.008 ú ê Dyú = ê - 0.008ú
ë ûë û ë û
Calculons Dx et Dy :

_______________________________________________________________________________
[Link] Boumediene -118-
Annexe

- 0.001 2 1 - 0.001
- 0.008 8.008 - 0.002 - 0.008
Dx = = 0.001 et Dy = = -0.001
1 2 1 2
- 0.002 8.008 - 0.002 8.008

La troisième solution approchée X (3 ) est :

(3 ) éx (2 ) + Dx = -0.001 + 0.001ù é0.000ù


X = ê (2 ) ú=ê ú
ë y + Dy = 1.001 - 0.001 û ë1.000 û
Et puisque :
(3 ) (2 ) é Dx = 0.001 £ 10 -2 ù
X -X =ê -2 ú
ëê Dy = 0.001 £ 10 úû

Nous avons abouti la précision fixée (10-2), la valeur approchée cherchée est :

é0.000 ù
X=ê ú
ë1.000 û

- Exercices Chapitre 8 -
Exercice 1 :

Soit le système de deux équations différentielles suivant:

ìy1¢ (t ) = y 2 (t ) ( y1 ( 0 ) = 2)
í
îy¢2 ( t ) = 2y 2 ( t ) - y1 ( t ) ( y 2 (0) = 1)

Dont la solution analytique est:

ìïy1 ( t ) = 2et - te t
í
ïîy 2 (t ) = e t - te t

Solutions de l’Exercice 1 :

On a alors:
ìf1 (t , y1 ( t ), y 2 ( t )) = y 2 ( t )
í
îf2 ( t, y1 ( t ), y 2 ( t )) = 2y 2 ( t ) - y1 ( t )

et la condition initiale (t0, y1,0, y2,0) = ( 0 , 2 , 1 ) . Si on prend par exemple


t = 0,1, on trouve:

k1,1 = 0.1( f1 (0, 2, 1)) = 0.1


k 2,1 = 0.1(f 2 (0, 2, 1)) = 0

_______________________________________________________________________________
[Link] Boumediene -119-
Annexe

k1,2 = 0.1( f1 (0.05, 2.05, 1.0)) = 0.1


k 2,2 = 0.1( f 2 (0.05, 2.05, 1.0)) = -0.005

k1,3 = 0.1( f1 (0.05, 2.05, 0.9975)) = 0.09975


k 2,3 = 0.1( f 2 (0.05, 2.05, 0.9975)) = -0.0055

k1,4 = 0.1( f1 (0.1, 2.09975, 0.9975)) = 0.09945


k 2,4 = 0.1( f 2 (0.1, 2.09975, 0.9975)) = -0.011075

1
y1,1 = y1,0 + (0.1 + 2(0.1) + 2(0.09975) + 0.09945)
6
= 2.099825
1
y 2,1 = y 2,0 + (0 + 2( -0.005) + 2( -0.0055) + ( -0.011075))
6
= 0.994651

Exercice 2 :

Soit l'équation différentielle :

y¢( t ) = -y( t ) + t + 1 et la condition initiale y(0) = 1.

Solutions de l’Exercice 2 :

On a donc t0 = 0 et y0 = 1, et on prend un pas de temps t = 0,1. De plus, on a:

f ( t, y ) = - y + t + 1

On peut donc utiliser la méthode d'Euler et obtenir successivement des


approximations de y(0,l), y(0,2), y(0,3) …, notées y1, y2, y3… La première
itération produit:

y1 = y 0 + t × f ( t 0 , y 0 ) = 1 + 0.1 × f (0,1) = 1 + 0.1(- 1 + 0 + 1) = 1

La deuxième itération fonctionne de manière similaire:

y2 = y1 + t × f ( t1 , y1 ) = 1 + 0.1 × f (0.1,1) = 1 + 0.1(- 1 + 0.1 + 1) = 1.01

On parvient à:

y3 = y2 + t × f ( t 2 , y2 ) = 1.01 + 0.1 × f (0.2,1.01) = 1.01 + 0.1(- 1.01 + 0.2 + 1) = 1.029

Le tableau suivant rassemble les résultats des dix premières itérations. On


peut montrer que la solution analytique de cette équation différentielle est:

y(t ) = e - t + t , ce qui permet de comparer les solutions numérique et


analytique, et de constater la croissance de l'erreur.

_______________________________________________________________________________
[Link] Boumediene -120-
Annexe

ti y( t i ) yi y(t i ) = y i

0,0 1,000000 1,000000 0,000000


0,1 1,004837 1,000000 0,004837
0,2 1,018731 1,010000 0,008731
0,3 1,040818 1,029000 0,011818
0,4 1,070302 1,056100 0,014220
0,5 1,106531 1,090490 0,016041
0,6 1,148812 1,131441 0,017371
0,7 1,196585 1,178297 0,018288
0,8 1,249329 1,230467 0,018862
0,9 1,306570 1,287420 0,019150
1,0 1,367879 1,348678 0,019201

Exercice 3 :

Soit l'équation différentielle déjà résolue par la méthode d'Euler:

y¢( t ) = -y( t ) + t + 1

et la condition initiale y(0) = 1.

Solutions de l’Exercice 3 :

Dans ce cas:

¶f ¶f
f ( t, y) = -y + t + 1 et = 1 et = -1
¶t ¶y

L'algorithme devient:

t2
y n +1 = y n + ( -y n + t n + 1)t + (1 + ( -1)( -y n + t n + 1))
2

La première itération de la méthode de Taylor d'ordre 2 donne (avec t = 0,1):

( 0 .1 ) 2
y1 = 1 + 0.1( -1 + 0 + 1) + (1 + ( -1)( -1 + 0 + 1)) = 1.005
2

Une deuxième itération donne:

(0.1)2
y 2 = 1.005 + 0.1( -1.005 + 0.1 + 1) + (1 + ( -1)( -1.005 + 0.1 + 1)) = 1.019025
2

Les résultats sont compiles dans le tableau qui suit :

_______________________________________________________________________________
[Link] Boumediene -121-
Annexe

ti y( t i ) yi y(t i ) = y i

0,0 1,000000 1,000000 0,000000


0,1 1,004837 1,005000 0,000163
0,2 1,018731 1,019025 0,000294
0,3 1,040818 1,041218 0,000400
0,4 1,070302 1,070802 0,000482
0,5 1,106531 1,107075 0,000544
0,6 1,148812 1,149404 0,000592
0,7 1,196585 1,197210 0,000625
0,8 1,249329 1,249975 0,000646
0,9 1,306570 1,307228 0,000658
1,0 1,367879 1,368541 0,000662

_______________________________________________________________________________
[Link] Boumediene -122-

Vous aimerez peut-être aussi