2 - INTERPOLATION SPLINE
J-P Croisille
Université Paul Verlaine-Metz
Semestre S7, master de mathématiques M1, année 2008/2009
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
1- INTRODUCTION
Fonctions splines: Fonctions interpolantes paticulièrement adaptées.
◮ Interpolation locale avec des polynômes de bas degré, mais
produisant des interpolations locales régulières.
◮ Principe: 2 paramètres:
• des points x0 < x1 < . . . < xn .
• un degré de régularité l tel que la spline soit dans C l−1 [a, b].
◮ En pratique, la spline cubique (qui est C 2 ) est très utilisée.
◮ Problèmes ouverts en théorie de l’interpolation spline en
plusieurs dimensions
◮ Quelques thèmes classiques: propriétés d’optimalité des splines,
zéros des splines.
◮ Importance des splines en CAO/CAD. Travaux de Paul de
Casteljau (1910-1999) chez Citroën et de Pierre Bézier
(1910-1999) chez Renault: design des pièces automobile.
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
2- INTERPOLATION POLYNOMIALE PAR
MORCEAUX
Idée générale: Il est plus efficace de réaliser une interpolation
polynômiale par morceaux, qu’une interpolation globale. Soit
s(x) , a ≤ x ≤ b, un polynôme défini par morceaux par
s(xj ) = f (xj ), j = 0, 1 . . . n (1)
où les valeurs x0 = a < x1 < . . . < xn = b , f (xj ), j = 0, 1 . . . n. Pour
x̄, s(x̄) est déterminé par quelques valeurs (xi , f (xi )) avec xi proche
de x̄.
Exemple 1: Interpolation linéaire par morceaux. Pour xj ≤ x ≤ xj+1 ,
x − xj
s(x) = f (xj ) + f (xj+1 − f (xj ) (2)
xj+1 − xj
Erreur: |f (x) − s(x)| ≤ 81 M2 h2 (cf Théorème B du ch.1). Attention à la
terminologie: “linéaire par morceaux” signifie que la fonction est
“d’ordre 1”, mais l’erreur est d’ordre 2.
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Spline cubique: s(x) est un polynôme cubique par morceaux.
s(x) = sj (x)|[xj ,xj+1 ] , j = 0, 1 . . . n − 1. (3)
Sur chaque intervalle, on les relations de collocation
sj (xj ) = f (xj ) , sj (xj+1 ) = f (xj+1 ) (4)
Il faut deux autres conditions pour déterminer un polynôme cubique
s(x) = a0 + a1 x + a2 x2 + a3 x3 ..
Méthode 1: On interpole deux valeurs de plus (élargissement du
stencil) (xj−1 , f (xj−1 )),(xj+2 , f (xj+2 )), 1 ≤ j ≤ n − 2.
Aux extrêmités de l’intervalle, on interpole
◮ f (x2 ), f (x3 ) pour le point j = 0
◮ f (xn−2 ), f (xn−3 ) pour le point j = n − 1.
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Méthode 2: On interpole de plus les dérivées en
(xj , s′ (xj )), (xj+1 , s′ (xj+1 )). On réalise donc une interpolation de type
Hermite. Les valeurs s′ (xj ) sont des valeurs approchées de la
dérivée, évaluées par ailleurs. On a
sj (x) = f (xj ) + s′ (xj )(x − xj ) + c2 (x − xj )2 + c3 (x − xj )3 (5)
Le système linéaire de 2 équations à 2 inconnues
sj (xj+1 ) = f (xj+1 )
(6)
s′j (xj+1 ) = s′ (xj+1 )
donne la solution
3(f (xj+1 ) − f (xj )) 2s′ (xj ) + s′ (xj+1 )
c2 =
−
(xj+1 − xj )2 xj+1 − xj
′ ′ (7)
2(f (xj ) − f (xj+1 )) s (xj ) + s (xj+1 )
c3 =
+
(xj+1 − xj )3 (xj+1 − xj )2
Pour la méthode 1, en appliquant le Théorème B du Ch.1, on obtient
une erreur d’ordre 4
j+2
1 Y
|f (x) − s(x)| ≤ |x − xi | max |f (4) (ξ)| (8)
24 i=j−1 xj−1 ≤ξ≤xj+2
4
splineée par 2 = 16.
Si le nombre de points
Jean-Pierre CROISILLE -double, l’erreur
Laboratoire LMAM est divis
2-Interpolation
Méthode 3: Spline cubique On effectue une interpolation par
morceaux, de sorte que s′′ , a ≤ x ≤ b soit continue. Sur [xj , xj+1 ], on
a
sj (x) = f (xj ) + s′ (xj )(x − xj ) + c2 (x − xj )2 + c3 (x − xj )3 (9)
donc
s′′ (x) = 2c2 + 6c3 (xj+1 − xj ) (10)
f (xj ) − f (xj+1 ) 2s′ (xj ) + 4s′ (xj+1 )
= 6 + (11)
(xj+1 − xj )2 xj+1 − xj
Ceci doit être identique à s′′j+1 (xj+1 ) On obtient l’équation
s′ (xk−1 ) + 2s′ (xk ) 2s′ (xk ) + s′ (xk+1 )
+ = (12)
xk − xk+1 xk+1 − xk
f (xk ) − f (xk−1 ) f (xk+1 ) − f (xk−1 )
3 +3 (13)
(xk − xk−1 )2 (xk+1 − xk )2
La dépendance de sj (x) = s|[xj ,xj+1 ] est non locale.
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Dans le cas où xk+1 = xk + h, ceci devient pour 1 ≤ k ≤ n − 1
3
s′ (xk−1 ) + 4s′ (xk ) + s′ (xk+1 ) = (f (xk+1 − f (xk )) (14)
h
Il faut fixer s′ (x0 ), s′ (xn ) par une procédure indépendante, par
exemple
f (x1 ) − f (x0 ) f (xn ) − f (xn−1 )
s′ (x0 ) ≃ , s′ (xn ) ≃ (15)
h h
On obtient un système linéaire en les inconnues s′ (xj ) , j = 0, . . . , n,
qui donne après résolution s′ (xj ) en fonction de f (xj ), ce qui permet
ensuite de calculer s(x) sur chaque [xj , xj+1 ].
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
2- SPLINES
Espace des splines: Ωn = (xν )ν=1,...,n avec
a = x0 < x1 < ... < xn = b (16)
réalise une partition d’un intervalle [a, b]. L’espace Sl (x0 , x1 , ..., xn )
des splines polynômiales de degré l est l’espace vectoriel des
fonctions s(x) ∈ Vect(1, x, . . . , xl , (x − x1 )l+ , . . . , (x − xn−1 )l+ ) ,c-a-d.
l
X n−1
X
s(x) = aλ xλ + bν (x − xν )l+ (17)
λ=0 ν=1
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
On a dim Sl = n + l. Les fonctions
(x − xν )l , x ≥ xν
ql,ν (x) = (x − xν )l+ = (18)
0
sont dans C l−1 [a, b] et
dim Sl = dim Pl (a, b) + dim Vect((x − xν )l+ ) = l + 1 + n− 1 = n+ l (19)
| {z } | {z }
l+1 n−1
Les fonctions (x − xν )l+ permettent d’assurer un raccord C l−1 aux
points xν des différents morceaux de polynômes.
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Theorème A (définition des splines)
L’espace des splines Sl (x0 , x1 , ..., xn ) est exactement l’espace Sl des
fonctions s : [a, b] → R tel que
◮ s ∈ C l−1 [a, b]
◮ s ∈ Pl [xν , xν+1 ]
Démonstration:
On montre que toute fonction s ∈ Sl possède une représentation de
la forme
X l n−1
X
s(x) = aλ xλ + bν (x − xν )l+ (20)
λ=0 ν=1
c-a-d., que Vect(xλ , (x − xν )l+ ) est génératrice pour Sl .
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Soit s(x) ∈ Sl . Montrons par récurrence la propriété (Pk ) suivante:
(Pk ): Il existe b1 , . . . , bk−1 t.q.
l
X k−1
X
s(x) = aλ xλ + bν (x − xν )l+ , sur [x0 , xk ] (21)
λ=0 ν=1
• k = 1. On regarde s(x) sur [x0 , x1 ]. C’est un polynôme de degré l,
donc
s(x) = a0 + a1 x + . . . + al xl (22)
Donc (21) est vrai pour k = 1.
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
• (Pk ) ⇒ (Pk+1 ) Soit
l
X k−1
X
λ
ρ(x) = s(x) − aλ x − bν (x − xν )l+ (23)
λ=0 ν=0
On a ρ(x) ∈ C l−1 [x0 , xk+1 ] et ρ(x) ≡ 0 sur [x0 , xk ], par hypothèse de
récurrence. De plus ρ(x) est solution de l’équation différentielle
y (l+1) (x) = 0 , xk ≤ x ≤ xk=1
(24)
y(xk ) = y ′ (xk ) = . . . y (l−1) (xk ) = 0
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
(24) est une [Link]. linéaire en l’inconnue
x 7→ Y (x) = [y(x), y ′ (x), . . . , y l (x)]T ∈ Rl+1 (25)
d
dx Y(x) = AY (x)
(26)
Y (xk ) = [0, . . . , α]T
α
Solution unique ρ(x) = l! (x − xk )l+ . En notant bk = −α/l!, s(x) s’écrit
sur [xk , xk+1 ]
l
X k
X
s(x) = aλ xλ + bν (x − xν )l+ (27)
λ=0 ν=1
ce qui prouve (Pk ) ⇒ (Pk+1 ). On a donc montré que la famille de
fonctions (xλ , (x − xν )l+ ) est génératrice pour Sl , c-a-d,
Sl ⊂ Sl (x0 , x1 , . . . , xn ). Réciproquement toute spline
s ∈ Sl (x0 , . . . , xn ) est dans Pl [xν , xν+1 ] pour tout ν et est dans
C l−1 [a, b], donc Sl (x0 , . . . , xn ) ⊂ Sl . Donc Sl = Sl (x0 , . . . , xn ).
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Conditions de bord pour les splines: On a vu que la spline cubique
s(x) ∈ S3 (x0 , x1 , . . . , xn ) est déterminée par les valeurs des dérivées
aux noeuds s′ (xj ) solution du système linéaire
3
s′ (xk−1 ) + 4s′ (xk ) + s′ (xk+1 ) =(f (xk+1 − f (xk−1 )) , 1 ≤ k ≤ n − 1
h
(28)
Il faut rajouter 2 conditions de bord en a = x0 et b = xn . Trois types
de conditions de bord sont usuelles
◮ Les conditions de Hermite
s′ (a) = f ′ (a) , s′ (b) = f ′ (b) (29)
Cela suppose que f ′ (a), f ′ (b) sont connues par ailleurs
◮ Les conditions “naturelles”
s′′ (a) = s′′ (b) = 0 (30)
◮ Les conditions périodiques
s′ (a) = s′ (b) , s′′ (a) = s′′ (b) (31)
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Theorème B (relation intégrale pour les splines cubiques)
Soit d(x) = f (x) − s(x), l’erreur entre la fonction f (x) et la spline s(x).
On a la condition limite
s′′ (a)d′ (a) = s′′ (b)d′ (b) (32)
pour la spline cubique s(x) si et seulement si on a la relation intégrale
Z b Z b Z b
f ′′ (x)2 dx = (f ′′ (x) − s′′ (x))2 dx + s′′ (x)2 dx (33)
a a a
Démonstration:
La relation intégrale (33) équivaut à
Z b
(f ′′ (x)s′′ (x) − s′′ (x)2 )dx = 0 (34)
a
soit Z b
s′′ (x)d′′ (x)dx = 0 (35)
a
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
En intégrant deux fois par parties, on obtient
Z b Z b
′′ ′′
s (x)d (x)dx = s′′ d′ |ba − s′′′ (x)d′ (x)dx
a a
n−1
X Z xν+1
= s′′ (b)d′ (b) − s′′ (a)d′ (a) − s′′′ (x)d′ (x)dx
ν=0 xν
n−1
X Z xν+1
= s′′ (b)d′ (b) − s′′ (a)d′ (a) + ( s(4) (x)d(x)dx − s(3) (x)d(x)|xxν+1
ν
)
ν=0 xν
Mais on a
◮ d(xν ) = 0, donc s(3) (x)d(x)|xxν+1
ν = 0.
◮ s(4) (x)|[xν−1 ,xν ] ≡ 0 car s(x) ∈ P3 [a, b].
On a finalement
Z b
s′′ (x)d′′ (x)dx = s′′ (b)d′′ (b) − s′′ (a)d′′ (a) (36)
a
ce qui donne le résultat.
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Chacune des 3 conditions de bord
◮ Hermite
◮ Naturelle
◮ Périodique
vérifie cette relation.
Corollaire (Interprétation mécanique de la spline cubique)
Soit s(x) ∈ S3 (x0 , x1 , . . . , xn ), la spline cubique interpolant
f (x) ∈ C 2 [a, b]. Alors s(x) est la solution du problème de minimisation
Z b
min g ′′ (x)2 dx (37)
g∈G a
où G est l’ensemble des fonctions C 2 [a, b] vérifiant les mêmes
conditions d’interpolation que s(x)
◮ g(xν ) = f (xν ), ; , ν = 0, 1, . . . , n
◮ l’une des trois conditions limite
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Démonstration:
La relation (33) est valable non seulement pour f , mais aussi avec
toute fonction g ∈ G. Donc
Z b Z b Z b
(g ′′ (x))2 dx = (g ′′ (x) − s′′ (x))2 dx + s′′ (x)2 dx (38)
a a a
Z b
≥ s′′ (x)2 dx (39)
a
(40)
donc Z Z
b b
′′ 2
min g (x) dx = s′′ (x)2 dx (41)
g∈G a a
Corollaire
Le problème (M) possède une unique solution. Autrement dit, il existe
une unique spline cubique s(x) ∈ S3 [a, b] solution du problème
d’interpolation (avec l’un des trois choix de conditions limite).
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Démonstration:
Si s1 (x), s2 (x) ∈ G alors
Z b
(s′′1 (x) − s′′2 (x))2 dx = 0 (42)
a
Donc s′′1 (x) = s′′2 (x), ce qui entraı̂ne que s1 et s2 diffèrent par une
fonction affine. Les conditions limite entraı̂nent que s1 ≡ s2 .
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Exemples:
bleu: Log(x), vert: interpolee spline cubique
1
−1
−2
−3
−4
−5
0 0.2 0.4 0.6 0.8 1 1.2 1.4
Interpolation spline avec 4 points de:
x ∈ [0.01, 1.01] 7→ Log(x) (43)
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
bleu: Log(x), vert: interpolee spline cubique
1
−1
−2
−3
−4
−5
0 0.2 0.4 0.6 0.8 1 1.2 1.4
Interpolation spline avec 10 points de:
x ∈ [0.01, 1.01] 7→ Log(x) (44)
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
bleu: Log(x), vert: interpolee spline cubique
1
−1
−2
−3
−4
−5
0 0.2 0.4 0.6 0.8 1 1.2 1.4
Interpolation spline avec 20 points de:
x ∈ [0.01, 1.01] 7→ Log(x) (45)
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
bleu: sin(2*pi*x), vert: interpolee spline cubique
1
0.8
0.6
0.4
0.2
−0.2
−0.4
−0.6
−0.8
−1
0 1 2 3 4 5 6 7
Interpolation spline avec 3 points de:
x ∈ [0., 2π] 7→ sin(4x) (46)
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
bleu: sin(2*pi*x), vert: interpolee spline cubique
1.5
0.5
−0.5
−1
−1.5
0 1 2 3 4 5 6 7
Interpolation spline avec 7 points de:
x ∈ [0., 2π] 7→ sin(4x) (47)
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
bleu: Log(x), vert: interpolee spline cubique
1.5
0.5
−0.5
−1
−1.5
0 1 2 3 4 5 6 7
Interpolation spline avec 15 points de:
x ∈ [0., 2π] 7→ sin(4x) (48)
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
bleu: Log(x), vert: interpolee spline cubique
1
0.8
0.6
0.4
0.2
−0.2
−0.4
−0.6
−0.8
−1
0 1 2 3 4 5 6 7
Interpolation spline avec 31 points de:
x ∈ [0., 2π] 7→ sin(4x) (49)
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
1
0.8
0.6
0.4
0.2
−0.2
−0.4
−5 −4 −3 −2 −1 0 1 2 3 4 5
Interpolation spline avec 5 points de la fonction de Runge:
1
x ∈ [−5, 5] 7→ (50)
1 + x2
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
1
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
−5 −4 −3 −2 −1 0 1 2 3 4 5
Interpolation spline avec 10 points de la fonction de Runge:
1
x ∈ [−5, 5] 7→ (51)
1 + x2
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
1
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
−5 −4 −3 −2 −1 0 1 2 3 4 5
Interpolation spline avec 20 points de la fonction de Runge:
1
x ∈ [−5, 5] 7→ (52)
1 + x2
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
Tableau de mesures
ν E(ν) ν E(ν) ν E(ν)
>0.00 -3.10 3.65 3.10 >5.88 4.77
0.29 -2.94 >4.00 4.19 6.25 5.02
>0.45 -2.72 4.17 4.90 >6.71 5.05
>0.80 -2.23 >4.35 5.77 7.18 5.39
1.11 -1.60 >4.57 6.57 >8.00 6.55
>1.43 -0.78 >4.76 6.23 8.50 7.45
1.82 0.00 5.00 5.52 >9.00 8.45
>2.27 1.00 >5.26 4.90 9.50 9.80
2.50 1.30 5.56 4.65 >10.00 11.30
>2.91 1.80
On calcule l’interpolé spline (cubique) basée sur les points marqués
> et on compare pour les autres points valeur mesurée/valeur prédite
par la spline.
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
12
10
−2
−4
0 1 2 3 4 5 6 7 8 9 10
Interpolation spline avec les mesures marquées >, (15 mesures).
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline
12
10
−2
−4
0 1 2 3 4 5 6 7 8 9 10
Interpolation spline avec seulement 9 mesures .
Jean-Pierre CROISILLE - Laboratoire LMAM 2-Interpolation spline