!**********************************************************************!
!
! File: cromosc.f90 (06-Apr-2001)
!                   (27-May-2005)
!
! Harmonischer Oszillator mit Cromer-Verfahren
!
!   dv/dt = -x
!   dx/dt = v
!
! dt..........Zeitschritt
! ntpr........Ausgabe alle ntpr Zeitschritte
! ntmax.......maximale Anzahl Zeitschritte
! x(0),v(0)...Startwerte
!
! Ausgabe auf File "cromosc.dat"
!
!**********************************************************************!

program cromosc

!**********************************************************************!

  implicit none

  integer,parameter::iout=2

  integer::nt,ntpr,ntmax
  real::dt,t,v,vn,x,xn
  
  write(unit=*,fmt="("" dt,ntpr,ntmax="")",advance="no")
  read(unit=*,fmt=*) dt,ntpr,ntmax
  write(unit=*,fmt="("" x(0),v(0)="")",advance="no")
  read(unit=*,fmt=*) x,v
  
  open(unit=iout,file="cromosc.dat",status="replace",action="write", &
       form="formatted")
  
  t=0.0
  write(unit=iout,fmt="(tr1,f10.5,2(tr1,es12.5))") t,x,v
  
  do nt=1,ntmax
    vn=v-dt*x      ! Cromer
    xn=x+dt*vn     !
    t=nt*dt
    x=xn
    v=vn
    if(ntpr*(nt/ntpr)==nt) then
      write(unit=iout,fmt="(tr1,f10.5,2(tr1,es12.5))") t,x,v
    end if
  end do
  
  close(unit=iout)
  
end program cromosc

!**********************************************************************!