Bonjour,
j'applique la méthode de runge kutta à l'ordre 2 pour résoudre l'équation différentielle du pendule simple sans frottement. Je la réduit à un système d'ordre 1 comme ceci:
u=d(theta)/dt
du/dt=(g/l)cos(theta)=F(theta,g,l)
avec F la fonction suivante:
Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10 long double F(long double theta, long double g, long double a) { long double r; r=(g/a)*cosl(theta); return r; }
je calcule les coefficients qui sont dans un tableau k[4][4] en prenant comme vecteur : (dt/dt=1,d^2 t/dt^2=0,dtheta_dt,F(theta,g,l)=d^2theta/dt^2)
et puis la formule de récurrence :
Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25 k[0][0]=h; k[0][1]=0; k[0][2]=h*dtheta_dt[i]; k[0][3]=h*F(theta[i],g,l); k[1][0]=h; k[1][1]=h*k[0][1]/2; k[1][2]=h*(dtheta_dt[i]+(k[0][2]/2)); k[1][3]=h*F((theta[i]+(k[0][3]/2)),g,l); k[2][0]=h; k[2][1]=h*k[1][1]/2; k[2][2]=h*(dtheta_dt[i]+(k[1][2]/2)); k[2][3]=h*F((theta[i]+(k[1][3]/2)),g,l); k[3][0]=h; k[3][1]=h*k[2][1]; k[3][2]=h*(dtheta_dt[i]+k[2][2]); k[3][3]=h*F((theta[i]+k[2][3]),g,l);
Les résultats sont représentés sur les deux graphiques joints. pend.pdf montre le résultat sur 300s (abcisse) et pend2.pdf sur 30s.
Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3 t[i+1]=t[i]+h*((k[0][0]+2*(k[1][0]+k[2][0])+k[3][0])/6); theta[i+1]=theta[i]+h*((k[0][2]+2*(k[1][2]+k[2][2])+k[3][2])/6);
Je ne comprends pas pourquoi l'amplitude des oscillations augmente avec le
temps comme c'est montré sur pend.pdf. Est ce lié à l'erreur de la méthode numérique utilisée ? je prends comme conditions initiales t[0]=0,theta[0]=0 et dtheta_dt[0]=0 et h le pas d'incrémentation égal à 0.01. je pense que le calcul des coef est bon.
Je vous remercie d'avance pour toute aide.
Partager