Cours M ethodes Num eriques et Programmation

Pr. Souad EL BERNOUSSI

Laboratoire de recherche :
Math ematiques, Informatique et Applications MIA
[Link] : [Link]@[Link]
Table des mati` eres
Chapitre 1. M ethodes num eriques et programmation v
1. Introduction v
2. Arithm etique des calculateurs et Sources derreurs v
3. Propagation des erreurs. ix
4. Conditionnement et stabilit e num erique. x
5. Instabilit e num erique : x
6. S erie de Travaux Dirig es N I xii
Chapitre 2. R esolution dune equation f(x)=0 xiii
1. Introduction xiii
2. M ethode de la bissection. xiii
3. M ethode de Newton-Raphson : xiv
4. M ethode de la s ecante xiv
5. Algorithme de la s ecante : xv
6. M ethode du point xe xv
7. S erie de travaux dirig es N2 xvii
Chapitre 3. R esolution de syst` emes lin eaires xix
1. Introduction xix
2. Rappels sur les syst emes lin eaires xix
3. M ethode Gauss xx
4. Factorisation LU xxiii
5. Mesure des erreurs xxiv
6. M ethodes it eratives xxv
Chapitre 4. Interpolation xxix
1. Introduction xxix
2. Une m ethode directe bas ee sur la r esolution dun syst eme lin eaire : xxx
3. Une m ethode it erative : M ethode de Lagrange xxx
4. Interpolation It er ee de Newton-C otes xxxiv
5. Erreur dInterpolation polynomiale : xxxv
6. Interpolation dHermite xxxvi
Chapitre 5. D erivation et Integration xxxix
1. Introduction : xxxix
2. D erivation. xxxix
3. M ethodes num eriques dint egration. xliii
iii
CHAPITRE 1
M ethodes num eriques et programmation
1. Introduction
1.1. Analyse num erique ? Etude et la construction dalgorithmes (du nom du math ematicien Al Khawarizmi)
de r esolution num erique dun probl` eme donn e.
En pratique, lAnalyse Num erique se propose d etudier les propri et es math ematiques des algorithmes et leur mise
en oeuvre (programmation).
1.2. Objectifs ? Lobjectif de lanalyse num erique est de :
concevoir et d etudier
des m ethodes de r esolution de certains probl` emes math ematiques (en g en eral issus de la mod elisation de probl` emes
r eels), et dont on cherche ` a calculer la solution ou son approximation ` a laide dun ordinateur.
1.3. Enjeux de lanalyse num erique ? R esoudre des probl` emes :
que lon ne sait pas r esoudre autrement
mieux quon ne le faisait avant :
plus pr ecis ement,
moins cher...
1.4. Etique ( Objectifs ) de lanalyse num erique.
Plus vite :
complexit e des algorithmes
complexit e des probl` emes
Plus pr ecis :
erreur darrondi (li ees ` a la machine)
erreur dapproximation (li ees ` a lalgorithme)
Plus able :
stabilit e dun algorithme
Facile ` a programmer :
comprendre pour mieux r eutiliser
1.5. 1` ere Conclusion. La conance aveugle dans ce que lon appelle les r esultats fournis par lordinateur peut
etre la cause derreurs qui peuvent co uter tr` es ch` eres
Alors que faire ?
1.6. Sources dErreurs. Il y a en 3 cat egories :
Erreurs li es ` a la machine,
Erreurs ` a la m ethode (algorithme),
Erreurs sur les donn es (r esultat dun calcul approch e, dune mesure physique,...)
2. Arithm etique des calculateurs et Sources derreurs
Si sophistiqu e quil soit , un calculateur ne peut fournir que des r eponses approximatives.
v
Les approximations utilis ees d ependent ` a la fois des contraintes physiques (espace m emoire, vitesse de lhor-
loge...) et du choix des m ethodes retenues par le concepteur du programme .
Le but de ce chapitre est de prendre connaissance de limpact de ces contraintes et de ces choix m ethologiques.
Dans certains cas il doit etre pris en compte dans lanalyse des r esultats dont une utilisation erron ee pourrait etre
co uteuse.
La premi` ere contrainte est que le syst` eme num erique de lordinateur est discret,
cest ` a dire quil ne comporte quun nombre ni de nombres ;
Il en d ecoule que tous les calculs sont entach es derreurs.
2.1. Evaluation de lerreur. Rappelons dabord quelques notion de base ;
Si X est une quantit e ` a calculer et X

la valeur calcul ee, on dit que :

X X

est lerreur et [ E [=[ X X

[est lerreur absolue.

Exemple :
Si X = 2.224 et X

= 2.223 alors lerreur absolue [ E [=[ X X

[= 2.224 2.223 = 0.001

E
r
=

XX

Xr

est lerreur relative,

X
r
,= 0. X
r
est une valeur de r ef erence pour X. En g en eral ,on prend X
r
= X.
Exemple :
Si X = 2.224 et X

= 2.223
alors , si on prend X
r
= X , lerreur relative E
r
=

XX

Xr

=
|XX

|
|X|
=
0.001
2.224
= 4. 496 10
4
Cependant, si X est la valeur dune fonction F(t) avec a t b, on choisira parfois une valeur de r ef erence
globale pour toutes les valeurs de t.
Exemple :
Si X = sin(t) avec 0 t

4
, on pourra prendre X
r
=

2
2
= sup
0t

4
sin(t).
En g en eral , on ne connait pas le signe de lerreur de sorte que lon consid` ere les erreurs absolues et les erreurs
relatives absolues.
Les op erations el ementaires propagent des erreurs.
Dans la pratique, on consid` ere que :
1) Lerreur absolue sur une somme est la somme des erreurs absolues.
2) Lerreur relative sur un produit ou un quotient est la somme des ereurs relatives.
On peut estimer leffet dune erreur E sur largument x dune fonction f(x) au moyen de la d eriv ee de f(x).
En effet f(x + E) f(x) + Ef

(x)
Exemple :
Calculer la valeur de (11111111)
2
La valeur fournie par une petite calculatrice ` a cinq chiffres est 1, 2345x10
14
Mais la r eponse exacte est 123456787654321.
La machine a donc tronqu e le r esultat ` a 5 chiffres et lerreur absolue est de 6 10
9
.
vi
Lerreur relative est de 0.005% .
Cet exemple montre quil faut etablir clairement lobjectif vis e.
Cet objectif est double ;
1) Nous voulons un bon ordre de grandeur (ici 10
14
) et avoir le maximum de d ecimales exactes,
2) Ce maximum ne peut exc eder la longueur des mots permis par la machine et d epend donc de la machine
2.2. La m emoire de lordinateur : le stockage des nombres. La m emoire dun ordinateur est form ee dun
certain nombre dunit es adessables appel ees OCTETS .
Un ordinateur moderne contient des millions voir des milliards doctets.
Les nombres sont stock es dans un ordinateur comme ENTIERS ou REELS.
2.3. Les nombres entiers : Les nombres entiers sont ceux que lon utilise dhabitude sauf que le plus grand
nombre repr esentable d epend du nombre doctets utilis es :
- avec deux (2) octets, on peut repr esenter les entiers compris entre 32768 et 32767
- avec quatre (4) octets on peut repr esenterr les entiers compris entre 2147483648 et 2147483647
2.4. Les nombres r eels. Dans la m emoire dun ordinateur, les nombres r eels sont repr esent es en notation ot-
tante.
Cette notation a et e introduite pour garder une erreur relative ` a peu pr es constante ; quelque soit lordre de gandeur
du nombre quon manipule.
En notation ottante, un nombre a la forme :
x = Y b
e
b est la base du syst` eme num erique utilis e
Y est la mantisse : une suite de s entier y
1
y
2
...y
s
avec y
1
,= 0 si x ,= 0 et 0 y
i
(b 1)
e est lexposant(un nombre entier relatif)
La norme choisie est celle o` u la mantisse est comprise entre 0 et 1 et o` u le premier chiffre apr` es la virgule est
diff erent de z ero.
Calcul de lerreur
Nous terminons ce chapitre en d enissant les notions de troncature et darrondie.
Exemple :
En base 10, x = 1/15 = 0.066666666......
Dans le cas dune repr esentation tronqu ee nous aurons, pour s = 5, fl(x) = 0.66666 10
1
.
Remarquez comment nous avons modi e lexposant an de respecter la r` egle qui veut que le premier chiffre de
la mantisse ne soit pas nul .
Dans ce cas, lerreur absolue X fl(X) est de 6 10
7.
. Lerreur relative est de lordre de 10
5
Dans une repr esentation tronqu ee ` a s chiffres, lerreur relative maximale est de lordre de 10
s
Dans une repr esentation arrondie, lorsque la premi` ere d ecimale n eglig ee est sup erieure ` a 5, on ajoute 1 ` a la
derni` ere d ecimale conserv ee.
Exemple :
x = 1/15 = 0.066666666.
Nous ecrirons fl(x) = 0.66667 10
1
Lerreur absolue serait alors 3.333 10
7
et lerreur relative serait 5 10
6
En g en eral, lerreur relative dans une repr esentation arrondie ` a s chiffres est de 5 10
(s+1)
soit la moiti e de
celle dune repr esentation tronqu ee.
vii
2.5. Les r egles de base du mod` ele. Pour effectuer une op eration sur deux nombres r eels, on effectue lop eration
sur leurs repr esentations ottantes et on prend ensuite la repr esentation ottante du r esultat.
2.6. laddition ottante. x y = fl(fl(x) + fl(y))
2.7. la soustraction ottante. x y = fl(fl(x) fl(y))
2.8. la multiplication ottante. x y = fl(fl(x) fl(y))
2.9. la division ottante. x y = fl(fl(x)/fl(y))
Chaque op eration interm ediaire dans un calcul introduit une nouvelle erreur darrrondi ou de troncature.
Dans la pratique, il faudra se souvenir du fait que deux expressions alg ebriquement equivalentes peuvent fournir
des r esultats diff erents et que lordre des op erations peut changer les r esultats.
Pour laddition et la soustraction on ne peut effectuer ces 2 op erations que si les exposants sont les m emes. On
transforme le plus petit exposant et donc on ne respecte plus la r egle voulant que le premier chiffre de la mantisse ne
soit pas nul.
2.10. Quelques remarques sur ce mod` ele : On constate une d eviation importante par rapport aux lois habi-
tuelles de larithm etique.
x + (y + z) peut etre diff erent de (x + y) + z.
Exemple :
Pour 4 chiffres signicatifs (s = 4) on a :
(1 + 0.0005) + 0.0005 = 1.000
car
0.1 10
1
+0.5. 10
3
= 0.1. 10
1
+0.00005. 10
1
= 0.1 10
1
+0.0000. 10
1
= 0.1 10
1
et
1 + (0.0005 + 0.0005) = 1.001
Ainsi, laddition ottante nest pas associative .(TD :Sommation dune s erie ` a termes positifs)
On constate aussi que si y est tr es petit par rapport ` a x, laddition de x et y donnera seulement x.
Exemple :
L equation 1 + x = 1 a x = 0 comme unique solution. Mais dans un syst` eme ` a 10 chiffres signicatifs, elle aura
une innit e de solutions (il suft de prendre [ x [ 5 10
11
)
La distributivit e de la multiplication par rapport ` a laddition.
Exemple :
Consid erons lop eration
122 (333 + 695) =
(122 333) + (122 695) =
125416
Si nous effectuons ces deux calculs en arithm etique ` a 3 chiffres (s = 3) et arrondi, nous obtenons :
122 (333 + 695) =
fl(122) fl(1028) =
122 103 10
1
=
(125660) = 126 10
3
(122 333) + (122 695) =
fl(40626) + fl(84790) =
406 10
2
+848 10
2
=
(406 + 848) 10
2
=
(1254 10
2
) =
125 10
3
Donc la distributivit e de la multiplication par rapport ` a laddition nest pas respect ee en arithm etique ottante.
viii
3. Propagation des erreurs.
Une etude de la propagation des erreurs darrondi permattra dexpliquer ce ph enom` ene.
Soit ` a calculer e
x
` a laide de son d eveloppement en s erie qui est convergent pour tout x :
e
x
= 1 +
x
1!
+
x
2
2!
+ ...
Il est evident que dans la pratique il est impossible deffectuer la sommation dune innit e de termes. On arr etera
donc lorsque le terme g en eral
x
k
k!
devient inf erieur ` a 10
t
(on a t digits). Pour x n egatif on sait que le reste de la serie
est inf erieur au premier terme n eglig e donc ` a 10
t
(puisque la serie est alt ern ee).
Les calculs suivant sont fait sur ordinateur pour t = 14.

x
10
15
20
25
30

e
x
4.54.10
5
3.06.10
7
2.06.10
9
1.39.10
11
9.36.10
14

S
4.54.10
5
3.05.10
7
1.55.10
7
1.87.10
5
6.25.10
4

On voit que pour x 20 les r esultats obtenus sont d epourvus de sens. Lexplication de ce ph enom` ene est
la suivante : pour x = 30 les termes de la serie vont en croissant jusqu` a
x
30
30!
= 8.10
11
puis ils d ecroissent et
x
107
107!
9.19.10
15
.
Lerreur absolue sur le terme maximal est de 8.10
11
.10
15
= 8.10
4
. Ainsi le r esultat obtenu pour S repr esente
uniquement laccumulation des erreurs darrondi sur les termes de plus grand module de d eveloppement en serie.
La propagation des erreurs est un des principaux probl` emes en calcul num erique.
Consid erons le cas dune somme :
Dans laddition, les erreurs absolues sadditionnent. Soit en effet
1
et
2
les erreurs absolues sur x
1
et x
2
.
On peut ecrire :
(x
1

1
) + (x
2

2
) = (x
1
+ x
2
) (
1
+
2
)
En arithm etique ottante, lerreur relative est ` a peu pr` es constante et les erreurs absolues peuvent etre approxi-
mativement explicit es par :

1
=[ x
1
[ ,
2
=[ x
2
[
Si les nombres en pr esence ont le m eme signe, lerreur relative reste la m eme que celle quon avait x
1
et x
2
. En
effet

1
+
2
[x
1
+ x
2
[
=
([ x
1
[ + [ x
2
[)
[x
1
+ x
2
[
=
Si par contre les nombres sont de signes diff erents, lerreur relative peut etre ampli ee de facon spectaculaire.
Dans la multiplication, les erreurs relatives sadditionnent.
En effet,soient x
1
et x
2
,on a :
(x
1

1
).(x
2

2
) = x
1
x
2
x
1

2
x
2

2
Si de plus les erreurs s ecrivent

1
= [x
1
[

2
= [x
2
[
on a alors en n egligeant certains termes :
[(x
1

1
).(x
2

2
) x
1
x
2
[
[x
1
x
2
[
=
[x
1

2
x
2

1
[
[x
1
x
2
[
=

1
x
1

2
x
2

2
Des formules equivalentes peuvent donner des r esultats diff erents ; on peut am eliorer le r esultat en utilisant une
formule math ematique equivalente n ecessitant des op erations diff erentes.
ix
3.1. Exemple : Si on consid` ere les nombres

7001 et

7000.
En arithm etique ottante ` a 8 chiffres , on a :

7001= 0.83671979 10
2

7000= 0.83666003 10
2
Donc

7001

7000 = ((0.83671979 0.83666003) 10

2
) = 0.59760000 10
2
On peut obtenir un r esultat plus pr ecis en utilisant lidentit e suivante :

y= (

y)

x +

x +

y
=
x y

x +

y
On obtient alors
1

7001 +

7000
=
1
0.16733798 10
3
= 0.59759297 10
2
4. Conditionnement et stabilit e num erique.
Le fait que certains nombres ne soient pas repr esent es de facon exacte dans un ordinateur entraine que lintro-
duction m eme de donn ee dun probl` eme en machine modie quelque peu le probl` eme initial ; Il se peut que cette
petite variation des donn ees entraine une variation importante des r esultats. Cest la notion de conditionnement dun
probl` eme.
On dit quun probl` eme est bien (ou mal) conditionn e, si une petite variation des donn ees entraine une petite (une
grande) variation sur les r esultats.
Cette notion de conditionnement est li ee au probl` eme math ematique lui m eme et est ind ependante de la m ethode
utilis ee pour le r esoudre.
Une autre notion importante en pratique est celle de stabilit e num erique. Un probl` eme peut etre bien conditionn e
et la m ethode utilis ee pour le r esoudre peut etre sujette ` a une propagation importante des erreurs num eriques.
Ces notions de conditionnement dun probl` eme et de stabilit e num erique dune m ethode de r esolution sont fonda-
mentales en analyse num erique. Si un probl` eme est mal conditionn e alors la solution exacte du probl` eme tronqu e ou
arrondi ` a t digits pourra etre tr` es diff erente de la solution exacte du probl` eme initial. Aucune methode ne pourra rien ;
il faudra essayer de donner une autre formulation au probl` eme.
5. Instabilit e num erique :
Si les erreurs introduites dans les etapes interm ediaires ont un effet n egligeable sur le r esultat nal,on dira que
le calcul ou lalgorithme est num eriquement stable. Si des petits changements sur les donn ees entrainent des petits
changements sur le r esultat. Sinon , on dira que lalgorithme est num eriquement instable.
5.1. Exemple. On veut calculer la valeur de
I
n
=
_
1
0
x
n
a + x
dx
o` u a est une constante plus grande que 1, pour plusieurs valeurs de n. pour ce faire, nous allons exprimer I
n
r ecursivement, i.e. nous allons exprimer I
n
en fonction de n et I
n1
.
x
I
n
=
_
1
0
x
n1
(x + a a)
a + x
dx
=
_
1
0
x
n1
dx a
_
1
0
x
n1
a + x
dx
=
1
n
aI
n1
=
n1

i=0
(a)
i
n i
+ (a)
n
I
0
Comme
I
0
= ln(
1+a
a
)
On peut calculer I
n
pour toutes les valeurs de n.
Mais lalgorithme est num eriquement instable car toute erreur dans le calcul deI
0
= ln(
1+a
a
) va se propager.
En effet si on note par I

0
la valeur approch ee de I
0
et si I

0
= I
0
+ alors
I

n
=
n1

i=0
(a)
i
n i
+ (a)
n
I

0
=
n1

i=0
(a)
i
n i
+ (a)
n
(I
0
+ )
donc [I
n
I

n
[ a
n
.
5.2. Remarques : Il y a en fait diff erentes sources derreur. Nous pouvons les classer en 3 cat egories :
- les erreurs li ees ` a limpr ecision des mesures physiques ou au r esultat dun calcul approch e
- les erreurs li ees ` a lalgorithme utilis e
- les erreurs de calcul li ees ` a la machine
En g en eral, pour lobjet de notre cours, si le premier chapitre met laccent sur les erreurs li ees ` a la machine, nous
nous interresserons beaucoup plus aux erreurs li ees aux m ethodes ou encore aux algorithmes utilis es.
xi
6. S erie de Travaux Dirig es N I
Exercice 1 :Ecrire le nombre d ecimale 25,125 en une repr esentation binaire.
Exercice 2 : En arithm etique ottante avec 3 chiffres signicatifs et arrondi, illustrer la non-validit e des lois dasso-
ciativit e et de distributivit e.
(On pourra prendre : x = 854, y = 251 et z = 852).
Exercice 3 : En arithm etique ottante, avec s = 2, calculer
10

i=1
1
i
2
.
1) En calculant :
1
1
+
1
4
+
1
9
+ ...
1
100
2) En calculant :
1
100
+
1
81
+ ... +
1
1
Quel r esultat est le plus pr ecis et pourquoi ?
Exercice 4 : On consid` ere le polyn ome ax
2
+bx+c = 0 (a ,= 0). On suppose que le discriminant ~ 0. On sait que
(1) x
1
=
b+

b
2
4ac
2a
et x
2
=
b

b
2
4ac
2a
1) V erier que x
1
+ x
2
=
b
a
et x
1
x
2
=
c
a
.
2) Utiliser ce r esultat pour montrer que ces racines peuvent aussi s ecrire sous la forme
(2) x
1
=
2c
b+

b
2
4ac
et x
2
=
2c
b

b
2
4ac
3) Pour s = 4 et trouver les racines de x
2
+ 53.1x + 1 = 0 en calculant
i) x
1
` a partir de (1) et la relation x
1
x
2
=
c
a
.
ii) x
1
` a partir de (2) et la relation x
1
x
2
=
c
a
.
Quel calcul donne le meilleur r esultat et pourquoi ?
iii) Si on clacule dabord x
2
, laquelle des formules (1) et (2) serait-il pr ef erable de choisir et pourquoi ?
Exercice 5 : R esolvez ces deux syst` emes lin eaires :
(1)
_
x + y = 2
x + 1.01y = 2.01
(2)
_
x + y = 2
x + 1.01y = 2.02
Que remarquez-vous ?
Exercice 6 : On cherche les racines de
p(x) = (x 1)(x 2)(x 10)
Elles sont evidentes.
Si on d eveloppe ce polynome en une valeur approch ee pour lune des racines par exemple 10.1
p(x) devient
p(x) = (x 10.1)(x
2
+ bx + c)
Calculer b et c .
En d eduire les racines du polynomes du second degr e. Que remarquez-vous ?
xii
CHAPITRE 2
R esolution dune equation f(x)=0
1. Introduction
Soit f une fonction num erique dune variable r eelle.
On cherche les racines simples de l equation
(1) f(x) = 0
Isoler les racines, cest ` a dire trouver un intervalle [a, b] dans lequel est lunique racine r eelle de (1).
Trouver cet intervalle : th eor` eme des valeurs interm ediaires :
f(a) f(b) < 0 f admet un nombre impair de racines
Si f(a) f(b) > 0 f admet un nombre pair de racines
f

(x) = 1 +
2e
x
(e
x
2)
2
donc f

(x) > 0 pour tout x.

L equation a donc 2 racines simples situ ees de chaque c ot e de ln(2).
On v erie sans probl` eme quune premi` ere racine appartient ` a [1, 0] et la deuxi` eme ` a [1, 2]
On supposera donc d esormais avoir trouv e un intervalle [a, b] o` u f admet une unique racine simple et on supposera
que f est d enie, continue, et autant de fois continument d erivable que n ecessaire.
D EFINITION 1. Nous appellerons algoritnme toute m ethode de r esolution dun probl` eme donn e.
Pour tout probl` eme, nous avons des donn ees et des r esultats.
Les donn ees sont appel ees param` etres dentr ee(input)
les r esultats param` etres de sortie (output)
Ils constituent linterface de lalgorithme.
Les algorithmes classiques que nous allons etudier sont les suivants :
(1) M ethode de la bissection
(2) M ethode de Newton-Raphson
(3) M ethode de la s ecante
(4) M ethode du point xe.
2. M ethode de la bissection.
Soit f(x) continue et cherchons p tel que f(p) = 0.
Supposons quon a localis e un intervalle [a, b] dans lequel la fonction change de signe.
on pose c =
a+b
2
,
si f(a) f(c) 0 on remplace b par c
sinon on remplace a par c,
on continue cette operation jusqu` a ce quon trouve p avec la pr ecision demand ee.
BUT :
Trouver une approximation de la solution de f(x) = 0 dans [a, b]. en construisant une suite dintervalles ([a
n
, b
n
])
n
contenant la racine et tets que a
n
ou b
n
est le milieu de lintervalle [a
n1
, b
n1
].
Entr ees : a, b, , N
0
Sortie : la valeur approch ee p : f(p) = 0
ALGORITHME 2. (1) Si f(a) = 0 imprimer la solution est a. Si f(b) = 0 imprimer la solution est b, aller
` a 10
(2) si f(b) f(a) > 0, imprimer (pas de changement de signe). Aller ` a 10
xiii
(3) poser N = 1
(4) Tant que N N
0
,faire les etapes 5 ` a 8
(5) poser p =
a+b
2
(6) Si f(p) = 0 ou
ba
2
, imprimer p. Aller ` a 10
(7) poser N = N + 1
(8) Si f(a) f(p) > 0,alors poser a = p, sinon poser b = p
(9) Imprimer apr es N
0
it erations lapproximation obtenue est p et lerreur maximale est
ba
2
(10) Fin
3. M ethode de Newton-Raphson :
Le principe consiste ` a construire une suite (x
n
)
n
, telle que x
n+1
soit lintersection de la tangente ` a la courbe de
f au point (x
n
, f(x
n
)) avec laxe horizontal.
On a :
_
A = (x
0
, f(x
0
)), B = (x
1
, 0) axe(Ox)
A et B D : y = ax + b
donc
_
f(x
0
) = ax
0
+ b
0 = ax
1
+ b

_
a = f

(x
0
)
x
1
= x
0

f(x0)
f

(x0)
BUT
Entr ees : une approximation initiale p
0
(la pr ecision d esir ee)
N
0
(le nombre maximum dit erations)
Sortie : valeur approch ee de p ou un message d echec
ALGORITHME 3. (1) N = 1
(2) Tant que N N
0
, faire les etapes 3 ` a 6.
(3) Poser p = p
0

f(p0)
f

(p0)
(4) Si[ p p
0
[ alors imprimer p, aller ` a l etape 8.
(5) Poser N = N + 1.
(6) Poser p
0
= p.
(7) Imprimer la m ethode a echou e apr` es N it erations.
(8) Fin.
4. M ethode de la s ecante
La m ethode de Newton-Raphson suppose le calcul de f

(p). On remplace dans la m ethode de Newton f

(p
n
) par
f(p
n
) f(p
n1
)
p
n
p
n1
.
L equation de la s ecante s ecrit :
s(x) = f(p
n
) + (x p
n
)
f(pn)f(pn1)
pnpn1
Si s(p
n+1
) = 0 , on en d eduit :
p
n+1
= p
n
f(p
n
)
pnpn1
f(pn)f(pn1)
xiv
5. Algorithme de la s ecante :
BUT et Algorithme de la s ecante
Trouver une solution de f(x) = 0
Entr ees : deux approximations initiales p
0
et p
1
(la pr ecision d esir ee)
N
0
(le nombre maximum dit erations)
Sortie :la valeur approch ee de p ou un message d echec
(1) poser N = 1, q
0
= f(p
0
), q
1
= f(p
1
)
(2) Tant que N N
0
+ 1,faire les etapes 3 ` a 6
(3) poser p = p
1
q
1
(p1p0)
q1q0
(4) Si [ p p
1
[ alors imprimer p, aller ` a l etape 8
(5) Poser N = N + 1
(6) Poser p
0
= p
1
, q
0
= q
1
, p
1
= p,q
1
= f(p)
(7) Imprimer la m ethode a echou e apr es N
0
it erations
(8) Fin
6. M ethode du point xe
Nous pouvons observer que la m ethode de Newton peut sinterpr eter comme p
n+1
= g(p
n
) o` u
g(x) = x (
f(x)
f

(x)
).
Maintenant , si la fonction g(x) est continue et si lalgorithme converge (c.` a.d. p
n
p), on tire de p
n+1
= g(p
n
)
que p satisfait l equation p = g(p) ; on dit que p est un point xe de g.
BUT
trouver une solution de g(x) = x
Entr ees : une approximation initiale p
0
(la pr ecision d esir ee)
N
0
le nombre maximale dit erations
Sortie : valeur approch ee de p ou un message d echec
6.1. Convergence et ordre de convergence.
D EFINITION 4. Soit D une partie de R et F une application de D dans D. On dit que la fonction F est contrac-
tante si
x, y D , k [0, 1[ tel que
[ F(x) F(y) [ k [ x y [ .
k est le co efcient de contraction ou de Lipschitz de F.
THEOR` EME 5. Consid erons le segment S = [p
0
a, p
0
+ a] D; si F est contractante sur S et si [ F(p
0
)
p
0
[ (1 k) a, alors lit eration p
n+1
= F(p
n
) de point initial p
0
, converge vers lunique point xe p S de F.
THEOR` EME 6. Si F est diff erentiable au voisinage dun point xe p et si [ F

(p) [< 1 alors :

V voisinage de p tels que p
0
V et p
n+1
= F(p
n
) converge vers p.
xv
6.2. Ordre de convergence.
D EFINITION 7. Consid erons une suite p
n
convergeant vers p et posons e
n
= p
n
p. On dit dans le cas o` u
_

en
en1

_
converge, que la suite p
n
converge lin eairement vers p ou encore que la m ethode est du premier ordre.
Si on a
_

en
(en1)
k

_
converge, alors la convergence est dite dordre k.
Exemple :
La m ethode de Newton est une m ethode de type point xe avec
F(x) = x
f(x)
f

(x)
.
Si x

est racine simple de f(x) = 0, alors f

(x

) ,= 0 et il existe un voisinage V de x

tel que pour tout p

0
V,
la suite (p
n
)
n
converge vers x

et lordre de convergence est 2.

Pour d eterminer lordre de convergence on utilise la formule de Taylors en x

: F(x) = F(x

) + F

(x

)(x
x

) + F

(x)
(xx

)
2
2
).
xvi
7. S erie de travaux dirig es N2
Exercices 1
R esoudre ` a laide de la m ethode de bisection tanx x = 0 dans lintervalle [4; 4.7].
Exercice 2
On consid` ere l equation
(1) e
x
4x = 0
1) D eterminer le nombre et la position approximative des racines de (1) situ ees dans .x 0
2) Utiliser lalgorithme de bissection pour d eterminer la plus petite de ces racines ` a pr` es.(par exemple 10
7
)
3) Sans faire dit erations, d eterminer combien vous devriez en faire pour calculer la plus grande racine ` a laide
de la bissection avec une pr ecision de 10
8
, si lintervalle de d epart est [2; 2, 5]
Exercice 3

Ecrire un algorithme pour calculer par la m ethode de Newton la racine k-i` eme dun nombre.
Quelle est la valeur de s =
_
2 +
_
2 +

2 + ..... ?
Suggestion : ecrire .p
n+1
=G(p
n
),p
0
=0 Quel est lordre de convergence ?
Exercice 4

Ecrire 3 m ethodes it eratives pour la r esolution de x

3
x 1 =0 et v erier exp erimentalement leur convergence
avec x
0
= 1, 5. Trouver ` a 10
6
pr` es la racine comprise entre 1 et 2. Connaissant la valeur de cette racine, calculer
lordre de convergence de vos 3 m ethodes. Ce r esultat coincide-t-il avec lexp erience ?
Exercice 5
R esoudre x
2
-1=0 en utilisant la m ethode de la s ecante avec x
0
= 3 et x
1
= 5/3. Quarrivera-t-il si on choisit et
x
0
= 5/3 et x
1
= 3 ? Expliquez.
Exercice 6
On consid` ere la fonction f(x) = x
4
2x
3
+ 2x 1.
1) Calculer lordre de convergence de la m ethode de Newton pour la racine x

= 1.
2) Etudier la convergence et lordre de la m ethode de points xe suivante :
x
n+1
= x
n
3
f (x
n
)
f

(x
n
)
, avec x
0
pr` es de x

= 1 racine de f.
xvii
CHAPITRE 3
R esolution de syst` emes lin eaires
1. Introduction
Un syt` eme lin` eaire s ecrit sous la forme :
(1) Ax = b
A est une matrice n n ` a coefcients r eels,
b R
n
x R
n
Les m ethodes de r esolution sont de deux types :
(1) Les m ethodes directes :
obtenir la solution en un nombre ni dop erations.
(2) Les m ethodes it eratives :
construire une suite (x
n
)
n
qui converge vers la solution.
Dans ce chapitre nous allons :
(1) Rapeler des notions et notations de base relatives aux syst` emes lin eaires et aux matrices
(2) Etudier une m ethode directe : la m ethode de Gauss.
(3) Etudier la d ecomposition (factorisation) LU.
(4) Etudier des applications : Inverse de matrices,...
2. Rappels sur les syst emes lin eaires
Si A est inversible alors (1) admet une solution unique
x = A
1
b
Ainsi th eoriquement le probl` eme revient ` a calculer A
1
? Mais en pratique ce calcul est difcile.
M ethodes classiques pour r esoudre (1) sans calculer A
1
.
Exemple
(1)
_
x + 2y = 5
2x + y = 4
2.1. La m ethode de Cramer : x =

5 2
4 1

1 2
2 1

=
3
3
= 1 et y =

1 5
2 4

1 2
2 1

=
6
3
= 2
2.2. La m ethode de substitution (ou d elimination).
_
x + 2y = 5
2x + y = 4

_
x = 2y + 5
2x + y = 4

_
x = 2y + 5
2(2y + 5) + y = 4

_
x = 2y + 5
3y = 6

_
x = 2y + 5
y = 2

_
x = 1
y = 2
Peut-on g en eraliser ces m ethodes pour un syst eme de n equations avec n N?
Th eoriquement OUI mais en pratique cela va n ecessiter beaucoup de calculs et de techniques.
xix
3. M ethode Gauss
Lid ee de base consiste ` a transformer (1) en un probl` eme que lon sait r esoudre.
(1) Si la matrice A = D avec D une matrice diagonale, alors on sait r esoudre (1).
(2) Si la matrice A = U (ou L) avec U (ou T) triangulaire sup erieure ( ou inf erieure) alors on sait r esoudre (1).
Probl` eme : Comment tranformer une matrice en une matrice triangulaire inf erieure ou sup erieure ?
Cas n = 3 :
AX = b s ecrit :
8
<
:
a11x1 + a12x2 + a13x3 = b1
a21x1 + a22x2 + a23x3 = b2
a31x1 + a32x2 + a33x3 = b3
la formeaugment ee
(A b) =
0
@
a11 a12 a13
a21 a22 a23
a31 a32 a33
b1
b2
b3
1
A
On suppose quea11 ,= 0, par elimination, on obtient :
(A1 b1) =
0
@
a11 a12 a13
0 a

22
a

23
0 a

32
a

33
b1
b

2
b

3
1
A
(1)
8
<
:
a11x1 + a12x2 + a13x3 = b1 (l1)
a21x1 + a22x2 + a23x3 = b2 (l2)
a31x1 + a32x2 + a33x3 = b3 (l3)
On note par (li) la i
` eme
equation du syst` eme pr ecedent.
On suppose que a11 ,= 0,
On pose :
(l

2
) = a11(l2) a21(l

1
) et (l

3
) = a11(l3) a31(l

1
)
Alors (1) s ecrit
(2)
8
<
:
a11x1 + a12x2 + a13x3 = b1 (l1)
a

22
x2 + a

23
x3 = b

2
(l

2
)
a

32
x2 + a

33
x3 = b

3
(l

3
)
On suppose que a

22
,= 0. On pose :
(l

3
) = a

22
(l

3
) a

32
(l

2
) Alors (2) s ecrit
(2)
8
<
:
a11x1 + a12x2 + a13x3 = b1 (l1)
a

22
x2 + a

23
x3 = b

2
(l

2
)
a

33
x3 = b

3
(l

3
)
Remarque :
(1) Les termes diagonaux sont appel es les pivots
(2) Si un pivot aii est nul on change de ligne (on permute) de i ` a n (pivotage partiel)
(3) Cette m ethode se g en eralise assez facilement bien quil faut etre prudent avec le choix du pivot. En pratique, il faut eviter
de prendre des pivots trop petits.
3.1. Algorithme :

Elimination de Gauss.
3.1.1. Partie 1 : R eduction ` a la forme triangulaire. Entr ee A et b Sortie A = U (forme triangulaire), et b.
Pour j = 1, ..., (n 1)
Pour i = j + 1, ..., n
lij
a
ij
a
jj
Pour k = j + 1, ..., n
a
ik
a
ik
lija
jk
Fin
bj bj lijbj
Fin
Fin
Sortie A = U (forme triangulaire), et b
xx
(1) Poser j = 1
(2) Tant que j n 1 faire
(3) Si ajj = 0 afcher pivot nul aller ` a etape 14, sinon
(4) Poser i = j + 1
(5) Tant que i n faire
(6) lij =
a
ij
a
jj
(7) Si lij = 0 ; aller ` a l etape 12
(8) Poser k = j + 1
(9) Tant que k n faire
(10) a
ik
= a
ik
lija
jk
, k = k + 1 ; aller ` a l etape 9.
(11) bi = bi lijbj ;
(12) poser i = i + 1; Aller ` a l etape 5.
(13) j = j + 1 ; Aller ` a l etape 2.
(14) Fin
Exemple :
8
>
>
<
>
>
:
x + y + 3t = 4
2x + y z + t = 1
3x y z + 2t = 3
x + 2y + 3z t = 4
qui s ecrit encore :
0
B
B
@
1 1 0 3
2 1 1 1
3 1 1 2
1 2 3 1
1
C
C
A
0
B
B
@
x
y
z
t
1
C
C
A
=
0
B
B
@
4
1
3
4
1
C
C
A
Nous appliquons lalgorithme ` a notre exemple en travaillant sur la matrice augment ee.
2
6
6
6
6
4
A b
1 1 0 3 . 4
2 1 1 1 . 1
3 1 1 2 . 3
1 2 3 1 . 4
3
7
7
7
7
5
2
6
6
4
1 1 0 3 . 4
0 1 1 5 . 7
0 4 1 7 . 15
0 3 3 2 . 8
3
7
7
5
2
6
6
4
1 1 0 3 . 4
0 1 1 5 . 7
0 0 3 13 . 13
0 0 0 13 . 13
3
7
7
5
Que lon peut ecrire sous la forme :
8
>
>
<
>
>
:
x + y + 3t = 4
y z 5t = 7
3z + 13t = 13
13t = 13
Notons que l etape j = 3 nous donnerait l43 = 0.
Nous avons maintenant un syst` eme triangulaire ` a r esoudre.
xxi
3.1.2. Partie 2 : Remont ee triangulaire. Entr ee A, b avec A matrice triangulaire sup erieure
Sortie x solution du syt` eme Ax = b
(1) xn =
bn
ann
(2) Pour i = n 1, n 2, ..., 1 faire :
xi =
1
aii
(bi
n
X
j=i+1
aijxj)
En appliquant cet algorithme ` a notre exemple, nous obtenons x = (1, 2, 0, 1).
Remarque :
(1) Dans la pratique le test (3) de lalgorithme d elimination de Gauss ne conduit pas ` a larr et. En fait, si le pivot est nul, on
cherche, dans la m eme colonne, un el ement dindice plus grand non nul, puis on echange les lignes correspondantes. Si
ceci est impossible, le syst` eme est singulier.
(2) On est parfois amen e, pour des raisons de stabilit e num erique, ` a effectuer des echanges de lignes m eme si le test (3) est
n egatif (cest ` a dire que le pivot est non nul). Ceci conduit ` a des strat egies dites de pivot que nous n etudierons pas ici.
Exemple :
8
<
:
2x + 6y + 10z = 0
x + 3y + 3z = 2
3x + 14y + 28z = 8

0
@
2 6 10
1 3 3
3 14 28
1
A
0
@
x
y
z
1
A
=
0
@
0
2
8
1
A
0
@
2 6 10 0
1 3 3 2
3 14 28 8
1
A

0
@
2 6 10 0
0 0 4 4
0 5 13 8
1
A

0
@
2 6 10 0
0 5 13 8
0 0 4 4
1
A
En utilisant la remont e on trouve :
8
<
:
z =
4
4
= 1
y =
1
5
(8 13 (1)) = 1
x =
1
2
(6 1 10 (1)) = 2
x

=
0
@
2
2
1
1
A
M ethode de Gauss avec normalisation :Elle consiste ` a normaliser le pivot : On a :
(1)
8
<
:
a11x1 + a12x2 + a13x3 = b1 (l1)
a21x1 + a22x2 + a23x3 = b2 (l2)
a31x1 + a32x2 + a33x3 = b3 (l3)
On note par (li) la i
` eme
equation du syst` eme pr ecedent.
On suppose que a11 ,= 0,
(l1) s ecrit : x1 +
a
12
a
11
x2 +
a
13
a
11
x3 =
b
1
a
11
(l

1
)
Si on pose :
(l

2
) = (l2) a21(l

1
)
et
(l

3
) = (l3) a31(l

1
)
Alors (1) s ecrit
(2)
8
<
:
x1 +
a
12
a
11
x2 +
a
13
a
11
x3 =
b
1
a
11
(l

1
)
a

22
x2 + a

23
x3 = b

2
(l

2
)
a

32
x2 + a

33
x3 = b

3
(l

3
)
On suppose que a

22
,= 0,
(l

2
) s ecrit x2 + a

23
x3 = b

2
(l

2
)
Si on pose :
(l

3
) = (l

3
) a

32
(l

2
)
si a

33
,= 0 on pose (l

3
) x3 =
b

3
a

33
Alors (2) s ecrit
xxii
(2)
8
>
<
>
:
x1 +
a
12
a
11
x2 +
a
13
a
11
x3 =
b
1
a
11
(l

1
)
x2 + a

23
x3 = b

2
(l

2
)
x3 =
b

3
a

33
(l

3
)
4. Factorisation LU
Supposons que dans l elimination de Gauss on nutilise aucune strat egie de pivotage et que tous les pivots a
(k)
k,k
,= 0.
Dans ce cas le passage de A
(k)
A
(k+1)
(1 k N 1) revient ` a multiplier ` a gauche la matrice A
(k)
par la matrice
N N :
E
(k)
=
0
B
B
B
B
B
B
B
B
B
B
B
B
B
@
1 0 0
0 1 0
.
.
.
.
.
.
.
.
.
.
.
.
.
.
. 0 1 0
.
.
.
.
.
. l
k+1,k
.
.
.
.
.
.
.
.
.
.
.
.
.
.
. 0
.
.
.
.
.
.
0 0 l
N,k
0 1
1
C
C
C
C
C
C
C
C
C
C
C
C
C
A
Avec l
i,k
=
a
(k)
i,k
a
(k)
k,k
, pour k + 1 i N.
La matrice A
(N)
, qui est triangulaire sup erieure, est alors egale ` a
A
(N)
= MA
avec M = E
(N1)
E
(N2)
...E
(1)
M est le produit de matrices triangulaires inf erieures, donc M est aussi triangulaire inf erieure, on a det(M) =
N1
Q
i=1
det(E
(i)
) =
1 et linverse de M est aussi triangulaire inf erieure.
En posant U = A
(N)
et L = M
1
, on a A = LU.
(E
(k)
)
1
=
0
B
B
B
B
B
B
B
B
B
B
B
B
B
@
1 0 0
0 1 0
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
. 0 1
.
.
. 0
.
.
.
.
.
. l
k+1,k
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
. 0
.
.
. 0
0 0 l
N,k
0 ..0 1
1
C
C
C
C
C
C
C
C
C
C
C
C
C
A
L = M
1
= (E
(1)
)
1
(E
(2)
)
1
.....(E
(N1)
)
1
L =
0
B
B
B
B
B
B
B
B
B
B
B
B
B
@
1 0 0
l2,1 1 0
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
. 1 0
.
.
.
.
.
. l
k+1,k
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
l
N,1
l
N,k
l
N,N1
1
1
C
C
C
C
C
C
C
C
C
C
C
C
C
A
Exemple :
A =
0
@
5 2 1
5 6 2
4 2 1
1
A
= A
(1)
E
(1)
=
0
@
1 0 0
1 1 0
4
5
0 1
1
A
; A
(2)
=
0
@
5 2 1
0 8 1
0
18
5
9
5
1
A
E
(2)
=
0
@
1 0 0
0 1 0
0
9
20
1
1
A
; A
(3)
=
0
@
5 2 1
0 8 1
0 0
9
4
1
A
= U
xxiii
L =
0
@
1 0 0
1 1 0
4
5
9
20
1
1
A
et A = LU
THEOR` EME 8. Soit A = (ai,j)
1i,jN
une matrice carr ee dordre N telle que les N sous-matrices de A :
0
B
B
B
B
@
a11 a
1k
.
.
.
.
.
.
.
.
.
.
.
.
a
k1
a
kk
1
C
C
C
C
A
, 1 k N
soient inversibles,
alors il existe une matrice triangulaire inf erieure L = (lij)
1iN;1jN
, avec
lii = 1 (1 i N), et une matrice triangulaire sup erieure U telles que A = LU.
De plus cette factorisation est unique.
Si A = LU on peut r esoudre : Ax = b en r esolvant
(1) Lz = b
(2) Ux = z
(1) Ax = b

(2) Lz = b
(3) Ux = z
Dans ce cas Ax = LUx = L(Ux) = Lz = b .
Les syst` emes (2) et (3) etant triangulaires, la r esolution ne n ecessite que lex ecution dune remont ee et dune descente trian-
gulaire.
Exemple :
0
B
B
@
1 1 0 3
2 1 1 1
3 1 1 2
1 2 3 1
1
C
C
A
=
0
B
B
@
1 0 0 0
2 1 0 0
3 4 1 0
1 3 0 1
1
C
C
A
0
B
B
@
1 1 0 3
0 1 1 5
0 0 3 13
0 0 0 13
1
C
C
A
.
Exemple :
A = LU =
0
B
B
@
1 0 0 0
2 1 0 0
3 4 1 0
1 3 0 1
1
C
C
A
0
B
B
@
1 1 0 3
0 1 1 5
0 0 3 13
0 0 0 13
1
C
C
A
, on r esoud le syst` eme Ax =
0
B
B
@
4
1
3
4
1
C
C
A
.
Lz = b
8
>
>
<
>
>
:
z1 = 4
2z1 + z2 = 1
3z1 + 4z2 + z3 = 3
z1 3z2 + z4 = 4
z =
0
B
B
@
4
7
13
13
1
C
C
A
.
Ux = z
8
>
>
<
>
>
:
13x4 = 13
3x3 + 13x4 = 13
x2 x3 5x4 = 7
x1 + x2 + 3x4 = 4
x =
0
B
B
@
1
2
0
1
1
C
C
A
.
5. Mesure des erreurs
Lutilisation dun calculateur pour implanter les algorithmes etudi es conduira in evitablement ` a des erreurs. Pour mesurer
celles-ci, nous devons mesurer la distance entre le vecteur repr esentant la solution exacte x = (x1, ..., xn) et le vecteur x =
( x1, ..., xn) repr esentant la solution approch ee. Nous pouvons, pour ce faire, utiliser la longueur usuelle de R
n
i.e. :
|x|2 = |
n
X
1
x
2
i

1
2
pourtant, dans la pratique on lui pr ef` ere souvent la longueur
|x| = max
1in
[xi[
xxiv
Par exemple si x = (1, 7, 2, 4) alors |x| = 7 .
Exemple :
Si x = (1, 1, 1, 1) alors |x| = 1 si x = (1.01, 1.1, 1, 1), on a
|x x| = 0.1
Consid erons alors le syst` eme
0
B
B
@
10 7 8 7
7 8 6 5
8 6 10 9
7 5 9 10
1
C
C
A
0
B
B
@
x1
x2
x3
x4
1
C
C
A
=
0
B
B
@
32
23
33
31
1
C
C
A
dont la solution exacte est x = (1, 1, 1, 1) .
Si dans le membre de droite nous remplacons b par :

b = (32.06; 22.87; 33.07; 30.89)

nons obtenons
x = (9.19; 12.59, 4.49, 1.09)
Cest-` a-dire quune erreur relative de lordre de :
|b

b|
|b|
= 3 10
1
sur b a entran e une erreur relative de lordre de
|x x|
|x|
= 13.52
sur la solution.
Nous devons donc soupconner que lapplication de larithm etique nie ` a la r esolution dun tel syst` eme serait d esastreuse.
L etude de cette question d epasse le cadre de ce programme.
6. M ethodes it eratives
On va voir un type de m ethodes it eratives de r esolution du syst` eme lin eaire Ax = b sous la forme :
(3)

x
(0)
vecteur arbitraire,
x
(k+1)
= Bx
(k)
+ c, k 0
lorsque
Ax = b x = Bx + c,
la matrice B et le vecteur c sont en fonction de A et b.
D enition
La m ethode it erative (3) est convergente si
lim
k+
x
(k)
= x, x
(0)
.
Remarque :
Si on Pose e
(k)
= x
(k)
x , pour k = 0, 1, ...
Comme x = Bx + c et x
(k+1)
= Bx
(k)
+ c,on a
e
(k)
= Bx
(k1)
Bx = Be
(k1)
= ... = B
k
e
(0)
lim
k+
x
(k)
= x lim
k+
e
(k)
= 0 lim
k+
B
k
e
(0)
= 0
Donc
La m ethode it erative (3) est convergente si
lim
k+
B
k
v = 0, v
ce qui equivaut ` a
lim
k+
[[B
k
v[[ = 0, v
pour toute norme vectorielle [[.[[.
6.1. Convergence des m ethodes it eratives :
THEOR` EME 9. Les propositions suivantes sont equivalentes :
(1)la m ethode it erative (3) est convergente ;
(2)(B) < 1;
(3)[[B[[ < 1 pour au moins une norme matricielle subordonn ee [[.[[.
xxv
6.2. M ethode de Jacobi, de Gauss-Seidel. Ces m ethodes sont des cas particuliers de la m ethode suivante :
A = M N avec M inversible et assez simple. On aurait alors :
Ax = b Mx = Nx + b
x = M
1
Nx + M
1
b
x = Bx + c,
avec B = M
1
N et c = M
1
b.
6.2.1. M ethode de Jacobi : En posant A = D (E + F),
M = D est la diagonale de A et N = E + F.
Ax = b Dx = (E + F)x + b.
On suppose que D est inversible, cest ` a dire aii ,= 0, 1 i N.
La matrice
J = D
1
(E + F) = IN D
1
A
est appel ee la matrice de Jacobi.

x
(0)
donn e,
Dx
(k+1)
= (E + F)x
(k)
+ b, k 0
A chaque etape, on calcule les N composantes x
(k+1)
1
, .., x
(k+1)
N
du vecteur x
(k+1)
:
8
>
>
>
>
<
>
>
>
>
:
a1,1x
(k+1)
1
= a1,2x
(k)
2
... a1,Nx
(k)
N
+ b1
a2,2x
(k+1)
2
= a21x
(k)
1
a2,3x
(k)
3
... a2,Nx
(k)
N
+ b2
.
.
.
aN,Nx
(k+1)
N
= a
(k)
N,1
x
(k)
1
... a
(k)
N,N1
x
(k)
N1
+ bN
6.2.2. M ethode de Gauss-Seidel. M est la partie triangulaire inf erieure de A : M = D E et N = F.
On pourrait am eliorer la m ethode pr ec edente en utilisant les quantit es d ej` a calcul ees, on calcule successivent les N compo-
santes
x
(k+1)
1
, .., x
(k+1)
N
du vecteur x
(k+1)
:
8
>
>
>
>
<
>
>
>
>
:
a1,1x
(k+1)
1
= a1,2x
(k)
2
... a1,Nx
(k)
N
+ b1
a2,2x
(k+1)
2
= a21x
(k+1)
1
a2,3x
(k)
3
... a2,Nx
(k)
N
+ b2
.
.
.
aN,Nx
(k+1)
N
= aN,1x
(k+1)
1
. .. aN,N1x
(k+1)
N1
+ bN
qui revient ` a ecrire :
Dx
(k+1)
= Ex
(k+1)
+ Fx
(k)
+ b
ou encore
(D E)x
(k+1)
= Fx
(k)
+ b

x
(k+1)
= (D E)
1
Fx
(k)
+ (D E)
1
b
/1 = (D E)
1
F est la matrice de Gauss-Seidel.
Elle est inversible si aii ,= 0, 1 i N.
Exemples : Etude de la convergence des m ethodes de Jacobi et de Gauss-Seidel dans le cas o` u la matrice du syst` eme lin eaire
est :
1- A1 =
0
@
4 1 0
1 4 1
0 1 4
1
A
;
Les deux m ethodes convergent.
2- A2 =
0
@
2 1 1
1 2 1
1 1 2
1
A
;
la m ethode de Jacobi diverge et La m ethode de Gauss-Seidel converge.
3- A =
0
@
1 2 2
1 1 1
2 2 1
1
A
.
la m ethode de Jacobi converge et la m ethode de Gauss-Seidel diverge.
subsectionconvergence des m ethodes de Jacobi, de Gauss-Seidel
THEOR` EME 10. Soit A une matrice sym etrique d enie positve.
On suppose que
A = M N, avec M inversible.
Si la matrice sym etrique
M
T
+ N
xxvi
est d enie positive,
alors
(M
1
N) < 1.
Exemple :
A =
0
B
B
B
B
@
2 1 0 0
1
.
.
.
.
.
. 0
0
.
.
.
.
.
. 1
0 0 1 2
1
C
C
C
C
A
1- Pour tout u R
N
,
u
T
Au = u
2
1
+ u
2
N
+
N
P
i=2
(ui ui1)
2
=la matrice sym etrique A est d enie positive.
2- Pour la m ethode de Jacobi :
M = D et N = E + F :
A = M N est sym etrique d enie positive, car
u
T
(M
T
+ N)u = u
2
1
+ u
2
N
+
N
P
i=2
(ui + ui1)
2
> 0, si u ,= 0.
M
T
+ N est d enie positive =(M
1
N) < 1.
xxvii
CHAPITRE 4
Interpolation
1. Introduction
Nous abordons dans ce chapitre un nouveau type de probl` eme, faisant intervenir la notion dapproximation dune fonction.
Exemples :
1) Dapr es la Formule de Taylor ` a lordre 5 de la fonction sin(x), on a :
x V ois(0), sin(x) x
x
3
3!
+
x
5
5!
+ sin
(6)
()
x
6
6!
Lerreur commise serait de lordre de sin
(6)
()
x
6
6!
Ainsi :
Si N = 3, sin(0.1) = (0.1)
(0.1)
3
3!
= 9. 983 3 10
2
Si N = 5, sin(0.1) = 0.1
(0.1)
3
3!
+
(0.1)
5
5!
= 9. 983 3 10
2
Avec le logiciel Maple on a : sin(0.1) = 9. 983 3 10
2
2) Avec les cours danalyse I et II, on ne connait pas dexpression explicite de I =
R
1
0
e
x
2
dx
Cependant dapr es :
La formule du trap` eze I =
R
1
0
e
x
2
dx
f(0)+f(1)
2
=
1+e
1
2
= 0.683 94
La formule de Simpson : I =
R
1
0
e
x
2
dx
1
6

f(0) + 4f(
1
2
) + f(1)

=
1
6
(1 + 4e

1
4
+ e
1
) = 0.747 18
Avec le Logiciel Maple, on a :
R
1
0
e
x
2
dx = 0.746 82
On ne connait pas ` a ce niveau du cours lexpression explicite de lerreur.
La notion dapproximation dune fonction consiste ` a remplacer un probl` eme donn e par un probl` eme voisin.
La question fondamentale serais de savoir la qualit e de cette approximation.
Remarque : En pratique la fonction f est connue explicitement, ou seulement par ses valeurs en quelques points.
La notion dinterpolation polynomiale est la facon la plus simple dobtenir une telle approximation.
Th eor` eme : Soit f une fonction continue dans [a, b] IR, alors pour tout > 0 donn e, il existe un polyn ome Pn de degr e n
tel que
Max
x[a,b]
[f(x) Pn(x)[ <
(1) Linterpolation polyn omiale est un outil pour la construction des m ethodes dint egration num erique ou des m ethodes
dapproximation des equations diff erentielles.
(2) Linterpolation par les fonctions splines est largement utilis ee dans tous les programmes de dessin assist e par ordinateur,
conception assist ee par ordinateur ou plus g en eralement de graphisme.
(3) Les s eries de Fourier et leur analogue discret, la transformation de Fourier discr` ete : Elles sont un moyen tr` es utile pour
lapproximation des fonctions p eriodiques.
Remarque :
Pour les equations aux d eriv ees partielles, la m ethode des el ements nis, un des outils de base de ling enierie moderne,
utilise de facon essentielle linterpolation multi-dimensionnelle
Une facon naturelle dapprocher les fonctions p eriodiques est dutiliser les polyn omes trigonom etrique.
Nous allons nous limiter ` a lintroduction de linterpolation Polyn omiale. Elle consiste ` a d eterminer un polyn ome Pn(x) de
degr e n qui puisse remplacer lors des applications la fonction f(x).
De plus, cest un outil efcace pour :
Calculer, pour x donn e, une approximation de f(x) en calculant Pn(x)
Construire :
(1) des m ethodes dint egration num erique
xxix
(2) des m ethodes de diff erentation
(3) des m ethodes dapproximation des equations diff erentielles
(4) ...
Le principe est simple, le proc ed e est le suivant :
On choisit (ou on se donne) (n + 1) points x0, x1, ..., xn.
On calcule y0 = f(x0), ..., yn = f(xn)
ou on se donne (xi, yi), i = 0, ..., n.
On cherche un polyn ome de degr e n tel que Pn(xi) = yi, i = 0, ..., n .
Remarque :
(1) Les points (xi, yi)i=0,n sont appel es points dinterpolation.
(2) Si la fonction f est connue seulement par ses valeurs en quelques points, les (n + 1) points x0, x1, ..., xn sont x es..
(3) Si on veut que Pm(xi) = f(xi) et P

m
(xi) = f

(xi), i = 0, ..., n , on obtient linterpolation dite dHermite

Il existe plusieurs techniques pour calculer Pn(x). Les plus connues sont celles de Lagrange et de Newton-C otes. Nous allons
en fait le faire de deux facons :
(1) Une m ethode directe bas ee sur la r esolution dun syst eme lin eaire
(2) Une m ethode it erative due ` a Lagrange.
(3) Une m ethode it erative dHermite.
Nous terminerons ce chapitre par :
(1) Une br` eve discution sur lerreur dinterpolation polynomiale
(2) Une br` eve description du principe de la m ethode it er ee de Newton-C otes
2. Une m ethode directe bas ee sur la r esolution dun syst eme lin eaire :
On se donne (n + 1) points x0, x1, ..., xn .
On calcule y0 = f(x0), ..., yn = f(xn) .
On cherche un polyn ome de degr e n tel que Pn(xi) = yi, i = 0, ..., n .

Ecrivons explicitement Pn(xi) = yi.

anx
n
i
+ an1x
n1
i
+ ... + a1xi + a0 = yi, i = 0, ..., n
Sous forme matricielle :
0
B
B
B
@
x
n
0
x
n1
0
x0 1
x
n
1
x
n1
1
x1 1
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
x
n
n
x
n1
n
xn 1
1
C
C
C
A
0
B
B
B
@
an
an1
.
.
.
a0
1
C
C
C
A
=
0
B
B
B
@
y0
y1
.
.
.
yn
1
C
C
C
A
Matrice de type Vandermonde. Son d eterminant est
det =
Y
i<j
(xi xj)
On a det ,= 0 si tous les xi sont distincts. On peut donc trouver un unique vecteur de coefcients (an, ..., a0) r esolvant le
probl` eme.
3. Une m ethode it erative : M ethode de Lagrange
3.1. Interpolation Lin eaire : On consid` ere deux points (x0, y0), (x1, y1) avec :

x0 ,= x1
y0 = f(x0) et y1 = f(x1).
Pour d eterminer le polyn ome P1(x) = ax + b) qui passe par deux points distincts (x0, y0), (x1, y1) (x0 ,= x1). On peut :
1) R esoudre le syst` eme d equations :

ax0 + b = y0
ax1 + b = y1
do` u
xxx
(
a =
(y
1
y
0
)
(x
1
x
0
)
b = y0 ax0 =
x
1
y
0
x
0
y
1
x
1
x
0
On a :
P1(x) =
(y1 y0)
(x1 x0)
x + (
x1y0 x0y1
x1 x0
)
et
P1(x0) = y0 et P1(x1) = y1
2) Poser
L0(x) =
x x1
x0 x1
L1(x) =
x x0
x1 x0
On a :
L
k
(xi) =

0 si i ,= k
1 si i = k
Ainsi,
P1(x) = y0L0(x) + y1L1(x)
= y0
(x x1)
(x0 x1)
+ y1(
x x0
x1 x0
)
=
(y1 y0)
(x1 x0)
x + (
x1y0 x0y1
x1 x0
)
On a :
P1(x0) = y0 et P1(x1) = y1
car
L
k
(xi) =

0 si i ,= k
1 si i = k
Ces deux proc ed es d eterminent evidemment le m eme polyn ome de d egr e 1 (la m eme droite).
Si maintenant, on veut d eterminer le polynome de degr e 2 qui passe par trois (3) points distincts alors :
i) la premi` ere expression de P1(x) est inad equate (il faut refaire les calculs)
ii) la dexi` eme expression se pr ete assez facilement ` a une g en eralisation par r ecurrence.
Exemple :
D eterminer le polyn ome dinterpolation P1(x) de degr e 1 tel que
P1(xi) = f(xi), i = 0, 1
avec yi = f(xi) i = 0, 1 , (x0, y0) = (0, 1), (x1, y1) = (2, 5)
Dapr es la m ethode de Lagrange,
P1(x) = y0L0(x) + y1L1(x)
= y0
(x x1)
(x0 x1)
+ y1(
x x0
x1 x0
)
= 1
(x 2)
(0 2)
+ 5
(x 0)
(2 0)
= 1
(x 2)
(2)
+ 5
(x)
(2)
= 2x + 1
xxxi
3.2. Interpolation parabolique. On consid ere trois points (x0, y0), (x1, y1) et (x2, y2) avec :

x0 ,= x1 ,et x0 ,= x2 et x1 ,= x2
y0 = f(x0), y1 = f(x1) et y2 = f(x2).
Pour d eterminer le polyn ome P2(x) de d egr e 2, d equation y = ax
2
+bx+c qui passe par trois points distincts (x0, y0), (x1, y1)
et (x2, y2), il suft de poser :
L0(x) =
(x x1)(x x2)
(x0 x1)(x0 x2)
L1(x) =
(x x0)(x x2)
(x1 x0)(x1 x2)
L2(x) =
(x x0)(x x1)
(x2 x0)(x2 x1)
On a :
L
k
(xi) =

0 si i ,= k
1 si i = k
Ainsi
P2(x) = y0L0(x) + y1L1(x) + y2L2(x)
= y0
(x x1)(x x2)
(x0 x1)(x0 x2)
+
y1
(x x0)(x x2)
(x1 x0)(x1 x2)
+
y2
(x x0)(x x1)
(x2 x0)(x2 x1)
est le polyn ome dinterpolation polyn omiale associ e.
Exemple :
D eterminer le polyn ome dinterpolation P2(x) de degr e 2 tel que
P2(xi) = f(xi), i = 0, 1 et 2
avec yi = f(xi) i = 0, 1 et 2 , (x0, y0) = (0, 1), (x1, y1) = (1, 2) et (x2, y2) = (2, 5)
Dapr es la m ethode de Lagrange,
P2(x) = y0L0(x) + y1L1(x) + y2L2(x)
= y0
(x x1)(x x2)
(x0 x1)(x0 x2)
+ y1
(x x0)(x x2)
(x1 x0)(x1 x2)
+ y2
(x x0)(x x1)
(x2 x0)(x2 x1)
= 1
(x 1)(x 2)
(1)(2)
+ 2
(x)(x 2)
(1)(1)
+ 5
(x)(x 1)
(2)(1)
= 1
(x 1)(x 2)
2
+ 2
(x)(x 2)
1
+ 5
(x)(x 1)
(2)(1)
= x
2
+ 1
Remarque :
(1) Pour calculer P2(x) ,on na pas utilis e le polynome P1(x) calcul e dans lexemple pr ec edent et pourtant on avait deux
points communs.
(2) Li(x), i = 0, 1, 2 sont des polyn omes de degr e 2 :
L0(x) =
(x1)(x2)
(1)(2)
=
1
2
(x 1) (x 2) =
1
2
x
2

3
2
x + 1
L1(x) =
(x)(x2)
(1)(1)
= x(x 2) = x
2
+ 2x
L2(x) =
(x)(x1)
(2)(1)
=
1
2
x(x 1) =
1
2
x
2

1
2
x
xxxii
3.3. Interpolation de Lagrange.
On choisit n + 1 points x0, x1, ..., xn .
On calcule y0 = f(x0), ..., yn = f(xn) .
On cherche un polyn ome de degr e n tel que Pn(xi) = yi, i = 0, ..., n .
On introduit les coefcients dinterpolation de Lagrange.
L
k
(x) =
(x x0)...(x x
k1
)(x x
k+1
)...(x xn)
(x
k
x0)...(x
k
x
k1
)(x
k
x
k+1
...(x
k
xn)
L
k
(x) =
j=n
Y
j=0 j=k
(x xj)
(x
k
xj)
L
k
(x) est un polyn ome de degr e n,
L
k
(xi) =

0 si i ,= k
1 si i = k
Donc
P(x) = y0L0(x) + y1L1(x) + ... + ynLn(x) =
n
X
k=0
y
k
L
k
(x)
est un polyn ome de degr e n qui v erie bien P(xi) = yi
Propri et e : Le Polyn ome dinterpolation polyn omiale est unique.
En effet si P(x) et Q(x) sont deux polyn omes dinterpolation alors :
P(x) Q(x) est un polyn ome de degr e n pour lequel
P(xi) Q(xi) = 0, i = 0, ..., n.
Ce polyn ome de degr e n ayant n + 1 racines, il est identiquement nul.
Exemple :
On suppose que f(x) =
3

x et que (x0, y0) = (0, 0), (x1, y1) = (1, 1) et (x2, y2) = (8, 2)
1) D eterminer le polyn ome P2(x) dinterpolation polyn omiale qui passent par les points (xi, yi)i=0,2
Dapr es la m ethode de Lagrange,
P2(x) = y0L0(x) + y1L1(x) + y2L2(x)
= y0
(x x1)(x x2)
(x0 x1)(x0 x2)
+ y1
(x x0)(x x2)
(x1 x0)(x1 x2)
+ y2
(x x0)(x x1)
(x2 x0)(x2 x1)
= 0
(x 1)(x 2)
(0 1)(0 2)
+ 1
(x 0)(x 8)
(1 0)(1 8)
+ 2
(x 0)(x 1)
(8 0)(8 1)
P2(x) = 1
x(x 8)
7
+ 2
x(x 1)
56
=
3
28
x
2
+
31
28
x
On a bien P2(0) = 0, P2(1) = 1 et P2(8) =
3
28
(8)
2
+
31
28
8 = 2
2) Calculer P2(x) et f(x) =
3

x pour x = 0.5, 0.95, 1, 1.5 et 3. On a :

x f(x) P2(x) =
3
28
x
2
+
31
28
x
0.5 0.793 7 0.526 79
0.95 0.983 05 0.955 09
1 1 1
1.5 1. 144 7 1. 419 6
3 1. 442 2
33
14
= 2. 357 1
Remarque :
1) En pratique, on utilise linterpolation polyn omiale avec des polyn omes de d egr e n assez grand ou linterpolation po-
lyn omiale par morceaux. Ainsi dans lexemple pr ecedent, il faut augmenter le nombre de points dinterpolations.
2) Si les valeurs y
k
sont des valeurs exp erimentales. Linterpolation polynomiale est une technique peu appropri ee pour de
telles situations. Les polyn omes de degr e elev e sont sensibles ` a la perturbation des donn ees.
3) La m ethode de Lagrange sadapte mal au changement du nombre de points (xi, yi)i. On ne peut utiliser les coefcients de
Lagrange si on passe de n ` a (n + 1) points.
xxxiii
4) Ph enom` ene de RUNGE (fonction de Runge) : Linterpolation polyn omiale ne fournit pas une bonne approximation de la
fonction f(x) =
1
1+25x
2
. Si on augmente le nombre de points dinterpolation le resultat devient plus mauvais.
4. Interpolation It er ee de Newton-C otes
On choisit n + 1 points x0, x1, ..., xn .
On calcule y0 = f(x0), ..., yn = f(xn) .
On cherche un polyn ome de degr e n tel que Pn(xi) = yi, i = 0, ..., n .
LInterpolation It er ee de Newton-C otes est un proc ed e it eratif qui permet de calculer le polyn ome dinterpolation Pn(x) de
d egr e n bas e sur (n + 1) points (xi, yi)i=0,n
` a partir du polyn ome dinterpolation P
(n1)
(x) de d egr e (n 1) bas e sur n points
(xi, yi)
i=0,(n1)
, en posant :
Pn(x) = P
(n1)
(x) + C(x), n 1
avec
C(x) = an(x x0)(x x1)...(x x
(n1)
)
an =
n
X
k=0
f(x
k
)
(x
k
x0)...(x
k
x
(k1)
)(x
k
x
(k+1)
)...(x
k
xn)
Les co efcients an sont appel es diff erences divis ees dordre n de la fonction f, on note :
an = f [x0, x1, ..., xn]
On appelle diff erence divis ee dordre 0 de f en un point x :
f [x] = f(x)
Diff erence divis ee dordre 1 de f en x et y :
f [x, y] =
f [x] f [y]
x y
on a
f [x, y] =
f(x)
x y
+
f(y)
y x
Diff erence divis ee dordre 2 de f en x, y et z :
f [x, y, z] =
f [x, y] f [y, z]
x z
=
f(x)
(x y) (x z)
+
f(y)
(y x) (y z)
+
f(z)
(z x) (z y)
et plus g en eralement :
f [x1, x2, ..., xn] =
n
X
i=1
f(xi)
n
Q
k=1
k=i
(xi x
k
)
Remarque :
Les diff erences divis ees sont ind ependants de lordre des points.
Quel est le lien entre f(x) et les diff erences divis ees ?
Soit x un point autre que les n + 1 points xi, i = 1, ..., n. On a
f [x, x0] =
f(x) f [x0]
x x0
d

o` u
f(x) = f [x0] + (x x0) f [x, x0]
xxxiv
mais comme
f [x, x0, x1] =
f [x, x0] f [x0, x1]
x x1
alors
f(x) = f [x0] + (x x0) f [x0, x1] (x x0)(x x1)f [x, x0, x1]
en continuant ainsi de proche en proche on obtient :
f(x) = f [x0] + (x x0) f [x0, x1]
+ ... + (x x0) ... (x xn1) f [x0, ..., xn] +
(x x0)...(x xn)f [x, x0, ..., xn]
on v erie que
f(x) = Pn(x) + L(x)f [x, x0, ..., xn]
o` u Pn(x) est un polyn ome de degr e n tel que Pn(xi) = f(xi), pour i = 0, ..., n. Cest donc le polyn ome dinterpolation de
Lagrange, on lappelle le polyn ome de Newton.
5. Erreur dInterpolation polynomiale :
Lerreur commise lors dune interpolation est une question fondamentale en analyse num erique :
elle renseigne ` a priori sur la nature de cette erreur
elle fournit des informations sur les termes qui y participent
elle permet davoir un ordre de grandeur de lerreur commise.
Th eor` eme :
Soient f une fonction de classe C
n+1
dans I et , (xi)i=0,n (n + 1) points distincts dans I avec x0 < x1 < ... < xn
Alors pour tout x [x0, xn] , il existe = (x) tel que
f(x) Pn(x) =
f
(n+1)
()
(n + 1)!
(x x0)(x x1)...(x xn) =
f
(n+1)
()
(n + 1)!
L(x)
Pn(x) = y0L0(x) + y1L1(x) + ... + ynLn(x) =
n
X
k=0
y
k
L
k
(x)
L
k
(x) =
n
Y
j=0 j=k
(x xj)
(x
k
xj)
L(x) = (x x0)(x x1)...(x xn)
Remarque :
1) Cette formule montre que :
i) lerreur est nulle pour x = xi i.e. x est un point dinterpolation.
ii) lerreur d epend de la fonction consid er ee ( de f
(n+1)
) et des points dinterpolations (xi)i.
2) Cette formule derreur permet de trouver des formules derreur pour lint egration num erique et la differentiabilit e num erique.
Dans le cas de lerreur dinterpolation ` a partir de la forme de Newton, on a :
f(x) Pn(x) = L(x).f[x, x0, ..., xn].
Comme on a la m eme fonction f selon les m emes points xi pour i = 0, ..., n, il sagit de deux formes du m eme polyn ome, et
lerreur dinterpolation est la m eme, do` u
f(x) Pn(x) =
f
(n+1)
()
(n + 1)!
L(x) = L(x).f[x, x0, ..., xn].
Exemple :
D eterminer lerreur dinterpolation polynomiale dans le cas de linterpolation parabolique
On approche la fonction f(x) par la parabole passant par les points (x0, y0) = (0, 1), (x1, y1) = (1, 2) et (x2, y2) = (2, 5).
Le polyn ome dinterpolation P2(x) de degr e 2 tel que P2(xi) = f(xi), i = 0, 1 et 2
avec yi = f(xi) i = 0, 1 et 2 , (x0, y0) = (0, 1), (x1, y1) = (1, 2) et (x2, y2) = (2, 5)
Dapr es la m ethode de Lagrange,
xxxv
P2(x) = y0L0(x) + y1L1(x) + y2L2(x)
= 1
(x 1)(x 2)
(1)(2)
+ 2
(x)(x 2)
(1)(1)
+ 5
(x)(x 1)
(2)(1)
= x
2
+ 1
Dapr` es le th eor` eme pr ec edent,
f(x) P2(x) =
f
(3)
()
3!
(x x0)(x x1)(x x2)
=
f
(3)
()
3!
x(x 1)(x 2)
Si

f
(3)
(x)

M alors
x [0, 2] , [f(x) P2(x)[
[f(x) P2(x)[
M
6
[x(x 1)(x 2)[

M
6
x(x 1)(x 2)
6.4 10
2
M.
(le maximum de u(x) = x(x 1)(x 2) est atteint en x

=
3

3
3
; do` u
1
6
u(x

) =
1
6
3

3
3

3
3
1

3
3
2

= 0.0
641 5 6.4 10
2
).
6. Interpolation dHermite
Linterpolation de Lagrange, qui fournit facilement un polyn ome prenant des valeurs donn ees, pr esente linconv enient dans
certains cas de donner une qualit e m ediocre : entre les points dinterpolation la diff erence entre une fonction et son polyn ome
dinterpolation peut etre grande, m eme si le nombre de points tend vers linni (ph enom` ene dit de Runge).
Pour rem edier ` a cela on peut essayer dutiliser non seulement les valeurs dune fonction mais aussi celles de ses d eriv ees : cest
linterpolation dHermite.
On se donne (xi, f
(k)
(xi)), pour i = 0, ..., n, k = 0, ..., mi o` u mi N.
En posant N =
P
n
i=0
(mi +1), on peut montrer que si les noeuds xi sont distincts, il existe un unique polyn ome HN1 PN1,
appel e polyn ome dinterpolation dHermite, tel que
H
(k)
N1
(xi) = y
(k)
i
, i = 0, ..., n, k = 0, ..., mi.
Ce polyn ome s ecrit
HN1(x) =
P
n
i=0
P
m
i
k=0
y
(k)
i
L
ik
(x) o` u y
(k)
i
= f
(k)
(xi), i = 0, ..., n, k = 0, ..., mi. Les fonctions L
ik
(x)sont appel ees po-
lyn omes caract eristiques dHermite et sont d enies par les relations
d
p
x
p
(L
ik
)(xj) =

1 si i = j, k = p
0 sinon
6.1. Interpolation dHermite, cas o` u k = 1. En d enissant les polyn omes
li0(x) =
n
Q
s=0
s=i
(
xxs
x
i
xs
)
2
, i = 0, ..., n
li1(x) = (x xi)
n
Q
s=0
s=i
(
xxs
x
i
xs
)
2
, i = 0, ..., n
et en posant
Li1(x) = li1(x), pour i = 0, ..., n,
Li0(x) = li0(x) l

i0
(xi)Li1(x).
Soit f une fonction de classe C
2n+2
sur un intervalle I contenant les points deux ` a deux distincts xi; i = 0; ...; 2n + 1, rang es
dans lordre croissant. Alors pour tout x [x0; x2n+1] il existe (au moins) un r eel x dans ce m eme intervalle tel que :
Concernant lerreur dinterpolation, on a lestimation f(x) H2n+1(x) =
f
(2n+1)
(x)
2n+1!
N(x), o` u x I(x; x0, ..., xn) et N est
le polyn ome de degr e (2n+1) d eni par 2n+1(x) = (x x0)
2
(x x1)
2
...(x xn)
2
.
xxxvi
6.1.1. Interpolation dHermite : Exemple. Soit f une fonction de classe C
2
sur lintervalle I = [0 ; 1]. Soit H(x) un polyn ome
de degr e 3 tel que :
H(0) = f(0)
H(1) = f(1)
H

(0) = f

(0)
H

(1) = f

(1)
D eterminons H. En appliquant les formules , on trouve
H(x) = f(0)(1 x)
2
(1 + 2x) + f(1)x
2
(3 2x) + f

(0)x(1 x)
2
+ f

(1)x
2
(x 1)
z
i
x
i
w
i
y
i
d
i
=
w
i+1
w
i
z
i+1
z
i
dd
i
=
d
i+1
d
i
z
i+2
z
i
ddd
i
=
dd
i+1
dd
i
z
i+3
z
i
z
1
= x
1
w
1
= y
1
d
1
= y

1
z
2
= x
1
w
2
= y
1
dd
1
=
d
2
d
1
z
3
z
1
d
2
=
w
3
w
2
z
3
z
2
ddd
1
=
dd
2
dd
1
z
4
z
1
z
3
= x
2
w
3
= y
2
dd
2
=
d
3
d
2
z
4
z
2
d
3
= y

2
z
4
= x
2
w
4
= y
2
xxxvii
CHAPITRE 5
D erivation et Integration
1. Introduction :
Si f est une fonction d erivable sur [a, b], la d eriv ee en c ]a, b[ est d enie par :
f

(c) = lim
h0
f(c)
h
o` u f(c) = f(c + h) f(c)
Si f est une fonction continue sur [a, b], lintegrale de f sur [a, b] est d enie par
Z
b
a
f(x)dx = lim
h0
R(h)
o` u R(h) =
n
X
k=1
f(a + kh).h
R(h) est la somme de Riemann avec h =
ba
n
.
Il existe des fonctions simples comme
sin x
x
ou
p
cos
2
x + 3 sin
2
x qui nont pas de primitive connue.
f peut- etre connue seulement en quelques points et sa formule est inconnue (exp : r esultats experimentaux,...),
Comment peut-on integrer de telles fonctions entre a et b?
Si P(x) est une approximation de f dans lintervalle [a, b], nous nous proposons d etudier les approximations :
f

(y) P

(y) y [a, b]
et
Z
b
a
f(x)dx
Z
b
a
P(x)dx.
2. D erivation.
2.1. D eriv ee premi` ere. La d erivation num erique nous permet de trouver une estimation de la d eriv ee ou de la pente dune
fonction, en utilisant seulement un ensemble discret de points. Soit f une fonction connue seulement par sa valeur en (n + 1) points
donn es xi i = 0, 1, ..., n distincts.
On suppose connue la valeur de la fonction en xi1, xi et xi+1; on pose f(xi1) = yi1, f(xi) = yi et f(xi+1) = yi+1.
Si on suppose que lespace entre deux points successifs est constant, donc on pose h = xi xi1 = xi+1 xi.
Alors les formules standarts en deux points sont :
Formule de difference progressive :
f

(xi)
f(xi+1) f(xi)
xi+1 xi
=
yi+1 yi
xi+1 xi
.
Formule de difference r egressive :
f

(xi)
f(xi) f(xi1)
xi xi1
=
yi yi1
xi xi1
.
Formule de difference centrale :
f

(xi)
f(xi+1) f(xi1)
xi+1 xi1
=
yi+1 yi1
xi+1 xi1
.
Exemple :
Pour illustrer les trois formules, consid erons les donn ees suivantes :
(x0, y0) = (1, 2); (x1, y1) = (2, 4); (x2, y2) = (3, 8); (x3, y3) = (4, 16) et (x4, y4) = (5, 32).
Nous voulons estimer la valeur de f

(x2).
xxxix
(1) Progressive : f

(x)
f(x
3
)f(x
2
)
x
3
x
2
=
168
43
= 8.
(2) Regressive : f

(x)
f(x
2
)f(x
1
)
x
2
x1
=
84
32
= 4.
(3) Centrale : f

(x)
f(x
3
)f(x
1
)
x
3
x1
=
164
42
= 6.
Les donn ees ont et e calcul e pour la fonction f(x) = 2
x
. f

(x) = 2
x
ln(2) et pour x = 3 f

(3) = 2
3
ln(2) = 5.544.
Remarque :
En utilisant la formule de Taylor :
f(x + h) = f(x) + hf

(x) +
h
2
2
f

().
x x + h
Formule progressive :
h = xi+1 xi
f

(xi) =
f(xi+1) f(xi)
h

h
2
f

()
xi xi+1
lerreur est
h
2
f

() donc en O(h).
Cette formule peut etre trouv ee aussi en utilisant le polyn ome dinterpolation de Lagrange pour les points (xi, f(xi)) et
(xi+1, f(xi+1)). Formule regressive :
h = xi xi1
f

(xi) =
f(xi) f(xi1)
h
+
h
2
f

()
xi1 xi
La formule de diff erence centrale de la d eriv ee en xi peut etre trouv ee en utilisant la formule de Taylor dordre 3 avec
h = xi+1 xi = xi xi1
f(xi+1) = f(xi) + hf

(xi) +
h
2
2
f

(xi) +
h
3
3!
f

(1)
f(xi1) = f(xi) hf

(xi) +
h
2
2
f

(xi)
h
3
3!
f

(2)
xi 1 xi+1, xi1 2 xi
si on suppose que f

est continue sur [xi1, xi+1] on peut ecrire la formule suivante :

f

(xi) =
f(xi+1) f(xi1)
2h
+
h
2
6
f

()
xi1 xi+1
lerreur est
h
2
6
f

() donc en O(h
2
). La formule de diff erence centrale peut aussi etre trouv ee ` a partir du polyn ome dint erpolation
de Lagrange en 3 points.
On peut interpoler les donn ees par un polynome au lieu dutiliser la droite, nous obtenons alors les formules de diff erence qui
utilisent plus de deux points. On suppose que le pas h est constant.
Formule de diff erence progressive utilisant trois points :
f

(xi)
f(xi+2) + 4f(xi+1) 3f(xi)
xi+2 xi
Formule de diff erence r egressive utilisant trois points :
f

(xi)
3f(xi) 4f(xi1) + f(xi2)
xi xi2
Exemple : Formules de diff erence en trois points :
En utilisant les donn ees de lexemple pr ec edent, on trouve :
f

(xi)
32+4(16)3(8)
2
= 4 progressive.
f

(xi)
3(8)4(4)+2
2
= 5 regressive.
xl
2.2. Formule g en erale en trois points. La formule dapproximation en 3 points de la d eriv ee premi` ere, bas ee sur le polyn ome
dinterpolation de Lagrange, nutilise pas des points equidistants.
Etant donn e trois points (x1, y1); (x2, y2) et (x3, y3) avec x1 < x2 < x3, la formule suivante permet dapprocher la d eriv ee
en un point x [x1, x3]. Les d eriv ees aux points xi sont les suivantes :
f

(x1) =
2x1 x2 x3
(x1 x2)(x1 x3)
y1+
x1 x3
(x2 x1)(x2 x3)
y2 +
x1 x2
(x3 x1)(x3 x2)
y3;
f

(x2) =
x2 x3
(x1 x2)(x1 x3)
y1+
2x2 x1 x3
(x2 x1)(x2 x3)
y2 +
x2 x1
(x3 x1)(x3 x2)
y3;
f

(x3) =
x3 x2
(x1 x2)(x1 x3)
y1+
x3 x1
(x2 x1)(x2 x3)
y2 +
2x3 x2 x1
(x3 x1)(x3 x2)
y3;
Le polyn ome de Lagrange est donn ee par
P(x) = L1(x)y1 + L2(x)y2 + L3(x)y3
o` u
L1(x) =
(x x2)(x x3)
(x1 x2)(x1 x3)
L2(x) =
(x x1)(x x3)
(x2 x1)(x2 x3)
L3(x) =
(x x1)(x x2)
(x3 x1)(x3 x2)
Lapproximation de la d eriv ee premi` ere est donn ee par f

(x) P

(x), qui peut secrire

P

(x) = L

1
(x)y1 + L

2
(x)y2 + L

3
(x)y3
o` u
L

1
(x) =
2x x2 x3
(x1 x2)(x1 x3)
L

2
(x) =
2x x1 x3
(x2 x1)(x2 x3)
L

3
(x) =
2x x1 x2
(x3 x1)(x3 x2)
donc
f

(x) =
2x x2 x3
(x1 x2)(x1 x3)
y1 +
2x x1 x3
(x2 x1)(x2 x3)
y2 +
2x x1 x2
(x3 x1)(x3 x2)
y3.
2.3. D eriv ees dordre sup erieur. Les formules de d eriv ees dordre sup erieur, peuvent etre trouv ees ` a partir des d eriv ees du
polyn ome de Lagrange ou en utilisant les formules de Taylor.
Par exemple, etant donn e 3 points xi1, xi, xi+1
equidistants, la formule de la d eriv ee seconde est donn ee par :
f

(xi) =
1
h
2
[f(xi+1) 2f(xi) + f(xi1)]
lerreur est en O(h
2
).
D eriv ee seconde ` a partir du polyn ome de Taylor.
f(x + h) = f(x) + hf

(x) +
h
2
2
f

(x) +
h
3
3!
f

(x) +
h
4
4!
f
(4)
(1)
f(x h) = f(x) hf

(x) +
h
2
2
f

(x)
h
3
3!
f

(x) +
h
4
4!
f
(4)
(2)
x 1 x + h et x h 2 x.
f

(x)
f(x + h) + f(x h) 2f(x)
h
2
lerreur est en O(h
2
).
xli
Pour obtenir les formules de la troisi` eme et la quatri` eme d eriv ee, on prend une combinaison lin eaire des d eveloppement de
Taylor, pour f(x + 2h), f(x + h), f(x h) et f(x 2h).
La table suivante donne diff erentes formules centrales toutes en O(h
2
) :
f

(xi)
1
2h
[f(xi+1) f(xi1)]
f

(xi)
1
h
2
[f(xi+1) 2f(xi) + f(xi1)]
f

(xi)
1
2h
3
[f(xi+2) 2f(xi+1) + 2f(xi1) f(xi2)]
f
(4)
(xi)
1
h
4
[f(xi+2) 4f(xi+1) + 6f(xi) 4f(xi1) + f(xi2)] .
En utilisant les polyn omes dinterpolation de Lagrange les d eriv ees dordre p sont calcul ees par :
f
(p)
()
n
X
i=0
Ai()f(xi)
o` u
Ai() = L
(p)
i
() p n
n
X
i=0
Ai()x
k
i
= 0 0 k p 1
n
X
i=0
Ai()x
k
i
= k(k 1)...(k p + 1)
kp
p k n.
Remarque :
: La fomule est exacte pour les polyn omes de degr es n.
: Le syst` eme lin eaire donnant les Ai() a un d eterminant de type Vandermonde diff erent de z ero si les xi sont distincts.
: Les Ai() sont ind ependants de f et peuvent etre calcul es une fois pour toutes.
2.4. Etude de lerreur commise. Dapr es le chapitre pr ec edent, si f est connue en (n + 1) points xi, i = 0, ..., n alors
f(x) = Pn(x) + e(x), o` u e(x) est lerreur dinterpolation. En d erivant on obtient
f

(x) = P

n
(x) + e

(x)
=
i=n
X
i=0
Ai(x).f(xi) + e

(x)
et e

(x) =
d
dx

1
(n + 1)!
L(x).f
(n+1)
(x)

=
d
dx
(L(x).f[x0, ..., xn, x])
=
1
(n + 1)!
L

(x).f
(n+1)
(x) +
1
(n + 1)!
L(x).
d
dx

f
(n+1)
(x)

Remarque
Lerreur de d erivation est nulle si f est un polyn ome de degr e inf erieur ou egale ` a n.
Si on prend pour x un point xi, le second terme de la d erni` ere somme sannule, sinon il faut connatre
d
dx

f
(n+1)
(x)

, ce
qui est difcile car la fonction x x
etant inconnue.
On peut donner une forme si f est n + 2 fois d erivable en utilisant la notion de diff erence. En effet
d
dx

f
(n+1)
(x)

=
d
dx
(f[x0, ..., xn, x])
= lim
h0
f[x0, ..., xn, x + h] f[x0, ..., xn, x]
h
= lim
h0
f[x0, ..., xn, x, x + h]
= lim
h0
1
(n + 2)!
f
(n+2)
(
x,h
).
On constate quon devra se contenter dune estimation
[ e(x) [
1
(n + 1)!
[ L

(x) [ Mn+1 +
1
(n + 2)!
[ L(x) [ Mn+2.
xlii
3. M ethodes num eriques dint egration.
Soit f : [a, b] R, une fonction continue donn ee. On d esire approcher num eriquement la quantit e
R
b
a
f(x)dx.
3.1. Formules ferm ees. On appelle ainsi les formules quand la fonction f est continue sur lintervalle [a, b].
Les points dinterpolation xi verient a = x0 < x1 < ... < xn1 < xn = b.
La formule des rectangles est une formule dite ` a un point x0 = a. Le polyn ome dinterpolation associ e est P0(x) = f(a) et
L(x) = x a pour tout x appartenant ` a [a, b]. Do` u
I(f) I(P0) = f(a)(b a).
Linterpr etation graphique consiste donc ` a remplacer
R
b
a
f(x)dx par laire du rectangle de base [a, b] et de hauteur f(a).
La formule des trap` ezes est une formule ` a 2 points : x0 = a et x1 = b. Le polyn ome de Lagrange associ e ` a ces deux points
est P1(x) = f(a)

xb
ab

+ f(b)

xa
ba

.Do` u
I(f) I(P1) =
Z
b
a
P1(x)dx =
f(a) + f(b)
2
(b a).
La formule de Simpson est une formule ` a trois points x0 = a , x1 =
a+b
2
et x2 = b : . Le polyn ome associ e ` a ces trois points
est P2(x) = f(a)L0(x) + f(
a+b
2
)L1(x) + f(b)L3(x). Notons que
L0(x) =
(x x1)(x b)
(a x1)(a b)

Z
b
a
L0(x)dx =
(b a)
6
,
L1(x) =
(x a)(x b)
(x1 a)(x1 b)

Z
b
a
L1(x)dx =
4(b a)
6
,
L2(x) =
(x x1)(x a)
(b x1)(b a)

Z
b
a
L2(x)dx =
(b a)
6
,
On tire donc la formule suivante :
I(f) I(P2) =
(b a)
6

f(a) + 4f(
a + b
2
) + f(b)

.
On appelle ainsi les formules quand la fonction f est continue sur lintervalle ]a, b[. Les points dinterpolation xi verient
a < x0 < x1 < ... < xn1 < xn < b.
Il en existe une innit e.
Une ` a 1 points avec x0 =
a+b
2
:
I(f) (b a)f(
a + b
2
)
Cette formule est exacte pour tout polyn ome de degr e 1.
Une ` a 2 points avec x0 =
2a+b
3
et x1 =
a+2b
3
:
I(f)
b a
2

2a + b
3

+ f

2b + a
3

.
Cette formule est exacte pour tout polyn ome de degr e 1.
Une ` a 3 points avec x0 =
3a+b
4
et x1 =
a+b
2
et x2 =
3b+a
4
:
I(f)
b a
6

4f

3a + b
4

+ 2f

a + b
2

2f

a + 3b
4

.
Cette formule est exacte pour tous les polyn omes de degr e 2.
3.2. Etude g en erale de lerreur commise. Estimer lerreur E(f) = I(f) I(Pn) avec pr ecision ? ?.
Si f est sufsamment d erivable, on a
E(f) = I(f Pn) =
Z
b
a

1
(n + 1)!
f
(n+1)
(x)L(x)

dx.
Th eor` eme : Supposons que E(f) = 0 pour les polyn omes de degr e au plus n et que la fonction f C
n+1
([a, b]). On dit alors
que la m ethode est dordre n + 1. Si on pose Mn+1 = max
x[a,b]
[ f
(n+1)
(x) [, Une premi` ere estimation de lerreur est
[ E(f) [
1
(n + 1)!
Mn+1
Z
b
a
[ L(x) [ dx.
xliii
Th eor` eme : En plus des hypoth` eses du Th pr ec edent, on suppose que le polyn ome L(x) ne change pas de signe sur [a, b], alors en
utilisant le Th de la moy` enne pour E(f), on obtient
E(f) =
1
(n + 1)!
f
(n+1)
()
Z
b
a
L(x)dx.
[a, b]
En utilisant ce dernier Th eor` eme on peut estimer les erreurs des m ethodes vues ci-dessus.
Pour la formule du rectangle on a :
E(f) = f

()
Z
b
a
(x a)dx = f

()
(b a)
2
2
[a, b]
cette m ethode est dordre 1.
Pour la formule du trap` eze on a :
E(f) =
1
2
f

()
Z
b
a
(x a)(x b)dx =
f

()
12
(b a)
3
la m ethode de Trap` eze est dordre 2.
Pour la formule de Simpson on a :
E(f) =
f
(4)
()
90

b a
2

5
,
la m ethode de Simpson est dordre 4.
Exemple :
I =
R
1
0
e
x
2
dx, a = 0,
a+b
2
=
1
2
, b = 1, f(0) = 1, f(
1
2
) = .7788, f(1) = .36788.
(1) Rectangle : I f(0) = 1.
(2) Trap` eze : I
h
f(0)+f(1)
2
i
= .68393.
(3) Simpson : I
1
6

f(0) + 4f(
1
2
) + f(1)

= .74718.
(4) La valeur de I ` a 5 d ecimales est .74718.
3.3. Formules compos ees. Plut ot que daugmenter le degr e du polyn ome dinterpolation, on peut obtenir une formule din-
tegration en d ecoupant lintervalle dint egration en sous-intervalles et en appliquant des formules simples sur chacun des sous-
intervalles.
Si n est entier, posons
h =
b a
n
, x
k
= a + kh, k = 0, ..., n.
alors
I(f) =
Z
b
a
f(x)dx =
n1
X
k=0
Z
x
k+1
x
k
f(x)dx

=
n1
X
k=0

f(x
k
) + f(x
k+1
)
2

h
h
3
12
f

(
k
)

,
o` u
k
[x
k
, x
k+1
] , k = 0, ..., n 1
D eveloppant et regroupant les termes qui apparaissent 2 fois, on obtient
I(f) =
h
2
"
f(a) + 2
n1
X
k=1
f(a + kh) + f(b)
#

h
3
12
n1
X
k=0
f

(
k
)
En appliquant le Th des valeurs interm ediaires, on peut r e ecrire lerreur sous la forme
E(f) =
nh
3
12
f

() =
(b a)
12
f

()h
2
.
Ceci nous donne la formule du trap` eze compos ee pour laquelle lapproximation est donn ee par :
Tn(f) =
h
2
"
f(a) + 2
n1
X
k=1
f(a + kh) + f(b)
#
xliv
et lerreur par
ET(f) =
(b a)
12
f

()h
2
.
Supposons que n soit pair, groupant les intervalles 2 ` a 2 et appliquant la formule de Simpson sur [xi, xi+2], on obtient
I(f) =
h
3
2
4
f(a) + 4
X
k impair
f(a + kh) + 2
X
k pair
f(a + kh) + f(b)
3
5

n
2
f
(4)
()
90
h
5
.
Ceci nous conduit ` a la formule de Simpson compos ee pour laquelle lapproximation est donn ee par
Sn(f) =
h
3
2
4
f(a) + 4
X
k impair
f(a + kh) + 2
X
k pair
f(a + kh) + f(b)
3
5
et lerreur par
ES(f) = f
(4)
()
(b a)
180
h
4
.
Exemple : D eterminer
R
1
0
e
x
2
dx.
Si n d esigne le nombre des intervalles utilis es.
n Tn(f) ET(f)
2 .73137 .015
4 .74298 3.84 10
3
8 .74658 9.58 10
4
16 .74676 1.39 10
4
32 .74680 5.98 10
5
Si nous d esirons obtenir 6 d ecimales exactes, il nous faut d eterminer h tel que
(1) max
01
[ f

() [
h
2
12
5 10
7
,
Pour une partition r eguli` ere x
k
= kh, h =
1
n
; donc nous cherchons n tel que
n
2

1
12
max
01
[ f

() [
1
5 10
7
.
or f

(x) = e
x
2
(4x
2
2) et f

(x) = e
x
2
4x(3 2x
2
). Puisque f

(x) ne change pas de signe sur [0; 1] ,

max
01
[ f

() [= max

[ f

(0) [, [ f

(1) [

= 2.
On voit que (1) sera satisfaite si
n
2

10
6
3
, n > 578.
Remarque Dans le choix de la pr ecision demand ee, il faut tenir compte des erreurs darrondi et de laccumutation des erreurs
xlv