TP Cinématique et trajectoires¶

Marc Buffat, dpt mécanique, UCB Lyon 1

Attention il faut exécuter toutes les cellules suivantes pour afficher son nom !!

In [1]:
%matplotlib inline
import numpy as np
from matplotlib import pyplot as plt
import matplotlib.pyplot as plt
plt.rc('font', family='serif', size='14')
In [2]:
from validation.valide_markdown import test_markdown, test_code
from validation.validation import info_etudiant, bib_validation
from IPython.display import display, Markdown, clear_output# nom etudiant
def printmd(string):
    display(Markdown(string))
# test si numero étudiant spécifier
try: NUMERO_ETUDIANT
except NameError: NUMERO_ETUDIANT = None 
if NUMERO_ETUDIANT is None:
    NOM,PRENOM,NUMERO_ETUDIANT = info_etudiant()
printmd("**Etudiant** {} {}  id={}".format(NOM,PRENOM,NUMERO_ETUDIANT))
# generation des parametres
notebook_name = "TP_trajectoire.ipynb"
bib_validation('cours','MGC2028L')
from Circuit import Circuit,test_vitesse,test_trajectoire,test_acceleration
np.random.seed(NUMERO_ETUDIANT)
Circuit("circuit.txt");
np.random.seed(NUMERO_ETUDIANT)
_H_ = 200.0 + np.random.randint(50)

Etudiant NBGrader Course id=2592889824

Objectif¶

L'objectif du TP est d'apprendre à manipuler des tableaux numpy représentant des fonctions, pour pouvoir numériquement analyser ces fonctions, calculer numériquement la dérivée de ces fonctions ou pouvoir numériquement les intégrer.

Les 2 cas traités sont:

  • une étude de cinématique 1D: chute libre
  • une étude de cinématique 2D: cycliste sur piste

Etude cinématique en 1D: chute libre¶

Dans un fichier "vchute.dat" sont stockées les mesures de vitesse d'un objet en chute libre en fonction du temps, sous la forme de 2 colonnes: le temps en seconde et la vitesse en km/h séparée par un espace.

Le contenu du fichier est de la forme suivante

# temps et vitesse pour une chute libre
0.000 0.000
0.100 -3.529
0.200 -7.045
0.300 -10.533

  • Utiliser la fonction numpy loadtxt pour lire les données et les mettre dans les 2 variables T et U.

  • Vérifier que les 2 variables T et U sont bien des tableaux numpy de même dimension et tracer la courbe de vitesse en fonction du temps avec un titre et des labels sur les axes pour vérifier.

  • Dans la suite, on utilisera comme unité : s. pour le temps, km/h pour la vitesse, et m pour la distance.

In [3]:
# lecture des données
T = None
U = None
## BEGIN SOLUTION
T,U = np.loadtxt("vchute.dat",delimiter=' ',unpack=True)
print("Vérification ", T.size==U.size)
plt.plot(T,U,label="U")
plt.title("vitesse de Chute libre")
plt.xlabel('t en sec.')
plt.ylabel('vitesse en km/h')
plt.legend()
## END SOLUTION

Vérification  True
Out[3]:
<matplotlib.legend.Legend at 0x7f4050ee6be0>
In [4]:
# ne rien ecrire
## BEGIN HIDDEN TESTS
assert(T.size == U.size)
assert(test_code(notebook_name,"cell-verif01",
                 ["plt.plot","plt.title","plt.xlabel","plt.ylabel"]))

## END HIDDEN TESTS

calcul de la vitesse limite¶

Sur le graphe précédent, on constate que la vitesse tend vers une valeur limite $u_{lim}$ que l'on va calculer.

  • Pour cela calculer dans un tableau dU, l'écart des vitesses sur chaque intervalle de temps: $U[k+1]-U[k]$

On obtient une estimation du temps limite $t_{lim}$ en déterminant où dans le tableau DU, la valeur devient, en valeur absolue, inférieure à une précision $\epsilon$ fixée: i.e. l'indice k tq. $|DU_k| < \epsilon$ avec une précision $\epsilon=10^{-2}$.

