Mouvement plan d'un solide

    Scilab est également bien adapté pour la résolution des équations de mouvement d'un solide ou d'un ensemble de solides. Ces équations sont des systèmes d'équations différentielles, linéaires ou non linéaires, du deuxième ordre par rapport au temps.

    Comme exemple nous avons pris le mouvement d'un solide sur un plan incliné (figure ci-dessous). Ce problème est traité au chapitre 23 (Section 23.2) de l'ouvrage Mécanique des Solides Rigides (2006). Ce problème est également traité dans l'ouvrage Mécanique des Solides (2012), Partie V

 
 

    Les équations de mouvement peuvent être écrites (27.30) sous la forme:                            

 

 

x et y sont les coordonnées du centre de masse G, ψ est l'angle de rotation du solide (S), ft et fr les coefficients de frottement en translation et rotation.

   Nous donnons ci-dessous un exemple du fichier d'éxécution de Scilab et les résultats obtenus pour les trajectoires des mouvements:

    -- du centre de masse (courbe en rouge);

    -- d'un point lié au solide (courbe en bleu);

pour diverses valeurs de l'inclinaison du plan et des coefficients de frottement.

    Les équations différentielles sont résolues en utilisant une méthode de Runge-Kutta d'ordre 4.

    Pour une analyse plus approfondie, le lecteur pourra se rapporter au chapître 27 de l'ouvrage Mécanique des Solides (2012), Partie VI.

 

Fichier.sce

//Mouvement plan sur plan d'un solide

ta=5 // temps final
alpha=20;//inclinaison du plan
alpha=alpha*%pi/180;
ft=0.5
fr=0.8
x0=0;y0=0;psi0=0//coordonnées et angle initiaux
vx0=0;vy0=5;vpsi0=360;//vitesses initiales
g=9.807;
n=200;//nombre de pas de calcul
t0=0;
pas=(ta-t0)/n;
t(1)=t0;x(1)=x0;y(1)=y0;psi(1)=psi0*%pi/180;vx(1)=vx0;vy(1)=vy0;vpsi(1)=vpsi0*%pi/180;
for i=1:n
    ti=t(i);xi=x(i);yi=y(i);psii=psi(i);vxi=vx(i);vyi=vy(i);vpsii=vpsi(i);
    Fvxi=g*sin(alpha)-ft*vxi;
    Fvyi=-ft*vyi;
    Fvpsii=-fr*vpsii;
    phi1x=pas*vxi;phi1y=pas*vyi;phi1psi=pas*vpsii;
    phi1vx=pas*Fvxi;phi1vy=pas*Fvyi;phi1vpsi=pas*Fvpsii;
    ti=ti+pas/2;
    xi=x(i)+phi1x/2;yi=y(i)+phi1y/2;psii=psi(i)+phi1psi/2;
    vxi=vx(i)+phi1vx/2;vyi=vy(i)+phi1vy/2;vpsii=vpsi(i)+phi1vpsi/2;
    Fvxi=g*sin(alpha)-ft*vxi;
    Fvyi=-ft*vyi;
    Fvpsii=-fr*vpsii;
    phi2x=pas*vxi;phi2y=pas*vyi;phi2psi=pas*vpsii;
    phi2vx=pas*Fvxi;phi2vy=pas*Fvyi;phi2vpsi=pas*Fvpsii;
    xi=x(i)+phi2x/2;yi=y(i)+phi2y/2;psii=psi(i)+phi2psi/2;
    vxi=vx(i)+phi2vx/2;vyi=vy(i)+phi2vy/2;vpsii=vpsi(i)+phi2vpsi/2;
    Fvxi=g*sin(alpha)-ft*vxi;
    Fvyi=-ft*vyi;
    Fvpsii=-fr*vpsii;
    phi3x=pas*vxi;phi3y=pas*vyi;phi3psi=pas*vpsii;
    phi3vx=pas*Fvxi;phi3vy=pas*Fvyi;phi3vpsi=pas*Fvpsii;
    ti=ti+pas/2;
    xi=x(i)+phi3x;yi=y(i)+phi3y;psii=psi(i)+phi3psi;
    vxi=vx(i)+phi3vx;vyi=vy(i)+phi3vy;vpsii=vpsi(i)+phi3vpsi;
    Fvxi=g*sin(alpha)-ft*vxi;
    Fvyi=-ft*vyi;
    Fvpsii=-fr*vpsii;
    phi4x=pas*vxi;phi4y=pas*vyi;phi4psi=pas*vpsii;
    phi4vx=pas*Fvxi;phi4vy=pas*Fvyi;phi4vpsi=pas*Fvpsii;

    t(i+1)=ti;
    x(i+1)=x(i)+(phi1x+2*phi2x+2*phi3x+phi4x)/6;
    y(i+1)=y(i)+(phi1y+2*phi2y+2*phi3y+phi4y)/6;
    psi(i+1)=psi(i)+(phi1psi+2*phi2psi+2*phi3psi+phi4psi)/6;
    vx(i+1)=vx(i)+(phi1vx+2*phi2vx+2*phi3vx+phi4vx)/6;
    vy(i+1)=vy(i)+(phi1vy+2*phi2vy+2*phi3vy+phi4vy)/6;
    vpsi(i+1)=vpsi(i)+(phi1vpsi+2*phi2vpsi+2*phi3vpsi+phi4vpsi)/6;
end
//mouvement du point M(xs,ys)
xs=4;ys=0;
for i=1:n
    ci=cos(psi(i));si=sin(psi(i));
    xm(i)=x(i)+xs*ci-ys*si;
    ym(i)=y(i)+xs*si+ys*ci;
end
clf

plot(xm,ym,"b",x,y,"r")
xtitle("Inclinaison = 20°  fr = 0.5 /s  ft = 0.8 /s","m","m");

 

Trajectoires obtenues