program cell_omp

  implicit none

  integer :: Ndim,Niter
  integer :: x,n,maxval=0
  real*8 :: sum=0
  integer :: myleft,myright
  character(len=80) :: paramin
  integer,pointer :: buffer(:)
  integer,pointer :: nextbuffer(:)
  integer,pointer :: tmp(:)

  ! Exit if # of arguments < 2
  if (command_argument_count().lt.2) then
     write(*,*) "Usage: cell_omp <Ndim> <Niter>"
     stop
  endif

  call getarg(1,paramin)
  read(paramin,*) Ndim
  call getarg(2,paramin)
  read(paramin,*) Niter

  write(*,*)
  write(*,'(A,T24,I8)') "Size of the array: ",Ndim
  write(*,'(A,T24,I8)') "Number of iterations: ",Niter

  allocate(buffer(Ndim))
  allocate(nextbuffer(Ndim))
  allocate(tmp(Ndim))

  ! Initialize the array.
!$omp parallel do
  do x=1,Ndim
     buffer(x)=mod(x-1,10)
  enddo

  ! Main loop
!$omp parallel do
  do n=1,Niter
     do x=1,Ndim
        myleft=x-1
        if (myleft.lt.1) myleft=myleft+Ndim
        myright=x+1
        if (myright.gt.Ndim) myright=myright-Ndim
        nextbuffer(x)=mod(buffer(myleft)+buffer(myright),10)
     enddo
     ! Swap the pointers
     tmp => buffer
     buffer => nextbuffer
     nextbuffer => tmp
  enddo

  ! Find the maximum and average.
!$omp parallel do reduction(+:sum) reduction(max:maxval)
  do x=1,Ndim
     sum=sum+buffer(x)
     if (buffer(x).gt.maxval) maxval=buffer(x)
  enddo

  write(*,*)
  write(*,'(A,T24,I4)') "Maximum: ",maxval
  write(*,'(A,T24,F10.6)') "Average: ",sum/float(Ndim)
  write(*,*)

  stop
end program cell_omp