Pour calculer cet indice k, on utilise la fonction numpy where:

la fonction np.where(condition) (où condition est une condition sur un tableau X (p.e. np.where(X>0)) renvoie les indices du tableau où cette condition est réalisée

  • Mettre dans la variable k le numéro de cet intervalle. Calculer les valeurs moyennes de t et U sur cet intervalle k pour en déduire $t_{lim}$ et $u_{lim}$ que l'on mettra dans Tlim et Ulim.

  • Affichez la valeur de Tlim et Ulim.

In [5]:
dU = None
k = None
Tlim = None
Ulim = None
### BEGIN SOLUTION
dU = U[1:]-U[:-1]
k = np.where(np.abs(dU)<0.01)[0][0]
Tlim = (T[k]+T[k+1])/2.
Ulim = (U[k+1]+U[k])/2.
print(Tlim,Ulim,max(T),max(np.abs(U)))
### END SOLUTION

8.25 -79.6235 10.0 79.707
In [6]:
# ne rien écrire dans cette cellule
### BEGIN HIDDEN TESTS
assert(np.allclose(dU,U[1:]-U[:-1]))
assert(k == np.where(np.abs(dU)<0.01)[0][0])
assert(np.isclose(Tlim,(T[k]+T[k+1])/2.))
assert(np.isclose(Ulim,(U[k]+U[k+1])/2.))
### END HIDDEN TESTS

Calcul de la position¶

D'après l'étude précédente sur chaque intervalle de temps $[t_i,t_{i+1}]$, on connaît la vitesse $u_i$ en $t_i$ et la vitesse $u_{i+1}$ en $t_{i+1}$. On peut donc en déduire la distance parcourue $d_i$ par l'objet entre $t_i$ et $t_{i+1}$ en utilisant la vitesse moyenne avec la relation :

$$ d_i = \frac{u_{i}+u_{i+1}}{2} (t_{i+1} - t_{i}) $$

A partir de la position $x_i$ à $t_i$, on peut alors en déduire la position $x_{i+1}$ à $t_{i+1}$

$$ x_{i+1} = x_i + d_i $$

Connaissant la position $x_0$ à l'instant initiale, on a donc une méthode itérative pour calculer la position du solide au cours du temps.

  • Écrire une fonction python Trajectoire1D(t,u) , qui à partir du tableau des temps t (en sec.) et des vitesses u (en km/h), calcule la distance parcourue en m dans un tableau x et renvoie ce tableau.
  • Appliquer la fonction aux données précédente pour en déduire la trajectoire dans un tableau X
  • Tracer la trajectoire en fonction du temps pour vérifier

Faire attention aux unités

In [7]:
# definition de la fonction Trajectoire
## BEGIN SOLUTION
def Trajectoire1D(t,u):
    """ calcule et retourne la trajectoire 1D dans x en m 
    les arguments sont le tableau des temps t en s. et des vitesses u en km/h"""
    n = t.size
    x = np.zeros(n)
    for i in range(1,n):
        um = (u[i]+u[i-1])/2./3.6
        x[i] = x[i-1] + (t[i]-t[i-1])*um
    return x
# calcul 
X = Trajectoire1D(T,U)
print(min(X),max(X))
plt.plot(T,X)
plt.title("Trajectoire")
plt.xlabel("t sen sec.")
plt.ylabel("x en m")
## END SOLUTION

-186.81406944444447 0.0
Out[7]:
Text(0, 0.5, 'x en m')
In [8]:
# ne rien écrire dans cette cellule
### BEGIN HIDDEN TESTS
assert (X.size == T.size)
assert (np.allclose(X , Trajectoire1D(T,U)))
### END HIDDEN TESTS

Calcul du temps d'impact¶

On suppose que l'objet est lâché sans vitesse initiale à partir d'une hauteur H donnée ci-dessous.

On veut calculer le temps $t_{fin}$ où l'objet atteint la surface de la terre. On connaît la position à la fin des mesures, et on sait que l'objet en fin de trajectoire atteint une vitesse limite.

  • En déduire alors la valeur du temps d'impact $t_{fin}$ que l'on stockera dans la variable tfin

  • Pour vérifier, tracer la position de l'objet par rapport à la surface de la terre en fonction du temps. Tracer aussi la position du point d'impact dans la cellule suivante avec des titres, labels et légendes.

Vous devriez obtenir une figure similaire à la figure ci-dessous:

In [9]:
print("hauteur H = ",_H_)

hauteur H =  239.0
In [10]:
H = 0
tfin = 0
## BEGIN SOLUTION
H = _H_
# distance restante à parcourir
dr = H + X[-1]
tr = dr/np.abs(Ulim/3.6)
tfin = T[-1]+tr
print("Temps impacte: ", tfin,T[-1],dr,tr)
## END SOLUTION

Temps impacte:  12.359471136033958 10.0 52.18593055555553 2.3594711360339584
In [11]:
# tracer résultat et pt d'impact
## BEGIN SOLUTION
plt.plot(T,H+X,label='h(t)')
plt.plot([tfin],[0],'o',label='impact')
plt.title('Altitude h(t)')
plt.xlabel('temps [s]')
plt.ylabel('hauteur [en m]')
plt.legend()
plt.grid();
#plt.savefig('Hauteur_Temps.png')
## END SOLUTION

Commentaires sur la figure

%%% BEGIN SOLUTION

%%% END SOLUTION

In [12]:
# ne rien écrire dans cette cellule
### BEGIN HIDDEN TESTS
assert (tfin > T[-1])
### END HIDDEN TESTS
In [13]:
# ne rien écrire dans cette cellule
### BEGIN HIDDEN TESTS
assert(test_code(notebook_name,"cell-verif02",
                 ["plt.plot","plt.title","plt.xlabel","plt.ylabel"]))

### END HIDDEN TESTS

Etude Cinématique en 2D: cycliste sur un circuit elliptique¶

On veut étudier la cinématique d'un cycliste sur piste à partir de mesures de sa vitesse en fonction du temps.

Les données du circuit sont dans le fichier circuit.txt qui contient 3 colonnes correspondants au temps (en s.), et aux 2 composantes de la vitesse suivant x et y (en km/h).

lecture des données de vitesse¶

Dans la cellule suivante:

  1. lire les données à partir du fichier circuit.txt en utilisant la fonction loadtxt pour les mettre dans 3 tableaux T,U,V. .
  2. à partir des composantes U et V, calculer le tableau Unorm qui contient la norme de la vitesse en km/h à chaque instant.
  3. déterminer Umean la moyenne de Umorm , son minimum Umin et son maximum Umax. On pourra utiliser les fonctions (méthodes) numpy: mean(), max() et min()appliquées aux tableaux.
  4. Dans la cellule suivante, tracer U,V et Unorm en fonction du temps avec un titre, des labels et des légendes
  5. Commenter l'ordre de grandeur des vitesses calculées. C'est un cycliste chevronné !
  • Dans la suite, on utilisera comme unité : s. pour le temps, km/h pour la vitesse, et m pour la distance.
In [14]:
T = None
U = None
V = None
Unorm = None
Umean = None
Umin  = None
Umax  = None
### BEGIN SOLUTION
T,U,V = np.loadtxt('circuit.txt',unpack=True)
Unorm = np.sqrt(U**2 + V**2)
Umean = Unorm.mean()
Umax = Unorm.max()
Umin = Unorm.min()
print("Vitesse moyenne={:.1f}km/h min={:.1f} km/h max={:.1f}km/h".format(Umean,Umin,Umax))
### END SOLUTION

Vitesse moyenne=42.0km/h min=31.1 km/h max=51.8km/h
In [15]:
# tracer des vitesses
### BEGIN SOLUTION
plt.figure(figsize=(10,7))
plt.plot(T,U,label="U")
plt.plot(T,V,label="V")
plt.plot(T,Unorm,label="amplitude")
plt.xlabel('t [sec]')
plt.ylabel('U [km/h]')
plt.title("Vitesse du cycliste en km/h")
plt.legend();
#plt.savefig("vitesse.png")
### END SOLUTION
In [16]:
# ne rien écrire dans cette cellule
### BEGIN HIDDEN TESTS
# test sur les résultats (a executer)
assert((U.size == T.size) and (V.size == T.size))
assert(test_vitesse("circuit.txt",Umean,Umin,Umax))
### END HIDDEN TESTS
In [17]:
# ne rien écrire dans cette cellule
### BEGIN HIDDEN TESTS
assert(test_code("TP_trajectoire.ipynb","cell-verif00",["plt.plot","plt.title","plt.legend"]))

### END HIDDEN TESTS

Calcul de la trajectoire¶

A partir des composantes de vitesse U et V , on va calculer la trajectoire en fonction du temps, i.e. la position $X_i$ et $Y_i$ du point $P_i$ de la trajectoire pour chacun des temps $T_i$.

On va généraliser l'approche précédente, en calculant la distance parcourue $\mathbf{D}_{i}$ entre $T_{i-1}$ et $T_{i}$, mais de façon vectorielle:

$$ \mathbf{D}_{i} = \frac{1}{2}(\vec{U}_{i-1}+\vec{U}_{i-1}) dt \mbox { avec } dt = T_i - T_{i-1} $$

d'où la position $P_i$ du cycliste donnée dans le repère $(O,x,y)$ par le vecteur $\mathbf{OP}_{i}$

$$ \mathbf{OP}_{i} = \mathbf{OP}_{i-1} + \mathbf{D}_{i} \leadsto X_{i} = X_{i-1} + \frac{1}{2}(U_{i-1}+U_i) dt , Y_{i} = Y_{i-1} + \frac{1}{2}(V_{i-1}+V_i) dt $$

On constate que pour calculer la trajectoire X et Y on peut utiliser la fonction Trajectoire1D précédente.

Dans la cellule suivante:

  1. calculer la position des points de la trajectoire dans 2 vecteurs numpy X et Y en utilisant la fonction Trajectoire1D
  2. calculer la longueur de la trajectoire dans la variable L
  3. pour vérifier comparer la valeur de L avec la longueur moyenne Lmean que l'on calculera en utilisant la vitesse moyenne précédente Umean
  4. dans la cellule suivante, tracer la trajectoire (X,Y) ainsi que la vitesse sur le même graphe et sauvegarder la figure dans un fichier, que vous devez inclure ensuite dans votre compte rendu.

On utilisera la fonction plt.quiver(X,Y,U,V) qui permet de tracer des flèches de longueur U,V aux points X,Y. Pour éviter de tracer trop de flèches, on pourra prendre par exemple un point sur 10 en remplaçant X par X[::10] et idem pour tous les autres arguments.

In [18]:
X,Y = None,None
L = None
Lmean = None
### BEGIN SOLUTION
X = Trajectoire1D(T,U)
Y = Trajectoire1D(T,V)
# calcul longueur
L = 0
for i in range(1,X.size):
    L += np.sqrt((X[i]-X[i-1])**2 + (Y[i]-Y[i-1])**2)
Lmean = Umean*(T[-1]-T[0])/3.6
print("Longueur du circuit {:.2f} m".format(L))
print("Verification / longueur moyenne :",Lmean)
### END SOLUTION

Longueur du circuit 470.56 m
Verification / longueur moyenne : 469.9997864019623
In [19]:
# tracer
### BEGIN SOLUTION
plt.figure(figsize=(12,8))
plt.plot(X,Y)
plt.quiver(X[::10],Y[::10],U[::10],V[::10])
plt.title("Trajectoire")
plt.axis('equal');
#plt.savefig("trajectoire.png")
### END SOLUTION
In [20]:
# ne rien écrire dans cette cellule
### BEGIN HIDDEN TESTS
# test sur les résultats (a executer)
assert((X.size == T.size) and (Y.size == T.size))
print("Longueur du circuit {:.2f}m".format(L))
assert(test_trajectoire("circuit.txt",X,Y,L))
assert(np.abs(L-Lmean)<L*0.01)
### END HIDDEN TESTS

Longueur du circuit 470.56m
In [21]:
# ne rien écrire dans cette cellule
### BEGIN HIDDEN TESTS
assert(test_code("TP_trajectoire.ipynb","cell-verif1",["plt.plot","plt.quiver","plt.title"]))
### END HIDDEN TESTS

Calcul de l'accélération¶

Pour le calcul de l'accélération, on se place sur chaque intervalle k de temps $[T_{k},T_{k+1}]$. Sur cet intervalle, la vitesse varie de ${U_{k},V_{k}}$ de ${U_{k+1},V_{k+1}}$ et l'accélération mesure cette variation de vitesse suivant x dans $Gx_k$ et suivant y dans $Gy_k$

$$ Gx_k = \frac{\Delta U}{\Delta t} = \frac{U_{k+1} - U_{k}}{T_{k+1} - T_{k}}, Gy_k = \frac{\Delta V}{\Delta t} = \frac{V_{k+1} - V_{k}}{T_{k+1} - T_{k}}$$

Cette accélération est calculée pour un temps $Tm_k$ correspondant à la moyenne de $[T_{k}$ et $T_{k+1}]$

$$ Tm_k = \frac{T_{k+1}+T_{k}}{2} $$

Dans la cellule suivante:

  1. Calculer les composantes suivant x et y de l'accélération au milieu de chacun des segments $[T_{k},T_{k+1}]$ de la trajectoire, en utilisant l'accroissement de vitesse suivant x et suivant y. On mettra le résultat dans les vecteurs numpy Gx et Gy respectivement. Attention le calcul de l'accélération se fait en unité SI: $ms^{-2}$
    1. Calculer le temps au milieu du kème dans le vecteur numpy Tm

A partir du calcul précédent, on se propose de calculer l'accélération tangentielle et normale, i.e. le vecteur accélération projeté dans le référentiel local direct formé par la tangente et la normale en chaque point de la trajectoire. En se plaçant sur chaque intervalle k, on peut calculer la tangente $\vec{t}$ à partir de la vitesse moyenne $\vec{u}_m$

$$ \vec{t} = \frac{\vec{u}_m}{ |{\vec{u}_m}|}$$

puis en déduire la normale $\vec{n}$. Il suffit alors de projeter l'accélération calculée suivant x et y pour obtenir l'accélération tangentielle et normale.

  1. Dans une boucle sur chacun des intervalles k, calculer les composantes suivant x et y de la tangente dans tx, ty et de la normale dans nx, ny. En déduire la composante d'accélération tangentielle Gt_k et normale Gn_k. On mettra ces composantes dans les tableaux numpy Gn et Gt.

  2. Tracer sur une figure l'accélération suivant x, suivant y, l'accélération normale et l'accélération tangentielle avec comme unité le g ($g=9.81 ms^{-2}$)

  3. Commenter le résultat dans la cellule suivante, et en particulier la variation de l'accélération tangentielle.

Faire attention aux unités

In [22]:
g = 9.81
Gx = None
Gy = None
Gt = None
Gn = None
### BEGIN SOLUTION
N = X.size
# calcul acceleration 
Gx = np.zeros(N-1)
Gy = np.zeros(N-1)
for i in range(0,N-1):
    Gx[i] = (U[i+1]-U[i])/(T[i+1]-T[i])/3.6
    Gy[i] = (V[i+1]-V[i])/(T[i+1]-T[i])/3.6
# temps
Tm = (T[:-1]+T[1:])/2.
# projection
Gt = np.zeros(N-1)
Gn = np.zeros(N-1)
for i in range(0,N-1):
    um = (U[i]+U[i+1])/2.
    vm = (V[i]+V[i+1])/2.
    tx = um/np.sqrt(um**2+vm**2)
    ty = vm/np.sqrt(um**2+vm**2)
    Gt[i] = Gx[i]*tx + Gy[i]*ty
    Gn[i] =-Gx[i]*ty + Gy[i]*tx
# tracer
plt.figure(figsize=(14,8))
plt.plot(Tm,Gx/g,label="Gx")
plt.plot(Tm,Gy/g,label='Gy')
plt.plot(Tm,Gn/g,label="Gn")
plt.plot(Tm,Gt/g,label='Gt')
plt.legend()
plt.xlabel("t en sec.")
plt.ylabel("$\gamma$ en g")
plt.title("Composantes de l'accélération");
#plt.savefig("acceleration.png")
### END SOLUTION

Commentaires sur la figure

%%% BEGIN SOLUTION

%%% END SOLUTION

In [23]:
# ne rien écrire dans cette cellule
### BEGIN HIDDEN TESTS
# test sur les résultats (a executer)
assert((Gx.size == T.size-1) and (Gy.size == T.size-1))
assert((Gn.size == T.size-1) and (Gt.size == T.size-1))
assert(test_acceleration("circuit.txt",Gt/g,Gn/g,None))
assert(test_code("TP_trajectoire.ipynb","cell-verif2",["plt.plot","plt.title"]))
### END HIDDEN TESTSS

Calcul de l'angle d'inclinaison (optionnel)¶

L'objectif est de calculer l'angle $\theta$ d'inclinaison du cycliste, qui correspond à l'équilibre des forces dans le référentiel lié au cycliste dans le plan de la figure:

$$ \vec{F_c} + \vec{P} + \vec{R} = \vec{0}$$

La force centrifuge $\vec{F_c} = - M \vec{\gamma_n}$ est une force fictive qui est opposée à l'accélération dans la direction normale, et qui traduit cet équilibre statique dans le référentiel lié au cycliste. Dans le référentiel fixe, on retrouve le PFD dans le plan de la figure:

$$ \vec{F_c} + \vec{P} = M \vec{\gamma_n}$$

Remarque dans la direction tangentielle (i.e. la direction perpendiculaire au plan de la figure) on a un équilibre entre l'accélération tangentielle, la force de traction, et les forces de frottements.

Pour calculer $\theta$, il faut donc utiliser l'accélération tangentielle et normale, i.e. le vecteur accélération projeté dans le référentiel local directe formé par la tangente et la normale en chaque point de la trajectoire.

  1. Calculer 'angle d'inclinaison $\theta$ en degré du cycliste pour composer la force centrifuge en prenant $g=9.81 \ m/s^{-2}$. On mettra le résultat dans le tableau numpy THETA. On pourra utiliser la fonction numpy arctan2(a,b) qui calcule l'angle en radian dont la tangente vaut $a/b$.

  2. Tracer l'angle d'inclinaison du cycliste.

  3. Commenter l'ordre de grandeur de l'angle d'inclinaison $\theta$, et discuter sa variation.

In [24]:
# calcul de l'angle
### BEGIN SOLUTION
THETA = np.rad2deg(np.arctan2(Gn,g))
# tracer
plt.figure(figsize=(14,8))
plt.plot(Tm,THETA)
plt.ylim(-20,20)
plt.xlabel("t en sec.")
plt.ylabel("$\\theta$ en degré")
plt.title("inclinaison en degré");
#plt.savefig("acceleration.png")
### END SOLUTION

Commentaires sur la figure

%%% BEGIN SOLUTION

%%% END SOLUTION

In [25]:
# ne rien écrire dans cette cellule
### BEGIN HIDDEN TESTS
# test sur les résultats (a executer)
assert(THETA.size == T.size-1)
assert(test_acceleration("circuit.txt",Gt/g,Gn/g,THETA))
assert(test_code("TP_trajectoire.ipynb","cell-verif3",["plt.plot","plt.title"]))
### END HIDDEN TESTSS

FIN¶