ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/do_Forces.F90
(Generate patch)

Comparing trunk/OOPSE/libmdtools/do_Forces.F90 (file contents):
Revision 843 by mmeineke, Wed Oct 29 20:41:39 2003 UTC vs.
Revision 895 by chuckv, Mon Jan 5 22:12:11 2004 UTC

# Line 4 | Line 4
4  
5   !! @author Charles F. Vardeman II
6   !! @author Matthew Meineke
7 < !! @version $Id: do_Forces.F90,v 1.35 2003-10-29 20:41:39 mmeineke Exp $, $Date: 2003-10-29 20:41:39 $, $Name: not supported by cvs2svn $, $Revision: 1.35 $
7 > !! @version $Id: do_Forces.F90,v 1.42 2004-01-05 22:12:11 chuckv Exp $, $Date: 2004-01-05 22:12:11 $, $Name: not supported by cvs2svn $, $Revision: 1.42 $
8  
9   module do_Forces
10    use force_globals
# Line 46 | Line 46 | module do_Forces
46    public :: setRlistDF
47  
48   #ifdef PROFILE
49 <  real(kind = dp) :: forceTime
50 <  real(kind = dp) :: forceTimeInitial, forceTimeFinal
51 <  real(kind = dp) :: globalForceTime
52 <  real(kind = dp) :: maxForceTime
53 <  integer, save :: nloops = 0
49 >  public :: getforcetime
50 >  real, save :: forceTime = 0
51 >  real :: forceTimeInitial, forceTimeFinal
52 >  integer :: nLoops
53   #endif
54  
55 +  logical, allocatable :: propertyMapI(:,:)
56 +  logical, allocatable :: propertyMapJ(:,:)
57 +
58   contains
59  
60    subroutine setRlistDF( this_rlist )
# Line 205 | Line 207 | contains
207      real ( kind = dp ), dimension(3,getNlocal()) :: f
208      !! Torsion array provided by C, dimensioned by getNlocal
209      real( kind = dp ), dimension(3,getNlocal()) :: t    
210 +
211      !! Stress Tensor
212      real( kind = dp), dimension(9) :: tau  
213      real ( kind = dp ) :: pot
# Line 231 | Line 234 | contains
234      integer :: listerror, error
235      integer :: localError
236  
237 <    real(kind=dp) :: listSkin = 1.0
235 <    
237 >    real(kind=dp) :: listSkin = 1.0  
238  
239      !! initialize local variables  
240  
# Line 257 | Line 259 | contains
259      do_pot = do_pot_c
260      do_stress = do_stress_c
261  
262 <
262 >
263 > #ifdef IS_MPI
264 >    if (.not.allocated(propertyMapI)) then
265 >       allocate(propertyMapI(5,getNrow())
266 >    endif
267 >
268 >    do i = 1, nrow
269 >       me_i = atid_row(i)
270 > #else
271 >    if (.not.allocated(propertyMapI)) then
272 >       allocate(propertyMapI(5,getNlocal())
273 >    endif
274 >
275 >    do i = 1, natoms
276 >       me_i = atid(i)
277 > #endif
278 >      
279 >       propertyMapI(1:5,i) = .false.
280 >
281 >       call getElementProperty(atypes, me_i, "propertyPack", propPack_i)
282 >    
283 >       ! unpack the properties
284 >      
285 >       if (iand(propPack_i, LJ_PROPERTY_MASK) .eq. LJ_PROPERTY_MASK) &
286 >            propertyMapI(1, i) = .true.
287 >       if (iand(propPack_i, DP_PROPERTY_MASK) .eq. DP_PROPERTY_MASK) &
288 >            propertyMapI(2, i) = .true.
289 >       if (iand(propPack_i, STICKY_PROPERTY_MASK) .eq. STICKY_PROPERTY_MASK) &
290 >            propertyMapI(3, i) = .true.
291 >       if (iand(propPack_i, GB_PROPERTY_MASK) .eq. GB_PROPERTY_MASK) &
292 >            propertyMapI(4, i) = .true.
293 >       if (iand(propPack_i, EAM_PROPERTY_MASK) .eq. EAM_PROPERTY_MASK) &
294 >            propertyMapI(5, i) = .true.
295 >
296 >    end do
297 >
298 > #ifdef IS_MPI
299 >    if (.not.allocated(propertyMapJ)) then
300 >       allocate(propertyMapJ(5,getNcol())
301 >    endif
302 >
303 >    do j = 1, ncol
304 >       me_j = atid_col(j)
305 > #else
306 >    if (.not.allocated(propertyMapJ)) then
307 >       allocate(propertyMapJ(5,getNlocal())
308 >    endif
309 >
310 >    do j = 1, natoms
311 >       me_j = atid(j)
312 > #endif
313 >      
314 >       propertyMapJ(1:5,j) = .false.
315 >
316 >       call getElementProperty(atypes, me_j, "propertyPack", propPack_j)
317 >    
318 >       ! unpack the properties
319 >      
320 >       if (iand(propPack_j, LJ_PROPERTY_MASK) .eq. LJ_PROPERTY_MASK) &
321 >            propertyMapJ(1, j) = .true.
322 >       if (iand(propPack_j, DP_PROPERTY_MASK) .eq. DP_PROPERTY_MASK) &
323 >            propertyMapJ(2, j) = .true.
324 >       if (iand(propPack_j, STICKY_PROPERTY_MASK) .eq. STICKY_PROPERTY_MASK) &
325 >            propertyMapJ(3, j) = .true.
326 >       if (iand(propPack_j, GB_PROPERTY_MASK) .eq. GB_PROPERTY_MASK) &
327 >            propertyMapJ(4, j) = .true.
328 >       if (iand(propPack_j, EAM_PROPERTY_MASK) .eq. EAM_PROPERTY_MASK) &
329 >            propertyMapJ(5, j) = .true.
330 >
331 >    end do
332 >
333      ! Gather all information needed by all force loops:
334      
335   #ifdef IS_MPI    
# Line 445 | Line 517 | contains
517   #ifdef IS_MPI
518      
519      if (update_nlist) then
448 –      
520         !! save current configuration, construct neighbor list,
521         !! and calculate forces
522         call saveNeighborList(nlocal, q)
# Line 454 | Line 525 | contains
525         nlist = 0      
526        
527         do i = 1, nrow
528 +
529            point(i) = nlist + 1
530            
531            inner: do j = 1, ncol
# Line 511 | Line 583 | contains
583   #else
584      
585      if (update_nlist) then
586 <      
586 >
587         ! save current configuration, contruct neighbor list,
588         ! and calculate forces
589         call saveNeighborList(natoms, q)
# Line 698 | Line 770 | contains
770         tau = tau_Temp
771         virial = virial_Temp
772      endif
773 <
773 >    
774   #endif
775 <
776 < #ifdef PROFILE
777 <    if (do_pot) then
706 <
707 < #ifdef IS_MPI
708 <
709 <      
710 <       call printCommTime()
711 <
712 <       call mpi_allreduce(forceTime,globalForceTime,1,MPI_DOUBLE_PRECISION, &
713 <            mpi_sum,mpi_comm_world,mpi_err)
714 <
715 <       call mpi_allreduce(forceTime,maxForceTime,1,MPI_DOUBLE_PRECISION, &
716 <            MPI_MAX,mpi_comm_world,mpi_err)
717 <      
718 <       call mpi_comm_size( MPI_COMM_WORLD, nprocs,mpi_err)
719 <      
720 <       if (getMyNode() == 0) then
721 <          write(*,*) "Total processor time spent in force calculations is: ", globalForceTime
722 <          write(*,*) "Total Time spent in force loop per processor is: ", globalforceTime/nprocs
723 <          write(*,*) "Maximum force time on any processor is: ", maxForceTime
724 <       end if
725 < #else
726 <       write(*,*) "Time spent in force loop is: ", forceTime
727 < #endif
728 <
729 <    
730 <    endif
731 <
732 < #endif
733 <
775 >    
776 >    
777 >    
778    end subroutine do_force_loop
779  
780    subroutine do_pair(i, j, rijsq, d, do_pot, do_stress, u_l, A, f, t, pot)
# Line 752 | Line 796 | contains
796      logical :: is_EAM_i,is_EAM_j
797      logical :: is_Sticky_i, is_Sticky_j
798      integer :: me_i, me_j
799 <
799 >    integer :: propPack_i
800 >    integer :: propPack_j
801      r = sqrt(rijsq)
802  
803   #ifdef IS_MPI
# Line 769 | Line 814 | contains
814      me_j = atid(j)
815  
816   #endif
817 <
817 >    
818      if (FF_uses_LJ .and. SimUsesLJ()) then
819 <       call getElementProperty(atypes, me_i, "is_LJ", is_LJ_i)
820 <       call getElementProperty(atypes, me_j, "is_LJ", is_LJ_j)
776 <
777 <       if ( is_LJ_i .and. is_LJ_j ) &
819 >
820 >       if ( propertyMapI(1, me_i) .and. propertyMapJ(1, me_j) ) &
821              call do_lj_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
822 +
823      endif
824  
825      if (FF_uses_dipoles .and. SimUsesDipoles()) then
782 –       call getElementProperty(atypes, me_i, "is_DP", is_DP_i)
783 –       call getElementProperty(atypes, me_j, "is_DP", is_DP_j)
826        
827 <       if ( is_DP_i .and. is_DP_j ) then
827 >       if ( propertyMapI(2, me_i) .and. propertyMapJ(2, me_j)) then
828            call do_dipole_pair(i, j, d, r, rijsq, pot, u_l, f, t, &
829                 do_pot, do_stress)
830            if (FF_uses_RF .and. SimUsesRF()) then
# Line 795 | Line 837 | contains
837  
838      if (FF_uses_Sticky .and. SimUsesSticky()) then
839  
840 <       call getElementProperty(atypes, me_i, "is_Sticky", is_Sticky_i)
799 <       call getElementProperty(atypes, me_j, "is_Sticky", is_Sticky_j)
800 <
801 <       if ( is_Sticky_i .and. is_Sticky_j ) then
840 >       if ( propertyMapI(3, me_i) .and. propertyMapJ(3, me_j)) then
841            call do_sticky_pair(i, j, d, r, rijsq, A, pot, f, t, &
842                 do_pot, do_stress)
843         endif
# Line 806 | Line 845 | contains
845  
846  
847      if (FF_uses_GB .and. SimUsesGB()) then
809 –
810 –
811 –       call getElementProperty(atypes, me_i, "is_GB", is_GB_i)
812 –       call getElementProperty(atypes, me_j, "is_GB", is_GB_j)
848        
849 <       if ( is_GB_i .and. is_GB_j ) then
849 >       if ( propertyMapI(4, me_i) .and. propertyMapJ(4, me_j)) then
850            call do_gb_pair(i, j, d, r, rijsq, u_l, pot, f, t, &
851                 do_pot, do_stress)          
852         endif
853 +
854      endif
855      
856  
857    
858     if (FF_uses_EAM .and. SimUsesEAM()) then
823 –      call getElementProperty(atypes, me_i, "is_EAM", is_EAM_i)
824 –      call getElementProperty(atypes, me_j, "is_EAM", is_EAM_j)
859        
860 <      if ( is_EAM_i .and. is_EAM_j ) &
861 <           call do_eam_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
862 <   endif
829 <
830 <
860 >      if ( propertyMapI(5, me_i) .and. propertyMapJ(5, me_j)) then
861 >         call do_eam_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
862 >      endif
863  
864 +   endif
865  
866    end subroutine do_pair
867  
# Line 1120 | Line 1153 | contains
1153      doesit = FF_uses_RF
1154    end function FF_RequiresPostpairCalc
1155    
1156 + #ifdef PROFILE
1157 +  function getforcetime() result(totalforcetime)
1158 +    real(kind=dp) :: totalforcetime
1159 +    totalforcetime = forcetime
1160 +  end function getforcetime
1161 + #endif
1162 +
1163   !! This cleans componets of force arrays belonging only to fortran
1164  
1165   end module do_Forces

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines