sylvain@sylbarth.com		Curriculum Vitae			My Blog

Interpolations polynomiales et cubic splines

Vous pouvez télécharger les fichiers interp.cc et interp.h qui vous permettrons d'expérimenter les méthodes présentées sur cette page.

 Interpolations de Lagrange (1735-1796)

Selon le théorême de Weiestrass' (et la preuve de Benstein) :

Si f(x) est continue sur l'intervalle [a,b] alors il quelque soit e>0 réel positif, il existe un entier n et un polynome pn(x) tels que |f(x)-pn(x)| < e

Supposons que vous disposez de N+1 points initiaux et que vous désirez effectuer une interpolation polynomiale sur(entre) ces points. La méthode de Lagrange propose de choisir N points parmis les N+1 et de résoudre le système de N équations à N inconnues formé par les équation polynomiales passant par ces points. Si vous disposez de 4 points, vous rechercherez donc un polynome de degré 3 passant par ces points.

double lagrange(const double& u, const double* x, const double* y, const int& n) { double interp = 0.0; for(int i=0; i < n; i++) { double p = 1.0; for(int j=0; j < n; j++) { if(i!=j) p *= (u-x[j])/(x[i]-x[j]); } interp += p * y[i]; } return interp; }

Le problème posé par cette méthode est qu'elle ne permet pas de choisir le degré du polynome. Un degré trop élevé engendrera des oscillations excessives et inversement. Les méthodes travaillant sur les différences et permettent de corriger ce défaut en choisissant un seuil de convergence "acceptable".

 Méthodes des différences

La méthode standard des différences peut donner le même résultat que la méthode de Lagrange, s'il n'est pas possible de trouver un degré polynomiale inférieur à N sans perte de précision excessive. Dans certains cas elle remplacera cependant avantageusement l'interpolation de Lagrange en ajustant un polynome de degré inférieur.

void divdiff_coef(double* dd, const double* x, const double* y, const int& n) { double temp1, temp2; for ( int i=0; i < n; i++ ) dd[i] = y[i]; for ( int j=1; j < n; j++ ) { temp1 = dd[j-1]; for ( int k=j; k < n; k++ ) { temp2 = dd[k]; dd[k] = (dd[k] - temp1)/(x[k] - x[k-j]); temp1 = temp2; } } } double divdiff(const double& u, const double* dd, const double* x, const int& n) { double interp = 0.0; for ( int i=(n-1); i>=1; i-- ) interp = (interp + dd[i]) * (u - x[i-1]); interp += dd[0]; return interp; }

 Cubic Splines

A la différence de la méthode précédente, le cubic spline n'est pas une méthode d'interpolation polynomiale à proprement parler. Nous cherchons plutot à lisser une série de points (un peu comme on le ferai avec une moyenne mobile) et non pas à trouver quel polynome représente le mieux ces points. Effectuer une interpolation par un cubic spline revient à résoudre un système tridiagonal. Comme pour la méthode des différences, on recherche tout d'abord les valeurs des coefficients avant de calculer les points interpolés.

L'algorithme utilisé pour rechercher un cubic spline sur une série de points est un peu plus compliqué que les deux précédents. Si vous être intéréssé par une implémentation en C++, vous pouvez downloader enl_interp.cc et enl_interp.h. Elle utilise la résolution d'une matrice tridiagonale proposée par les fameux Numerical Recipies in C.

 Comparaison

Pour démontrer les différences qui peuvent exister entre une interpolation polynomiale et par les splines, nous avons généré 9 points non équidistants de facon (pseudo) aléatoire. Le cubic spline donne des résultats assez differents de l'interpolation de Lagrange, avec un aspect bien plus lissé (courbe bleue).

Comme vous pourez le constater sur le tableau suivant, la méthode de Lagrange donne bien les mêmes résultats que la méthode des différences.