|
|
import numpy as np
|
|
|
np.set_printoptions(suppress=True) # (pour mieux arrondir à 0 lors des print)
|
|
|
from numpy.polynomial.polynomial import polyval
|
|
|
|
|
|
import math
|
|
|
|
|
|
import matplotlib.pyplot as plt
|
|
|
|
|
|
# On active le "mode interactif" de pyplot. Ainsi, les plt.show() ne sont plus nécessaires.
|
|
|
plt.ion()
|
|
|
|
|
|
# --------------------------------------------------------------------------------------
|
|
|
# Fonction de "fin temporaire" pendant l'avancée du TP. Attend une frappe Entrée, puis quitte le script.
|
|
|
|
|
|
def stop_ici():
|
|
|
plt.pause(0.1) # nécessaire pour que matplotlib ait le temps d'afficher les figures
|
|
|
input()
|
|
|
exit()
|
|
|
|
|
|
# --------------------------------------------------------------------------------------
|
|
|
# Tous les passages indiqués "TODO()" sont à remplacer par vos soins
|
|
|
|
|
|
def TODO():
|
|
|
print("à vous!")
|
|
|
stop_ici()
|
|
|
|
|
|
# --------------------------------------------------------------------------------------
|
|
|
# Cette petite fonction sert à tracer les flèches avec une syntaxe un peu plus intuitive que le code matplotlib de base.
|
|
|
# - xdepart, ydepart : coordonnées du point d'attache de la flèche
|
|
|
# - xfleche, yfleche : longueurs de la flèche en x et en y
|
|
|
|
|
|
def myfleche(xdepart, ydepart, xfleche, yfleche):
|
|
|
plt.annotate("", (xdepart+xfleche,ydepart+yfleche) , (xdepart,ydepart),
|
|
|
arrowprops={'color':'red','shrink': 0,'width':0.02}, annotation_clip=False)
|
|
|
|
|
|
|
|
|
|
|
|
#################################################################################################
|
|
|
###
|
|
|
### Rappels des TP2 et TP3 : courbe paramétrée et approximation affine
|
|
|
###
|
|
|
#################################################################################################
|
|
|
|
|
|
|
|
|
# Une *courbe paramétrée polynômiale* est une courbe paramétrée de la forme
|
|
|
# (P(t), Q(t))
|
|
|
# dont les 2 fonctions de coordonnées sont polynômiales.
|
|
|
#
|
|
|
# - La variable "t" s'appelle le **paramètre** de la courbe.
|
|
|
# - Chaque valeur de t définit UN point, dont l'abscisse vaut P(t) et l'ordonnée vaut Q(t).
|
|
|
# - La courbe est obtenue en représentant les points associés à TOUTES les valeurs de t.
|
|
|
|
|
|
|
|
|
#################################################################################
|
|
|
# Approximation affine pour une courbe polynômiale
|
|
|
#################################################################################
|
|
|
#
|
|
|
# L'approximation affine de la courbe polynômiale au point t s'écrit
|
|
|
#
|
|
|
# P(t+h) ~= P(t) + h.P'(t) [en abscisses]
|
|
|
# Q(t+h) ~= Q(t) + h.Q'(t) [en ordonnées]
|
|
|
#
|
|
|
# On représente cette approximation affine par une FLÈCHE.
|
|
|
#
|
|
|
# - La flèche démarre du point "t" sur le graphe
|
|
|
# Point(t) = ( P(t) , Q(t) )
|
|
|
#
|
|
|
# - La flèche a des longueurs (en abscisses et en ordonnées) de valeurs respectives
|
|
|
# fleche = ( h.P'(t) , h.Q'(t) )
|
|
|
# (Concrètement, la valeur choisie pour h détermine donc la mise à l'échelle de la flèche.)
|
|
|
|
|
|
# Avec ces notations, l'approximation affine ci-dessus peut se réécrire
|
|
|
# Point(t+h) ~= Point(t) + fleche
|
|
|
# c'est-à-dire que "fleche" représente l'évolution (approximative) de la courbe
|
|
|
# lorsque le point d'intérêt passe du paramètre "t" à "t+h".
|
|
|
#################################################################################
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
#################################################################################################
|
|
|
###
|
|
|
### EXERCICE 5 : Interpolation d'Hermite en 2D
|
|
|
###
|
|
|
#################################################################################################
|
|
|
|
|
|
|
|
|
# Il s'agit exactement du même problème qu'à l'exercice 4, mais pour une courbe polynômiale en 2D.
|
|
|
# On impose un certain nombre de points de passage en 2D, chacun associé à une flèche (= dérivée).
|
|
|
# Le but est de trouver une courbe polynômiale qui passe par ces points, et dont la dérivée à ces points
|
|
|
# soit donnée par les flèches : c'est-à-dire que la courbe doit être tangente à chaque flèche, et évoluer
|
|
|
# d'autant plus "fort" que la flèche est grande.
|
|
|
#
|
|
|
#
|
|
|
# Pour illustrer cela, on a implémenté une petite GUI (fenêtre interactive).
|
|
|
# - les points de passage eux-mêmes sont fixés à une position imposée
|
|
|
# - par contre, les clicks utilisateur permettent de définir la flèche (=dérivée) associée à chaque point
|
|
|
# - au dernier click, la courbe interpolatrice est calculée puis affichée.
|
|
|
|
|
|
|
|
|
print('''Complétez le code de la GUI ci-dessous, afin de calculer et tracer la courbe interpolatrice.''')
|
|
|
|
|
|
|
|
|
plt.close('all') # ferme les fenêtres précédentes, pour y voir plus clair
|
|
|
|
|
|
##############################################################
|
|
|
# Nombre de points à interpoler (laisser à 2 pour l'instant)
|
|
|
N = 2
|
|
|
##############################################################
|
|
|
|
|
|
# Position des N points cible (fixes, le long d'un polygone régulier à N sommets) :
|
|
|
|
|
|
# Elles sont stockées dans le np.array "point", sous le format suivant :
|
|
|
# - point[0,k] renvoie l'abscisse du point cible numéro k
|
|
|
# - point[1,k] renvoie l'ordonnée du point cible numéro k
|
|
|
# - point[0,:] renvoie toutes les abscisses des points cible
|
|
|
# - point[1,:] renvoie toutes les ordonnées des points cible
|
|
|
|
|
|
theta = np.pi - np.arange(0, 2*np.pi, 2*np.pi/N)
|
|
|
point = np.array([np.cos(theta),np.sin(theta)])
|
|
|
|
|
|
# La figure + les points cible :
|
|
|
fig = plt.figure()
|
|
|
limits = (-2,2) # limites de représentation (en x et en y)
|
|
|
ax = fig.add_subplot(111, xlim=limits, ylim=limits)
|
|
|
ax.plot(*point,linestyle='none',marker='o',color='blue')
|
|
|
|
|
|
# Dérivées imposées à chaque point cible (vont être définies par les clicks successifs) :
|
|
|
|
|
|
deriv = np.zeros((2,N))
|
|
|
# Une fois définies, elles seront stockées dans le np.array "deriv", sous le format suivant :
|
|
|
# - deriv[0,k] renvoie l'abscisse de la dérivée cible numéro k
|
|
|
# - deriv[1,k] renvoie l'ordonnée de la dérivée cible numéro k
|
|
|
# - deriv[0,:] renvoie toutes les abscisses des dérivées cible
|
|
|
# - deriv[1,:] renvoie toutes les ordonnées des dérivées cible
|
|
|
|
|
|
# Fonction d'interpolation (recopiée de l'exercice 4)
|
|
|
|
|
|
def interpol_hermite(ys,ds):
|
|
|
TODO()
|
|
|
|
|
|
# La fonction de callback de la GUI (sera appelée à chaque click)
|
|
|
|
|
|
nclick = -1
|
|
|
def onclick(event):
|
|
|
global clickpos, nclick, deriv, ax
|
|
|
clickpos = np.array([event.xdata, event.ydata]) # position du click
|
|
|
nclick = (nclick+1) % (N+2) # incrémente nclick (modulo N+2)*
|
|
|
|
|
|
# Récupère l'objet "axe" associé à la figure. Les commandes graphiques devront être
|
|
|
# appelées directement sur cet objet : "ax.ma_commande(...)" et non pas "plt.ma_commande(...)"
|
|
|
ax = fig.get_axes()[0]
|
|
|
|
|
|
if nclick < N:
|
|
|
### Trace la flèche associée à chaque point
|
|
|
fleche = clickpos - point[:,nclick]
|
|
|
myfleche(*(point[:,nclick]),*fleche)
|
|
|
ax.set(xlim=limits, ylim=limits)
|
|
|
# Écart h de l'approximation affine (fixe la proportionnalité gobale des flèches)
|
|
|
h = 0.2
|
|
|
# stocke la valeur imposée de dérivée correspondante (! en divisant par h !)
|
|
|
deriv[:,nclick] = fleche/h
|
|
|
|
|
|
elif nclick == N:
|
|
|
### au (N+1)ème click, on calcule la courbe interpolatrice et on l'affiche
|
|
|
# --> À VOUS !
|
|
|
# Calcule les poolynômes P(X) et Q(X) par interpolation
|
|
|
p = TODO()
|
|
|
q = TODO()
|
|
|
# Discrétisation de l'intervalle [1,N] en 200 points
|
|
|
t = TODO()
|
|
|
# Trace la courbe
|
|
|
ax.plot(TODO())
|
|
|
|
|
|
else:
|
|
|
### au (N+2)ème click, on nettoie pour pouvoir recommencer le cycle
|
|
|
ax.clear()
|
|
|
ax.plot(*point,linestyle='none',marker='o',color='blue')
|
|
|
ax.set(xlim=limits, ylim=limits)
|
|
|
|
|
|
# La commande ci-dessous "active" la GUI, en connectant la figure "fig" à la fonction de callback "onclick"
|
|
|
cid = fig.canvas.mpl_connect('button_press_event', onclick)
|
|
|
|
|
|
|
|
|
|
|
|
stop_ici() # ---------- Supprimez cette ligne une fois que le code avant fonctionne ------------------
|
|
|
|
|
|
# --------------------------------------------------------------------------------------
|
|
|
print("-"*80)
|
|
|
|
|
|
|
|
|
|
|
|
# Jusqu'à présent, la fonction "interpol_hermite" que vous avez écrite
|
|
|
# n'est capable d'interpoler que 2 points (avec leurs 2 dérivées associées).
|
|
|
#
|
|
|
# Dans ces dernières questions du TP, on vous guide pour réécrire la fonction afin
|
|
|
# qu'elle interpole N points (avec leurs N dérivées associées), pour n'importe
|
|
|
# quelle valeur de N.
|
|
|
|
|
|
print('''(Question papier).
|
|
|
Écrivez la matrice "A" et le second membre "b" du système si vous deviez
|
|
|
interpoler 3 points (et leurs 3 dérivées associées) :
|
|
|
P(1) = y1
|
|
|
P(2) = y2
|
|
|
P(3) = y3
|
|
|
P'(1) = d1
|
|
|
P'(2) = d2
|
|
|
P'(3) = d3
|
|
|
par un polynôme P de degré 5 (c'est-à-dire, avec 6 coefficients inconnus).''')
|
|
|
|
|
|
stop_ici() # ---------- Supprimez cette ligne une fois que le code avant fonctionne ------------------
|
|
|
|
|
|
# --------------------------------------------------------------------------------------
|
|
|
print("-"*80)
|
|
|
|
|
|
|
|
|
print('''(Question papier).
|
|
|
Même question dans le cas général où vous devez interpoler N points (et leurs
|
|
|
N dérivées associées). Trouvez la *formule générale* donnant la valeur de
|
|
|
l'élément de la matrice A en position (i,j).
|
|
|
Indice :
|
|
|
- pour les N premières lignes de A (contrainte de passer par un point yi),
|
|
|
vous avez déjà trouvé cette formule à la séance précédente (partie 1,
|
|
|
interpolation de Lagrange).
|
|
|
- pour les N dernières lignes (contrainte de vérifier une dérivée di), à
|
|
|
vous de jouer!
|
|
|
''')
|
|
|
|
|
|
stop_ici() # ---------- Supprimez cette ligne une fois que le code avant fonctionne ------------------
|
|
|
|
|
|
# --------------------------------------------------------------------------------------
|
|
|
print("-"*80)
|
|
|
|
|
|
|
|
|
print('''Complétez la nouvelle fonction "interpol_hermite(ys,ds)" (ébauche de code ci-dessous),
|
|
|
dans laquelle "ys" et "ds" peuvent maintenant être des tableaux de n'importe quelle taille N.
|
|
|
|
|
|
Enfin testez votre fonction en modifiant le paramètre N de la GUI plus haut (celui qui
|
|
|
valait 2 jusqu'à présent).
|
|
|
''')
|
|
|
|
|
|
### Fonction d'interpolation. Renvoie l'unique polynôme P de degré 2N-1 tel que
|
|
|
# P(1) = y1
|
|
|
# ...
|
|
|
# P(N) = yN
|
|
|
# P'(1) = d1
|
|
|
# ...
|
|
|
# P'(N) = dN
|
|
|
|
|
|
def interpol_hermite(ys,ds):
|
|
|
N = len(ys)
|
|
|
# les N premières lignes de la matrice A
|
|
|
A1 = [ [ TODO() for j in range(2*N) ] for i in range(1,N+1) ]
|
|
|
# les N dernières lignes de la matrice A
|
|
|
A2 = [ [ TODO() for j in range(2*N) ] for i in range(1,N+1) ]
|
|
|
# la matrice A entière
|
|
|
A = TODO()
|
|
|
# le second membre du système linéaire à résoudre
|
|
|
b = TODO()
|
|
|
# la solution du système
|
|
|
p = TODO()
|
|
|
return p
|
|
|
|
|
|
# Remarque : pas besoin d'effacer l'ancienne fonction 'interpol_hermite(ys,ds)'
|
|
|
# que vous avez écrite plus haut. Au moment où vous définissez la nouvelle version
|
|
|
# interpol_hermite(ys,ds), Python la "remplace" purement et simplement (exactement
|
|
|
# comme lorsqu'on redéfinit une variable, du genre v=v+4).
|
|
|
# En Python, les fonctions ne sont rien d'autre que des "variables qu'on peut appeler".
|
|
|
#
|
|
|
# Et puisque vous clickerez dans la GUI *après* avoir défini à Python la nouvelle
|
|
|
# version, la GUI fera naturellement appel à la nouvelle version :)
|
|
|
|
|
|
|
|
|
stop_ici() # ---------- Supprimez cette ligne une fois que le code avant fonctionne ------------------
|
|
|
|
|
|
|
|
|
########################################
|
|
|
# Vous avez tout fini ??? Bravo !!!
|
|
|
# Une fois arrivé ici, si vous vous ennuyez, modifiez le code de la GUI précédente pour
|
|
|
# que les clicks permettent également de définir la *position* des N points de passage.
|
|
|
# Il y aura donc un total de 2*N clicks : une position + une flèche pour chaque point de passage.
|
|
|
# Vous obtiendrez ainsi une GUI permettant de tester la version la plus générale de
|
|
|
# l'interpolation de Hermite (positions+flèches).
|
|
|
|
|
|
input()
|
|
|
|
|
|
|