
program md
    use inputmodule
	use initialconfignewmodule
	use accmodule
	!use averagepropmodule
    !use velrescalemodule
	use velverletmodule

	Implicit none
	integer :: Elec,Proton,N,Ne,t,itr,seed,i,j,p,rescale,total
	double precision dt,length,Temp,c,a,b
	double precision, allocatable, dimension(:) :: rx,ry,rz,vx,vy,vz,ax,ay,az
	double precision, allocatable, dimension(:) :: Potential,Etotal,caltemp,Kinetic,LJ   !average temperature at each instant time


    call input(Elec,Proton,N,Ne,length,dt,seed,Temp,itr,a,b)



	allocate (rx(N),ry(N),rz(N),vx(N),vy(N),vz(N),ax(N),ay(N),az(N))

	allocate (Caltemp(itr),Kinetic(itr),LJ(itr),Potential(itr),Etotal(itr))



	call initialconfignew(N,length,seed,rx,ry,rz,vx,vy,vz)  !that will give initial value of x,v at time t = 1


    open(unit=6,file="newmd.xyz")


    write(6,*) "330"
    write(6,*) "particle dynamics p & e"

   ! do j=1,N
   !    write(6,*) "Ar",rx(j),ry(j),rz(j)
   ! end do


    call acc(Elec,Proton,N,Ne,1,length,rx,ry,rz,ax,ay,az,LJ,a,b)  !to call acc. at initial value t=1

    i=2

    do i=2,itr

            call velverlet(Elec,Proton,N,Ne,dt,i,length,rx,ry,rz,vx,vy,vz,ax,ay,az,LJ,a,b)   !will give new values for x,v at t

            
	if (i>1995000) then
	    do j=1,N
                write(6,*) "Ar",rx(j),ry(j),rz(j)
            end do
	end if 

	end do

    close(unit=6)



end program md


		
