C 		Programme de simulation d'un pendule simple par la methode RK4
C		Dominique LEFEBVRE - mars 2006

	  PROGRAM PenduleRK4

	  IMPLICIT None
C 		Declaration des routines externes
C		ResolutionRK4 : routine d'implementation du schema RK4
C		Derivee       : description du systeme differentiel de l'oscillateur	  
      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 g, PI
      DATA g/9.81/, PI/ 3.141592654/
	 

C		Declaration des variables
	  INTEGER i, NbPas
	  REAL y(NbEquation), t, l, PasTemps, a0

C		Declaration des variables communes
	  COMMON /Pendule/l,g

C       Saisie des parametres du mouvement
      WRITE(*,'(a,$)') 'Longueur du pendule (en m) : '
      READ(*,*) l
      WRITE(*,'(a,$)') 'angle initial du pendule (en degres) : '
      READ(*,*) a0
      WRITE(*,'(a,$)') 'Pas temporel de calcul (en s) : '
      READ(*,*) PasTemps	  
      WRITE(*,'(a,$)') 'Nombre de pas : '
      READ(*,*) NbPas
      
C		Declaration des conditions initiales pour l'integration
	  y(1) = a0*PI/180.		! angle initial du pendule
	  y(2) = 0.0			! vitesse angulaire initiale
		
C		Ouverture du fichier de stockage des resultats de calcul
C		Ce fichier sera trace par GnuPlot
	  OPEN(10, FILE='PenduleRK4.dat')

C		Ecriture des conditions initiales
	  t = 0.0
	  WRITE(10,* ) t, y(1)

C		Boucle de calcul
	  DO i=1, NbPas
		t = i*PasTemps
		CALL ResolutionRK4(t, y, PasTemps, NbEquation, Derivee)
		WRITE(10,*) t, y(1)*180./PI		! notez la conversion en ° 	
	  ENDDO

C		Fermeture du fichier de donnees
	  CLOSE(10)

	  STOP
	  END