The program









To solve numerically the differential systems in the pendulum equation, we just need to discretize it. I had first chosen an easy Euler (explicit) scheme for each equations. It appears that this scheme is unstable for a damping factor equal to zero. With the recommendation of A.Merle, I changed the scheme, and the second equation is solved with an implicit scheme. The results seem to be very good.

The coding resulting is simple and it is like coding a serie.

This program is written in matlab language and is suitable for each case developed in this project. In the beginning, the variable type allow us to choose the type of pendulum we want to simulate. Then, several parameters can be modified, such as the damping factor and the excitation.



clear;

type=input('Type de pendule a modeliser (1=normal, 2=avec forcage): ');

nbcas=input('Nombre de cas a tester: ');

n=input('Nb de pas de tps simules: ');

Dt=0.01;
w=9.81;
 

for j=1:nbcas
 
    u(1,j)=input('Angle de depart: ');
    v(1,j)=input('Vitesse de depart: ');
 
    t(1)=0;
 
    if type == 1
 
      l=input('Valeur de l ammortissement: ');

      for i=1:n
         u(i+1,j)=v(i,j)*Dt + u(i,j);
         v(i+1,j)= 1/(1+Dt*l) * (-Dt*w*sin(u(i+1,j)) + v(i,j));
         t(i+1)=t(i)+Dt;
      end
 
    end
 
    if type == 2
 
      m=input('Valeur du forcage: ');
      k=input('Valeur k de l ammortissement (de forme k*teta^2): ');
 
      for i=1:n
         u(i+1,j)=v(i,j)*Dt + u(i,j);
         v(i+1,j)= 1/(1+Dt*(k*u(i+1,j)*u(i+1,j)-m)) * (-Dt*w*sin(u(i+1,j)) + v(i,j));

         t(i+1)=t(i)+Dt;
      end
 
    end
 
 
end
 

plot(u,v);
figure;
plot(t,u);


The outputs of the program are the phase portrait plot(u,v) and the evolution of the angle u during the time plot(t,u).