Le pendule pesant est en liaison pivot avec un chariot, lequel glisse en translation rectiligne sur un bâti. On suppose que le bâti constitue un référentiel galiléen. Le moteur lié au bâti entraîne en rotation une roue. Le chariot est lié à cette roue par une bielle. On a ainsi un système de type manivelle-bielle qui permet de convertir le mouvement de rotation du moteur en mouvement de translation du chariot.
Figure pleine pageOn note L la distance entre le centre de gravité du pendule et son axe de rotation, θ son angle avec la direction verticale. Soit un point O sur l'axe de rotation du moteur. On note a=OA la longueur de la manivelle et b=AB la longueur de la bielle. L'extrémité B de celle-ci (liaison pivot avec le chariot) a pour coordonnées (x,h).
L'angle que fait la manivelle avec l'axe horizontal Ox est α=ωt, où ω
est la vitesse angulaire de rotation du moteur. L'angle que fait la bielle avec cet axe est β. En décomposant
les vecteurs
Remarque : sur la figure ci-dessus, h est négatif.
En éliminant l'angle β de ces équations on obtient :
Voici deux exemples de courbe x(t) :
import numpy from matplotlib.pyplot import * def position(a,b,h,t): alpha = 2*numpy.pi*t return a*numpy.cos(alpha)+numpy.sqrt(b**2-numpy.power(h-a*numpy.sin(alpha),2)) figure() N = 500 t = numpy.arange(N)*2.0/N a=10 b=40 plot(t,position(a,b,0,t),label="h=0") plot(t,position(a,b,-10,t),label="h=-10") xlabel("t") ylabel("x") legend(loc="lower right") grid()Figure pleine page
Lorsque h=0 et
On utilisera cette expression pour la modélisation du pendule.
Le référentiel lié au chariot est non galiléen. L'accélération d'entraînement dans ce référentiel est :
Les forces de pesanteur s'appliquent au centre de gravité G du pendule. Il en est de même des forces d'inertie d'entraînement (pour un référentiel en translation). Soit m la masse du pendule et J son moment d'inertie par rapport à son axe de rotation. Le théorème du moment cinétique appliqué, dans le référentiel du chariot, par rapport à l'axe de rotation, conduit à l'équation différentielle suivante :
On a introduit un frottement proportionnel à la vitesse angulaire, avec un coefficient ν.
Il faut ajouter à cette équation deux conditions initiales :
L'équation différentielle est non linéaire. Il est toutefois possible d'obtenir une équation linéaire en supposant que l'angle reste petit :
La solution analytique exacte de cette équation est connue. En revanche, l'équation complète non linéaire ne peut être résolue que par des calculs numériques, par exemple avec la méthode d'Euler.
Pour simplifier l'étude on divise l'équation par J et on pose :
Il s'agit de la pulsation propre du pendule linéaire (petits angles), c'est-à-dire sa pulsation d'oscillation lorsque le chariot est immobile (ω=0). Soit T la période correspondante. Il est intéressant de faire le changement de variable suivant pour le temps :
Le temps t' est sans dimensions. Pour faire le changement de variable, on utilise les identités suivantes :
On obtient ainsi l'équation différentielle
Cette équation est adimensionnée, car les variables θ et t' sont sans dimensions. Lorsqu'on résout numériquement une équation du mouvement, on a toujours intérêt à la mettre sous forme adimensionnée. La période des petites oscillations libres est T'=1. On voit apparaître trois paramètres sans dimensions :
Voici finalement l'équation différentielle à étudier :
Pour résoudre numériquement cette équation avec les deux conditions initiales, il faut la mettre sous la forme d'un système différentiel du premier ordre, en introduisant la vitesse angule du pendule, que l'on notera va :
Le système est défini par une fonction (on étudie tout d'abord un système sans frottement) :
import math import scipy.integrate gamma = 0.0 A = 0.3 f = 1.0 w02 = (2*math.pi)**2 def systeme(Y,t): return [Y[1],-w02*math.sin(Y[0])-gamma*Y[1]-A*w02*f**2*math.cos(Y[0])*math.cos(2*math.pi*f*t)]
On fixe une condition initiale :
Yi = [math.pi*0.5,0.0]
On calcule les instants pour lesquels on veut les valeurs approchées :
T = 100.0 t = numpy.arange(start=0,stop=T,step=0.01)
On fixe la tolérance absolue et la tolérance relative :
atol=1e-8 rtol=1e-8
On fait l'intégration numérique :
tab_y = scipy.integrate.odeint(systeme,Yi,t,rtol=rtol,atol=atol) yt = tab_y.transpose() theta = yt[0] va = yt[1] figure() plot(t,theta) xlabel("t") ylabel("theta") axis([0,20,-2,2]) grid()figB.pdf
Pour analyser le mouvement, on peut faire une analyse spectrale :
import numpy.fft tfd = numpy.fft.fft(theta) N = theta.size spectre = numpy.absolute(tfd)/N freq = numpy.arange(N)*1.0/T figure() plot(freq,spectre) xlabel("freq") ylabel("theta") axis([0,5,0,spectre.max()])figC.pdf
On peut aussi représenter la trajectoire dans l'espace des phases. La dimension de l'espace des phases est égale au nombre de degrés de libertés du système (au sens de système du premier ordre). Ici, il y a trois degrés de liberté (θ,va,α).
On peut écrire le système en faisant apparaître explicitement le degré de liberté correspondant à la rotation de la manivelle :
La représentation graphique (à deux dimensions) de la trajectoire dans l'espace des phases à 3 dimensions n'est pas faisable. Une première solution est de représenter la projection sur le plan (θ,va) :
figure() plot(theta,va) xlabel("theta") ylabel("va") grid()figD.pdf
Une représentation simplifiée est possible en remarquant que la variable α est de période 2π. En trace alors dans le plan de projection seulement les points pour lesquels α est multiple de 2π. Cette représentation est appelée section de Poincaré.
t = numpy.arange(start=0,stop=T,step=1.0/f) tab_y = scipy.integrate.odeint(systeme,Yi,t,rtol=rtol,atol=atol) yt = tab_y.transpose() theta = yt[0] va = yt[1] figure() plot(theta,va,".") xlabel('theta') ylabel('va') grid()figE.pdf
Les points se répartissent sur une courbe fermée, et remplissent cette courbe lorsque t tend vers l'infini. Il s'agit d'un mouvement quasipériodique.
Voyons ce qu'il se passe si l'on augmente l'amplitude A. Pour cela, il faut augmenter la longueur de la manivelle.
A = 2.0 T = 100.0 t = numpy.arange(start=0,stop=T,step=0.01) Yi = [math.pi*0.5,0.0] tab_y = scipy.integrate.odeint(systeme,Yi,t,rtol=rtol,atol=atol) yt = tab_y.transpose() theta = yt[0] va = yt[1] figure() plot(t,theta) xlabel("t") ylabel("theta") grid()figF.pdf
On obtient un mouvement chaotique, au cours duquel le pendule fait des tours complets, dans un sens puis dans l'autre. Une caractéristique importante d'un mouvement chaotique est la très grande sensibilité aux conditions initiales. Refaisons le calcul avec un angle initial légèrement différent :
A = 2.0 T = 100.0 t = numpy.arange(start=0,stop=T,step=0.01) Yi = [math.pi*0.5+1e-3,0.0] tab_y = scipy.integrate.odeint(systeme,Yi,t,rtol=rtol,atol=atol) yt = tab_y.transpose() theta = yt[0] va = yt[1] plot(t,theta) xlabel("t") ylabel("theta") grid()figG.pdf
Au début, les deux mouvements sont très proches, mais l'écart grandit très vite. On obient finalement deux trajectoires complètement différentes. En ce sens, l'évolution d'un système chaotique n'est pas prévisible bien que son équation différentielle soit déterministe. Il faudrait une précision infinie sur la condition initiale pour prévoir l'état du système à un instant t quelconque, ce qui est impossible en pratique.
Voyons l'analyse spectrale de ce mouvement chaotique :
tfd = numpy.fft.fft(theta) N = theta.size spectre = numpy.absolute(tfd)/N freq = numpy.arange(N)*1.0/T figure() plot(freq,spectre) xlabel("freq") ylabel("theta") axis([0,1,0,spectre.max()])figH.pdf
Un mouvement chaotique a un spectre continu. Voyons la section de Poincaré :
t = numpy.arange(start=0,stop=T,step=1.0/f) tab_y = scipy.integrate.odeint(systeme,Yi,t,rtol=rtol,atol=atol) yt = tab_y.transpose() theta = yt[0]%(2*math.pi)-math.pi va = yt[1] figure() plot(theta,va,".") xlabel('theta') ylabel('va') grid()figI.pdf
Pour un mouvement chaotique, les points de la section de Poincaré se répartissent sur une surface, et non pas sur une courbe.
f=1.0
A=0.3
gamma=0.0