--- trunk/mdtools/md_code/simulation_module.F90	2002/12/04 21:19:38	194
+++ trunk/mdtools/md_code/simulation_module.F90	2003/01/02 21:45:45	222
@@ -6,10 +6,12 @@ module simulation
   PRIVATE
 
 
-  
-  real ( kind = dp ), public, dimension(3) :: box
+  integer :: nLR
 
 
+  real ( kind = dp ), dimension(3) :: box
+
+
   real ( kind = dp ), public :: rlist
   real ( kind = dp ), public :: rcut
   real ( kind = dp ), public :: rlistsq
@@ -19,6 +21,25 @@ module simulation
   integer,public, allocatable, dimension(:)   :: list
 
 
+
+
+
+  public :: check
+  public :: save_nlist
+  public :: wrap
+  public :: getBox
+
+  interface wrap
+     module procedure wrap_1d
+     module procedure wrap_3d
+  end interface
+
+  interface getBox
+     module procedure getBox_3d
+     module procedure getBox_dim
+  end interface
+
+
 !MPI dependent routines
 
 #ifdef IS_MPI
@@ -71,4 +92,97 @@ end module simulation
   end subroutine change_box_size
 
 
+  elemental function getBox_3d() result(thisBox)
+    real( kind = dp ), dimension(3) :: thisBox
+    thisBox = box
+  end function getBox_3d
+
+  function getBox_dim(dim) result(thisBox)
+    integer, intent(in) :: dim
+    real( kind = dp ) :: thisBox
+
+    thisBox = box(dim)
+  end function getBox_dim
+    
+  subroutine check(update_nlist)
+  
+    integer :: i
+    real( kind = DP ) :: dispmx 
+    logical, intent(out) :: update_nlist
+    real( kind = DP ) :: dispmx_tmp
+    
+    dispmx = 0.0E0_DP
+    
+    !! calculate the largest displacement of any atom in any direction
+    
+#ifdef MPI 
+    dispmx_tmp = 0.0E0_DP
+    do i = 1, nlocal
+       dispmx_tmp = max( abs ( q(1,i) - q0(1,i) ), dispmx )
+       dispmx_tmp = max( abs ( q(2,i) - q0(2,i) ), dispmx )
+       dispmx_tmp = max( abs ( q(3,i) - q0(3,i) ), dispmx )
+    end do
+    call mpi_allreduce(dispmx_tmp,dispmx,1,mpi_double_precision, & 
+       mpi_max,mpi_comm_world,mpi_err)
+#else
+    do i = 1, natoms
+       dispmx = max( abs ( q(1,i) - q0(1,i) ), dispmx )
+       dispmx = max( abs ( q(2,i) - q0(2,i) ), dispmx )
+       dispmx = max( abs ( q(3,i) - q0(3,i) ), dispmx )
+    end do
+#endif
+  
+    !! a conservative test of list skin crossings
+    dispmx = 2.0E0_DP * sqrt (3.0E0_DP * dispmx * dispmx)
+    
+    update_nlist = (dispmx.gt.(skin_thickness))
+    
+  end subroutine check
+  
+  subroutine save_nlist()
+    integer :: i
+#ifdef MPI
+    do i = 1, nlocal
+#else
+       do i = 1, natoms
+#endif
+          q0(1,i) = q(1,i)
+          q0(2,i) = q(2,i)
+          q0(3,i) = q(3,i)
+       end do
+
+  end subroutine save_nlist
+  
+
+  function wrap_1d(r,dim) result(this_wrap)
+    
+    
+    real( kind = DP ) :: r
+    real( kind = DP ) :: this_wrap
+    integer           :: dim
+    
+    if (use_pbc) then
+       !     this_wrap = r - box(dim)*dsign(1.0E0_DP,r)*int(abs(r/box(dim)) + 0.5E0_DP) 
+       this_wrap = r - box(dim)*nint(r/box(dim))
+    else
+       this_wrap = r
+    endif
+    
+    return
+  end function wrap_1d
+
+  elemental function wrap_3d(r) result(this_wrap)
+    real( kind = dp ), dimension(3), intent(in) :: r
+    real( kind = dp ), dimension(3) :: this_wrap
+
+    
+    if (use_pbc) then
+       !     this_wrap = r - box(dim)*dsign(1.0E0_DP,r)*int(abs(r/box(dim)) + 0.5E0_DP) 
+       this_wrap(1:3) = r(1:3) - box(1:3)*nint(r(1:3)/box(1:3))
+    else
+       this_wrap = r
+    endif
+  end function wrap_3d
+
+  
 end module simulation