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 891 by mmeineke, Fri Dec 19 20:36:35 2003 UTC vs.
Revision 897 by chuckv, Mon Jan 5 22:18:52 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.40 2003-12-19 20:36:35 mmeineke Exp $, $Date: 2003-12-19 20:36:35 $, $Name: not supported by cvs2svn $, $Revision: 1.40 $
7 > !! @version $Id: do_Forces.F90,v 1.43 2004-01-05 22:18:52 chuckv Exp $, $Date: 2004-01-05 22:18:52 $, $Name: not supported by cvs2svn $, $Revision: 1.43 $
8  
9   module do_Forces
10    use force_globals
# Line 52 | Line 52 | contains
52    integer :: nLoops
53   #endif
54  
55 +  logical, allocatable :: propertyMapI(:,:)
56 +  logical, allocatable :: propertyMapJ(:,:)
57 +
58   contains
59  
60    subroutine setRlistDF( this_rlist )
# Line 204 | 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 224 | Line 228 | contains
228      real( kind = DP ) ::  rijsq
229      real(kind=dp),dimension(3) :: d
230      real(kind=dp) :: rfpot, mu_i, virial
231 <    integer :: me_i
231 >    integer :: me_i, me_j
232      logical :: is_dp_i
233      integer :: neighborListSize
234      integer :: listerror, error
235      integer :: localError
236 +    integer :: propPack_i, propPack_j
237  
238      real(kind=dp) :: listSkin = 1.0  
239  
# Line 255 | Line 260 | contains
260      do_pot = do_pot_c
261      do_stress = do_stress_c
262  
263 <
263 >
264 > #ifdef IS_MPI
265 >    if (.not.allocated(propertyMapI)) then
266 >       allocate(propertyMapI(5,nrow))
267 >    endif
268 >
269 >    do i = 1, nrow
270 >       me_i = atid_row(i)
271 > #else
272 >    if (.not.allocated(propertyMapI)) then
273 >       allocate(propertyMapI(5,nlocal))
274 >    endif
275 >
276 >    do i = 1, natoms
277 >       me_i = atid(i)
278 > #endif
279 >      
280 >       propertyMapI(1:5,i) = .false.
281 >
282 >       call getElementProperty(atypes, me_i, "propertyPack", propPack_i)
283 >    
284 >       ! unpack the properties
285 >      
286 >       if (iand(propPack_i, LJ_PROPERTY_MASK) .eq. LJ_PROPERTY_MASK) &
287 >            propertyMapI(1, i) = .true.
288 >       if (iand(propPack_i, DP_PROPERTY_MASK) .eq. DP_PROPERTY_MASK) &
289 >            propertyMapI(2, i) = .true.
290 >       if (iand(propPack_i, STICKY_PROPERTY_MASK) .eq. STICKY_PROPERTY_MASK) &
291 >            propertyMapI(3, i) = .true.
292 >       if (iand(propPack_i, GB_PROPERTY_MASK) .eq. GB_PROPERTY_MASK) &
293 >            propertyMapI(4, i) = .true.
294 >       if (iand(propPack_i, EAM_PROPERTY_MASK) .eq. EAM_PROPERTY_MASK) &
295 >            propertyMapI(5, i) = .true.
296 >
297 >    end do
298 >
299 > #ifdef IS_MPI
300 >    if (.not.allocated(propertyMapJ)) then
301 >       allocate(propertyMapJ(5,ncol))
302 >    endif
303 >
304 >    do j = 1, ncol
305 >       me_j = atid_col(j)
306 > #else
307 >    if (.not.allocated(propertyMapJ)) then
308 >       allocate(propertyMapJ(5,nlocal))
309 >    endif
310 >
311 >    do j = 1, natoms
312 >       me_j = atid(j)
313 > #endif
314 >      
315 >       propertyMapJ(1:5,j) = .false.
316 >
317 >       call getElementProperty(atypes, me_j, "propertyPack", propPack_j)
318 >    
319 >       ! unpack the properties
320 >      
321 >       if (iand(propPack_j, LJ_PROPERTY_MASK) .eq. LJ_PROPERTY_MASK) &
322 >            propertyMapJ(1, j) = .true.
323 >       if (iand(propPack_j, DP_PROPERTY_MASK) .eq. DP_PROPERTY_MASK) &
324 >            propertyMapJ(2, j) = .true.
325 >       if (iand(propPack_j, STICKY_PROPERTY_MASK) .eq. STICKY_PROPERTY_MASK) &
326 >            propertyMapJ(3, j) = .true.
327 >       if (iand(propPack_j, GB_PROPERTY_MASK) .eq. GB_PROPERTY_MASK) &
328 >            propertyMapJ(4, j) = .true.
329 >       if (iand(propPack_j, EAM_PROPERTY_MASK) .eq. EAM_PROPERTY_MASK) &
330 >            propertyMapJ(5, j) = .true.
331 >
332 >    end do
333 >
334      ! Gather all information needed by all force loops:
335      
336   #ifdef IS_MPI    
# Line 451 | Line 526 | contains
526         nlist = 0      
527        
528         do i = 1, nrow
529 +
530            point(i) = nlist + 1
531            
532            inner: do j = 1, ncol
# Line 721 | Line 797 | contains
797      logical :: is_EAM_i,is_EAM_j
798      logical :: is_Sticky_i, is_Sticky_j
799      integer :: me_i, me_j
800 <
800 >    integer :: propPack_i
801 >    integer :: propPack_j
802      r = sqrt(rijsq)
803  
804   #ifdef IS_MPI
# Line 738 | Line 815 | contains
815      me_j = atid(j)
816  
817   #endif
818 <
818 >    
819      if (FF_uses_LJ .and. SimUsesLJ()) then
743       call getElementProperty(atypes, me_i, "is_LJ", is_LJ_i)
744       call getElementProperty(atypes, me_j, "is_LJ", is_LJ_j)
820  
821 <       if ( is_LJ_i .and. is_LJ_j ) &
821 >       if ( propertyMapI(1, me_i) .and. propertyMapJ(1, me_j) ) &
822              call do_lj_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
823 +
824      endif
825  
826      if (FF_uses_dipoles .and. SimUsesDipoles()) then
751       call getElementProperty(atypes, me_i, "is_DP", is_DP_i)
752       call getElementProperty(atypes, me_j, "is_DP", is_DP_j)
827        
828 <       if ( is_DP_i .and. is_DP_j ) then
828 >       if ( propertyMapI(2, me_i) .and. propertyMapJ(2, me_j)) then
829            call do_dipole_pair(i, j, d, r, rijsq, pot, u_l, f, t, &
830                 do_pot, do_stress)
831            if (FF_uses_RF .and. SimUsesRF()) then
# Line 764 | Line 838 | contains
838  
839      if (FF_uses_Sticky .and. SimUsesSticky()) then
840  
841 <       call getElementProperty(atypes, me_i, "is_Sticky", is_Sticky_i)
768 <       call getElementProperty(atypes, me_j, "is_Sticky", is_Sticky_j)
769 <
770 <       if ( is_Sticky_i .and. is_Sticky_j ) then
841 >       if ( propertyMapI(3, me_i) .and. propertyMapJ(3, me_j)) then
842            call do_sticky_pair(i, j, d, r, rijsq, A, pot, f, t, &
843                 do_pot, do_stress)
844         endif
# Line 775 | Line 846 | contains
846  
847  
848      if (FF_uses_GB .and. SimUsesGB()) then
778
779
780       call getElementProperty(atypes, me_i, "is_GB", is_GB_i)
781       call getElementProperty(atypes, me_j, "is_GB", is_GB_j)
849        
850 <       if ( is_GB_i .and. is_GB_j ) then
850 >       if ( propertyMapI(4, me_i) .and. propertyMapJ(4, me_j)) then
851            call do_gb_pair(i, j, d, r, rijsq, u_l, pot, f, t, &
852                 do_pot, do_stress)          
853         endif
854 +
855      endif
856      
857  
858    
859     if (FF_uses_EAM .and. SimUsesEAM()) then
792      call getElementProperty(atypes, me_i, "is_EAM", is_EAM_i)
793      call getElementProperty(atypes, me_j, "is_EAM", is_EAM_j)
860        
861 <      if ( is_EAM_i .and. is_EAM_j ) &
862 <           call do_eam_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
863 <   endif
798 <
799 <
861 >      if ( propertyMapI(5, me_i) .and. propertyMapJ(5, me_j)) then
862 >         call do_eam_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
863 >      endif
864  
865 +   endif
866  
867    end subroutine do_pair
868  

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines