Python for dummies : problème 2

Splines cubiques périodiques

L’interpolation par les splines cubiques usuelles sur l’intervalle [Xi−1 , Xi ] est définie par l’expression :

′′
Ui−1 U ′′
uhi (x) = (Xi − x)3 + i (x − Xi−1 )3
6hi 6hi

U ′′ hi
 
Ui−1
+ − i−1 (Xi − x)
hi 6

U ′′ hi
 
Ui
+ − i (x − Xi−1 )
hi 6

où Ui et Ui′′ sont les valeurs de la fonction et de la dérivée seconde en x = Xi . Finalement, on définit
les écarts hi = Xi − Xi−1 . Les dérivées secondes Ui′′ sont obtenues en résolvant une système linéaire
permettant d’obtenir une interpolation cubique par morceaux de classe C2 d’une fonction u(x) dont on
connaı̂t les ordonnées Ui aux abscisses Xi . Le système à résoudre s’écrit sous la forme

hi ′′ 2(hi + hi+1 ) ′′ hi+1 ′′ (Ui+1 − Ui ) (Ui − Ui−1 )

U + Ui + Ui+1 = −
6 i−1 6 6 hi+1 hi

i = 1, . . . , n − 1

Dans cette mission, nous allons calculer l’interpolation

cubique par morceaux d’une courbe fermée (x(t), y(t))
qui sera définie à partir de n points (Xk , Yk ) pour des
instants équidistants Tk = kh avec k = 0, . . . , n.

Attention, on ne dédouble pas le point initial de la

courbe que l’on souhaite construire !

Plus concrètement, il s’agit d’écrire une fonction

uh = splines(x,h,U) qui détermine les valeurs en x
de l’interpolation périodique cubique par morceaux pour
les valeurs Ui aux points Xi = ih. Le point initial et le
point final étant identiques pour une courbe fermée, il
faut évidemment légèrement modifier le système linéaire
pour obtenir une solution périodique, dont la dérivée
seconde, au point de départ, sera identique à celle au
point d’arrivée. Il est utile d’observer, par contre, qu’on
peut tirer profit du fait que les abscisses temporelles
sont toujours équidistantes. En d’autres mots, tous les
intervalles ont une longueur identique hi = h.
1. Comme python est très compliqué, nous vous donnons, exceptionnellement, un canevas contenant
déjà une grosse partie de la solution car elle calcule l’interpolation linéaire par morceaux....

def spline(x,h,U):
n = size(U); X = arange(0,n+1)*h

i = zeros(len(x),dtype=int)
for j in range(1,n):
i[X[j]<=x] = j

U = append(U,U[0])
return U[i]*(X[i+1]-x)/h + U[i+1]*(x-X[i])/h

Bien comprendre les instructions qui permettent d’identifier pour chaque élément x[index] du
tableau x l’intervalle i[index] à considérer !
Cela permet ensuite de bêtement calculer l’interpolation linéaire en prenant les deux points qui
encadrent cet intervalle. De manière vectorisée, tout cela se fait avec une unique ligne de code :-)
Et après, vous allez me dire que les petits programmes python sont trop compliqués à réaliser.
2. Ensuite pour tester, votre programme, on vous a vous fourni un tout petit programme simple
splineTest.py pour interpoler un cercle.

n = 4;
h = 3*pi/(2*(n+1));
T = arange(0,3*pi/2,h)
X = cos(T); Y = sin(T)

fig = plt.figure()
plt.plot(X,Y,’.r’,markersize=10)
t = linspace(0,2*pi,100)
plt.plot(cos(t),sin(t),’--r’)

t = linspace(0,3*pi/2,100)
plt.plot(spline(t,h,X),spline(t,h,Y),’-b’)
plt.axis("equal"); plt.axis("off")
plt.show()

3. Finalement, nous avons toujours le même programme un peu plus rigolo splineTestFun.py qui
permet de définir plein de points en cliquant sur la figure. Un double click permet d’obtenir la
courbe fermée obtenue par l’interpolation périodique cubique par morceaux.

4. Votre fonction (avec les éventuelles sous-fonctions que vous auriez créées) sera soumise via le site
web du cours.
Un petit programme interactif avec matplotlib...

Pour les petits curieux, voici l’implémentation de splineTestFun.py pour définir de manière magique
une courbe à interpoler...

import matplotlib
from matplotlib import pyplot as plt
from numpy import *
from splineTest import spline

# ====================== callback pour les événements avec la souris ======

#
# Observer la gestion distincte du clic simple et double :-)
# Apres un evenement, on redessine la figure avec draw()
#

def mouse(event):
global X,Y,n
if (event.dblclick):
t = arange(0,n+0.001,0.001)
x = spline(t,1.0,X)
y = spline(t,1.0,Y)
plt.plot(x,y,’-b’)
X,Y = [],[]; n = 0
else :
x = event.xdata
y = event.ydata
if (x != None and y != None) :
n = n + 1
X = append(X,[x])
Y = append(Y,[y])
print("New data : " + str(x) + "," + str(y))
plt.plot([x],[y],’.r’,markersize=10)
fig.canvas.draw()

# ============================= mainProgram ===============================

matplotlib.rcParams[’toolbar’] = ’None’
matplotlib.rcParams[’lines.linewidth’] = 1
plt.rcParams[’figure.facecolor’] = ’lavender’

X,Y = [],[]; n = 0
fig = plt.figure("Cubic spline interpolation")
fig.canvas.mpl_connect(’button_press_event’,mouse)
plt.ylim((0,1)); plt.xlim((0,1.3)); plt.axis("off")

plt.show()
Monsieur, j’arrive pas à utiliser le programme splineTestFun.py avec spyder

Tout d’abord, c’est mieux d’exécuter ton programme directement sur la console en introduisant
python splineTestFun.py ou d’utiliser l’environnement conseillé par les gentils informaticiens thonny...
Mais, comme nous sommes ouverts à toutes les options, j’ai donc essayé d’exécuter mon programme avec
spyder. Et cela marche, si, si !

Il y a, toutefois, une petite astuce : il faut modifier un paramètre dans les préférences de spyder. Plus
précisément, dans l’onglet Tools/Outils, il faut sélectionner Préférences IPython Console et modifier la
valeur de Backend = Automatic au lieu de Backend = Inline.
Et d’un coup, toutes vos difficultés disparaissent !
Et comment, j’ai trouvé cela, c’est simple : je me débrouille et je cherche la solution. C’est cela le job
d’un ingénieur, trouver tout seul des solutions et ne pas rester bêtement bloqué !

Là, j’ai mis 30 secondes à résoudre le problème via mon ami google et le site merveilleux stackoverflow:
so easy :-) Eh oui, il faut trouver comme un grand la solution de votre problème et pas rester bêtement
coincé : -)