Bonjour à tous j'ai récemment codé une résolution numérique d'équations différentielles du second degré avec la méthode de runge kutta 4 pour obtenir la trajectoire d'une particule.J'ai réussi à compiler le programme snas aucune erreur, cependant les résultats obtenus sont faux, j'en conclut que j'ai surement du faire une erreur dans mon programme.J'aurais voulu savoir si quelqu'un pourrait m'aider
à résoudre ce problème. Voici mon programme pour la trajectoire de la particule en fonction du temps.
Merci à tous d'avance.
PROGRAM ParticuleRK4
IMPLICIT None
c Declaration des routines externes
c ResolutionRK4:routine d'implementation du schema RK4
c Derivee :description du sysrtème differentiel du mouvement
EXTERNAL ResolutionRK4
EXTERNAL Derivee
c Declaration des parametres RK4
c Le systeme differentiel comporte deux equations
INTEGER NbEquation
PARAMETER (NbEquation=2)
c Declaration des constantes
REAL R,k,l,alfa
DATA R/20E-6/,k/440667./,l/-9.79487/,alfa/1.55235/
c Declaration des variables
INTEGER i, NbPAS
REAL z(NbEquation),t,PasTemps,z0,w0
c Saisie des paramètres du mouvement
WRITE(*,'(a,$)') 'Position initiale de la particule (en m) :'
READ(*,*) z0
WRITE(*,'(a,$)') 'Vitesse initiale de la particule (en m/s) :'
READ(*,*) w0
WRITE(*,'(a,$)') 'Pas temporel de calcul (en s) :'
READ(*,*) PasTemps
WRITE(*,'(a,$)') 'Nombre de pas :'
READ(*,*) NbPas
c Declaration des conditions initiales pour l'intégration
z(1) = z0 !position initiale de la paricule
z(2) = w0 !vitesse initiale de la particule
c Ouverture du fichier de stockage des resultats de calcul
c Ce fichier sera trace par Gnuplot
OPEN(10, FILE='resolutions.dat')
c Ecriture des conditios initiales
t = 0.0
c Boucle de calcul
DO i=1, NbPas
WRITE(10,*) t, z(i)
t = i*PasTemps
CALL ResolutionRK4(t,z,PasTemps,NbEquation)
ENDDO
WRITE(10,*) t, z(NbPas)
c Fermetrue du fichier de donnees
CLOSE(10)
STOP
END
c------------------------------------------------------------------------------------
c Routine de definition equations differentielles de
c la particule
c Par habitude et pour faciliter la reutilisation,j'appelle
c cette fonction Derivee.
SUBROUTINE Derivee(t,z,dz, NbEquation)
IMPLICIT None
REAL k,alfa
DATA k/440667./,alfa/1.55235/
INTEGER NbEquation
REAL t, z(NbEquation), dz(NbEquation)
c Description du systeme diferentiel de la particule
dz(1) = z(2)
dz(2) = -k*z(1)+alfa
RETURN
END
c Subroutine de resolution d'une equation differentielle
c par la méthode de Runge-Kutta d'ordre 4
SUBROUTINE ResolutionRK4(t,z,dt,NbEquation)
c t : abscisse
c z() : valeurs retournees par le systeme differentiel
c dt : pas de calcul
c NbEquation : nombre d'eqution dans le systeme differentiel
c Derivee : nom de la fonction decrivant le systeme differentiel
IMPLICIT NONE
c EXTERNAL Derivee !fonction de calcul de la derivee
c
c Typage de l'interface de ResolutionRK4
c
INTEGER NbEquation
REAL t,z(NbEquation), dt
c
c Declaration des variables locales
REAL pred1(NbEquation),pred2(NbEquation),pred3(NbEquation)
REAL pred4(NbEquation), ztemp(NbEquation), halfdt
INTEGER i
halfdt = dt/2
c Premiere prediction
CALL Derivee(t,z(NbEquation),pred1(NbEquation),NbEquation)
DO i=1, NbEquation
ztemp(i) = z(i) + pred1(i)*halfdt
ENDDO
c Seconde prediction
CALL Derivee(t+halfdt,ztemp(NbEquation),pred2(NbEquation)
& , NbEquation)
DO i=1, NbEquation
ztemp(i) = z(i) + pred2(i)*halfdt
ENDDO
c Troisieme prediction
CALL Derivee(t+halfdt,ztemp(NbEquation),pred3(NbEquation)
& , NbEquation)
DO i=1, NbEquation
ztemp(i) = z(i) + pred3(i)*dt
ENDDO
c Quatrieme prediction
CALL Derivee(t+halfdt,ztemp(NbEquation),pred4(NbEquation)
& , NbEquation)
c Estimation de y pour le pas suivant
DO i=1, NbEquation
z(i)=z(i)+(dt/6)*(pred1(i)+2*pred2(i)+2*pred3(i)+pred4(i))
ENDDO
END
c ----------------------------------------------------------------
__________________________
studientj
