ASI 3 Méthodes numériques pour l’ingénieur Interpolation f(x)
1
Approximation de fonctions • Soit une fonction f (inconnue explicitement) – connue seulement en certains points x0,x1…xn – ou évaluable par un calcul coûteux.
• Principe : – représenter f par une fonction simple, facile à évaluer
• Problème : – il existe une infinité de solutions ! 2 1 0 -1 -2 -3 -4 0
0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
1.8
2
Approximation de fonctions • Il faut se restreindre à une famille de fonctions – polynômes, – exponentielles, – fonctions trigonométriques…
2 1 0 -1 -2 -3 -4
0
0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
1.8
3
Quelques méthodes d'approximation • Interpolation polynomiale – polynômes de degré au plus n • polynômes de Lagrange • différences finies de Newton
• Interpolation par splines – polynômes par morceaux
• Interpolation d'Hermite – informations sur les dérivées de la fonction à approcher
• ...voir le groupe de TT… 4
Théorème d’approximation de Weierstrass soit f une fonction définie et continue sur l' intervalle [a, b] Alors, ∀ε > 0, il existe un polynôme P( x), définit sur [a, b] tel que : f ( x) − P ( x) < ε
∀ x ∈ [a, b]
ε )+ε ff(x)+(x P(x) P (x) f(x) f (x)
f (x) − ε
f(x)-ε
plus ε , est petit, plus l' ordre du polynome est grand
ab
ab
Interpolation : n + 1 points, n + 1 contraintes, n + 1 équations, n + 1 inconnues : ordre n5
Interpolation polynomiale • Le problème : les données, la solution recherchée
( x0 , y0 = f ( x0 ) ),..., ( xi , yi = f ( xi ) ),..., ( xi , yi = f ( xi ) ) P( x) tel que P ( xi ) = f ( xi ), i = 0, n
• mauvaise nsolution : résoudre le système linéaire P( x) = ∑ ai x i ⇒ Va = y V : matrice de Vandermonde i =0
• la combinaison linéaire de polynômes est un polynôme P( x) = y0 P0 ( x) + ... + yi Pi ( x) + yn Pn ( x) tel que Pi ( xi ) = 1 et Pi ( x j ) = 0
j≠i
ainsi P( xi ) = y0 P0 ( xi ) + ... + yi Pi ( xi ) + yn Pn ( xi ) ↓ ↓ ↓ 0 1 0
6
Interpolation polynomiale : Lagrange • Théorème – Soient n+1 points distincts xi réels et n+1 réels yi, il existe un unique polynôme p ∈ Pn tel que p(xi) = yi pour i = 0 à n démonstration
– Construction de p :
n
p( x ) = ∑ y i Li ( x ) i =0
(x − x ) L (x)=∏ (x − x ) n
avec Li polynôme de Lagrange – Propriétés de Li
j
i
j =0 j ≠i
i
j
L est un polynôme d’ordre n
• Li(xi)=1 • Li(xj)=0 (j ≠ i) 7
Lagrange : exemple n°1 • Exemple avec n=1 – on connaît 2 points (x0,y0) et (x1,y1) – on cherche la droite y=ax+b (polynôme de degré 1) qui passe par les 2 points : • y0 = a x0 + b • y1 = a x1 + b
–
a = (y0 - y1) / (x0 - x1) b = (x0 y1 - x1 y0) / (x0 - x1)
y1 y0
y0 − y1 x 0 y 1 − x1 y 0 y= x+ x0 − x 1 x0 − x1
x0
x1
x − x0 x − x0 x − x1 x − x1 − y1 = y0 + y1 y = y0 x0 − x1 x 0 − x1 x0 − x1 x1 − x0 L0(x)
L1(x) 8
Lagrange : exemple n°2 • Exemple avec n=2 – on connaît 3 points (0,1), (2,5) et (4,17) – polynômes de Lagrange associés : ( x( x − 2 ) x( x − 4 ) x − 2 )(x − 4 ) L2 ( x ) = L1 ( x ) = L0 ( x ) = 8 8 −4
-
0
2
4
0
2
4
0
2
4
9
Lagrange : exemple n°2 – calcul du polynôme d'interpolation 30
25
• p(x)=L0(x) + 5 L1(x) + 17 L2(x)
20
15
10
5
0 -1
0
1
2
3
4
5
• en simplifiant, on trouve p(x)=x2+1
10
Lagrange : l’algorithme Fonction y = lagrange(x,xi,yi) pour i = 1 jusqu' à n pour j = 1 jusqu' à n, j ≠ i; x − xi ( j ) l ←l* xi (i ) − xi ( j ) fait y ← y + yi * l fait
Complexité du calcul : n2 11
Lagrange : exemple n°3 • Exemple avec n=2
(fonction à approcher y=ex) – on connaît 3 points (0,1), (2,7.3891) et (4,54.5982) – Polynôme d'interpolation • p(x) =L0(x) + 7.3891 L1(x) + 54.5982 L2(x) 60
50
40
30
20
10
0
-10 -1
0
1
2
3
4
5
12
Lagrange : exemple n°3 – Erreur d'interpolation e(x) = f(x) - p(x) 60
50
40
30
20
10
0
-10
-20 -1
0
1
2
3
4
5
13
Lagrange : erreur d'interpolation • Théorème : – si f est n+1 dérivable sur [a,b], ∀ x ∈ [a,b], notons : • I le plus petit intervalle fermé contenant x et les xi • φ(x)=(x-x0)(x-x1)…(x-xn) ( n +1 )
– alors, il existe ξ ∈ I tel que
e( x ) =
(ξ ) φ (x ) (n + 1)!
f
– NB : ξ dépend de x
• Utilité = on contrôle l’erreur d’approximation donc la qualité de l’approximation
. 14
Lagrange : choix de n • Supposons que l'on possède un nb élevé de points pour approcher f … faut-il tous les utiliser ? – (calculs lourds)
• Méthode de Neville : – on augmente progressivement n – on calcule des Li de manière récursive – on arrête dès que l'erreur est inférieure à un seuil (d’ou l’utilité du calcul de l’erreur)
15
La méthode de Neuville • Définition Pm1 ,m2 ,...,mk ( x) polynôme de Lagrange calculé sur
(xm , ym ), (xm , ym ),..., (xm , ym ) Théorème (x − x j )P0,1,..., j −1, j +1,...,n ( x) − ( x − xi )P0,1,...,i−1,i+1,...,n ( x) P( x) = les k points
•
1
1
2
2
k
k
xi − x j
• Démonstration
P( xi ) = f ( xi ); P( x j ) = f ( x j ) et P( xk ) = f ( xk )
• Application systématique x0 x1 x2 x3
P0 = Q0,0 P1 = Q1,0 P2 = Q2,0 P3 = Q3,0
P0,1 = Q1,1 P1, 2 = Q2,1 P2,3 = Q3,1
Qi , j = Pi − j ,i − j +1,...,i −1,i
P0,1, 2 = Q2, 2 P1, 2,3 = Q3, 2
P0,1, 2, 4 = Q3,3
16
L’algorithme de Neuville Fonction y = Neuville(x,xi,yi) pour i = 1 jusqu' à n Q (i,0) ← yi (i ) fait pour i = 1 jusqu' à n pour j = 1 jusqu' à i ( x − xi (i − j ) )Q(i,j − 1 ) − ( x − xi (i ) )Q(i − 1,j − 1 ) Q(i,j) ← xi (i ) − xi (i − j ) fait y ← Q (n, n) fait Complexité du calcul : n2
17
Interpolation polynomiale : Newton • Polynômes de Newton : – base = {1, (x-x0), (x-x0)(x-x1), …, (x-x0)(x-x1)…(x-xn-1)} – on peut ré-écrire p(x) : p(x)=a0 + a1(x-x0) + a2(x-x0)(x-x1)+…+ an(x-x0)(x-x1)…(x-xn-1) – calcul des ak : méthode des différences divisées
18
Newton : différences divisées • Définition : – Soit une fonction f dont on connaît les valeurs en des points distincts a, b, c, … – On appelle différence divisée d’ordre 0, 1, 2,...,n les expressions définies par récurrence sur l’ordre k : – k=0 f [a] = f (a) – k=1 f [a,b] = ( f [b] - f [a] ) / ( b - a ) – k=2 f [a,b,c] = ( f [a,c] - f [a,b] ) / ( c - b ) – …
f [X,a,b] = ( f [X,b] - f [X,a] ) / ( b - a ) a∉X, b∉X, a≠b
19
Newton : différences divisées • Théorèmes : – détermination des coefficients de p(x) dans la base de Newton : f [x0, x1,…, xk] = ak
avec k = 0 … n
– erreur d'interpolation : e(x) = f [x0, x1,…, xn, x] φ(x)
20
Newton : différences divisées • Calcul pratique des coefficients : a0
x0
f [x0]
x1 f [x1] x2 f [x2] … … … … xn
f [xn]
a1
f [x0, x1] f [x1, x2] … …
f [x0, x1, x2] ... f [xn-3, xn-2, xn-1]
f [xn-1, xn]
f [xn-2, xn-1, xn]
a2
an
…
f [x0, ..., xn]
21
Newton : exemple • (ex. n°2) : n=2
(0,1), (2,5) et (4,17)
a0
0
f [x0]=1
2
f [x1]=5
4
f [x2]=17
a1 f [x0, x1] =(1-5)/(0-2)=2 f [x1, x2] =(5-17)/(2-4)=6
p(x)=1 + 2x + x(x-2)
f [x0, x1, x2] a2 = (2-6)/(0-4)=1
(et on retombe sur p(x) = 1 + x2) 22
Newton : l’algorithme Fonction a = Newton(xi,yi) pour i = 1 jusqu' à n F (i,0) ← yi (i ) fait pour i = 1 jusqu' à n pour j = 1 jusqu' à i F(i,j − 1 ) − F(i − 1,j − 1 ) F(i,j) ← xi (i ) − xi (i − j ) fait fait pour i = 1 jusqu' à n a(i) ← F (n, i ) fait Complexité du calcul : n2
23
A bas les polynômes • Ex : y=2(1+tanh(x)) - x/10
avec 9 points
– entre les points, le polynôme fait ce qu'il veut… et plus son degré est élevé plus il est apte à osciller ! 6
4
2
0
-2
-4
-6 -6
-4
-2
0
2
4
6
24
Interpolation par splines cubiques • Principe : – on approche la courbe par morceaux (localement) – on prend des polynômes de degré faible (3) pour éviter les oscillations
25
Splines cubiques : définition • Définition : – On appelle spline cubique d’interpolation une fonction notée g, qui vérifie les propriétés suivantes : • g ∈ C2[a;b] (g est deux fois continûment dérivable), • g coïncide sur chaque intervalle [xi; xi+1] avec un polynôme de degré inférieur ou égal à 3, • g(xi) = yi pour i = 0 … n
• Remarque : – Il faut des conditions supplémentaires pour définir la spline d’interpolation de façon unique – Ex. de conditions supplémentaires : • g"(a) = g"(b) = 0
spline naturelle. 26
Splines : illustration P2(x)=a2 (x-x2) 3+b2 (x-x2) 2+c2 (x-x2) +d2
P1(x)=α1x3+β1x2+χ1x+δ1 =a1 (x-x1)3+b1 (x-x1) 2+c1 (x-x1) +d1
27
Splines cubiques : détermination • Détermination de la spline d’interpolation – g coïncide sur chaque intervalle [xi; xi+1] avec un polynôme de degré inférieur ou égal à 3 g" est de degré 1 et est déterminé par 2 valeurs: • mi = g"(xi) et mi+1 = g"(xi+1)
(moment au noeud n°i)
– Notations : • hi = xi+1 - xi pour i = 0 … n-1 • δi= [xi; xi+1] • gi(x) le polynôme de degré 3 qui coïncide avec g sur l’intervalle δi
28
Splines cubiques : détermination – g"i(x) est linéaire : • ∀ x ∈ δi
• on intègre (ai constante) • on continue (bi constante)
– gi(xi) = yi – gi(xi+1) = yi+1
x − xi xi + 1 − x + mi hi hi ( ( xi +1 − x )2 x − xi )2 g i′ (x ) = mi +1 − mi + ai 2hi 2hi
g i′′(x ) = mi +1
g i (x ) = mi +1
(x − xi )3 + m (xi +1 − x )3 + a (x − x ) + b i
6 hi
mi hi 2 yi = + bi 6
1 yi +1
6 hi
i
mi +1 hi 2 = + ai hi + bi 6
i
i
2
29
Splines cubiques : détermination – g'(x) est continue : – 1 et 2
g i′ (xi ) = −mi
hi hi −1 + ai = mi + ai −1 = g i′−1 (xi ) 3 2 2
hi 1 ai = ( yi +1 − yi ) − (mi +1 − mi ) hi 6
– on remplace les ai dans : 3 4
hi −1mi −1 + 2(hi + hi −1 )mi + hi mi +1
1 1 ( yi − yi −1 ) = 6 ( y i +1 − yi ) − hi −1 hi
– Rappel : on cherche les mi (n+1 inconnues) on a seulement n-1 équations grâce à il faut rajouter 2 conditions, par exemple m0 = mn =0
4
(spline naturelle) 30
Splines cubiques : détermination 4
hi −1mi −1 + 2(hi + hi −1 )mi + hi mi +1
1 1 ( yi − yi −1 ) = 6 ( y i +1 − yi ) − hi −1 hi
– Ex de résolution avec hi = xi+1 - xi constant : • mi −1 + 4 mi + mi +1 = • forme matricielle :
Tm=f
1
h
2
( y i −1 − 2 y i + y i + 1 ) =
4 1
fi
1 0 m1 f 1 4 1 l = l p p p 1 4 1 1 4 mn −1 f n −1
• T inversible (diagonale strictement dominante) 31
Splines cubiques : l’algorithme pour i = 2; n − 1 T (i, i ) ← 2(hi + hi −1 ) T (i, i − 1) ← hi −1 T (i, i + 1) ← 2hi yi +1 − yi yi − yi −1 f (i − 1) ← 6 − hi hi −1 fait T ← T (2:n-1,2:n-1) m ←T \ f m ← [0, m,0] pour i = 1; n − 1 h 1 a(i ) ← ( yi +1 − yi ) − i (mi +1 − mi ) 6 hi mi hi b(i ) ← y (i ) − 6 fait Complexité du calcul : complexité du solveur
32
Splines cubiques : exemple • Ex : y=2(1+tanh(x)) - x/10
avec 9 points
6
4
spline
2
0
-2
-4
-6 -6
polynôme de degré 8 -4
-2
0
2
4
6
33
Conclusion • Interpolation polynomiale – évaluer la fonction en un point : Polynôme de Lagrange -> méthode de Neville – compiler la fonction : Polynôme de Newton
• Interpolation polynomiale par morceau : splines – – – –
spline cubique d’interpolation spline cubique d’approximation (on régularise) b spline spline généralisée : splines gausiènnes (multidimensionelle)
• approximation - apprentissage 34
