ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE-4/src/UseTheForce/doForces.F90
(Generate patch)

Comparing trunk/OOPSE-4/src/UseTheForce/doForces.F90 (file contents):
Revision 2262 by chuckv, Sun Jul 3 20:53:43 2005 UTC vs.
Revision 2268 by gezelter, Fri Jul 29 19:38:27 2005 UTC

# Line 45 | Line 45
45  
46   !! @author Charles F. Vardeman II
47   !! @author Matthew Meineke
48 < !! @version $Id: doForces.F90,v 1.23 2005-07-03 20:53:43 chuckv Exp $, $Date: 2005-07-03 20:53:43 $, $Name: not supported by cvs2svn $, $Revision: 1.23 $
48 > !! @version $Id: doForces.F90,v 1.26 2005-07-29 19:38:27 gezelter Exp $, $Date: 2005-07-29 19:38:27 $, $Name: not supported by cvs2svn $, $Revision: 1.26 $
49  
50  
51   module doForces
# Line 137 | Line 137 | module doForces
137  
138    type, public :: Interaction
139       integer :: InteractionHash
140 <     real(kind=dp) :: rList = 0.0_dp
140 >     real(kind=dp) :: rCut = 0.0_dp
141 >     real(kind=dp) :: rCutSq = 0.0_dp    
142       real(kind=dp) :: rListSq = 0.0_dp
143    end type Interaction
144    
# Line 150 | Line 151 | contains
151  
152    subroutine createInteractionMap(status)
153      integer :: nAtypes
154 <    integer :: status
154 >    integer, intent(out) :: status
155      integer :: i
156      integer :: j
157      integer :: ihash
158      real(kind=dp) :: myRcut
159 < ! Test Types
159 >    !! Test Types
160      logical :: i_is_LJ
161      logical :: i_is_Elect
162      logical :: i_is_Sticky
# Line 171 | Line 172 | contains
172      logical :: j_is_EAM
173      logical :: j_is_Shape
174      
175 <    
175 >    status = 0  
176 >
177      if (.not. associated(atypes)) then
178         call handleError("atype", "atypes was not present before call of createDefaultInteractionMap!")
179         status = -1
# Line 246 | Line 248 | contains
248         end do
249  
250      end do
251 +
252 +    haveInteractionMap = .true.
253    end subroutine createInteractionMap
254  
255 < ! Query each potential and return the cutoff for that potential. We build the neighbor list based on the largest cutoff value for that atype. Each potential can decide whether to calculate the force for that atype based upon it's own cutoff.
256 <  subroutine createRcuts(defaultRList)
257 <    real(kind=dp), intent(in), optional :: defaultRList
255 >  ! Query each potential and return the cutoff for that potential. We
256 >  ! build the neighbor list based on the largest cutoff value for that
257 >  ! atype. Each potential can decide whether to calculate the force for
258 >  ! that atype based upon it's own cutoff.
259 >  
260 >  subroutine createRcuts(defaultRcut, defaultSkinThickness, stat)
261 >
262 >    real(kind=dp), intent(in), optional :: defaultRCut, defaultSkinThickness
263      integer :: iMap
264      integer :: map_i,map_j
265      real(kind=dp) :: thisRCut = 0.0_dp
266      real(kind=dp) :: actualCutoff = 0.0_dp
267 +    integer, intent(out) :: stat
268      integer :: nAtypes
269 +    integer :: myStatus
270  
271 <    if(.not. allocated(InteractionMap)) return
271 >    stat = 0
272 >    if (.not. haveInteractionMap) then
273  
274 +       call createInteractionMap(myStatus)
275 +
276 +       if (myStatus .ne. 0) then
277 +          write(default_error, *) 'createInteractionMap failed in doForces!'
278 +          stat = -1
279 +          return
280 +       endif
281 +    endif
282 +
283      nAtypes = getSize(atypes)
284 < ! If we pass a default rcut, set all atypes to that cutoff distance
284 >    !! If we pass a default rcut, set all atypes to that cutoff distance
285      if(present(defaultRList)) then
286 <       InteractionMap(:,:)%rList = defaultRList
287 <       InteractionMap(:,:)%rListSq = defaultRList*defaultRList
286 >       InteractionMap(:,:)%rCut = defaultRCut
287 >       InteractionMap(:,:)%rCutSq = defaultRCut*defaultRCut
288 >       InteractionMap(:,:)%rListSq = (defaultRCut+defaultSkinThickness)**2
289         haveRlist = .true.
290         return
291      end if
# Line 273 | Line 295 | contains
295            iMap = InteractionMap(map_i, map_j)%InteractionHash
296            
297            if ( iand(iMap, LJ_PAIR).ne.0 ) then
298 < !            thisRCut = getLJCutOff(map_i,map_j)
298 >             ! thisRCut = getLJCutOff(map_i,map_j)
299               if (thisRcut > actualCutoff) actualCutoff = thisRcut
300            endif
301            
302            if ( iand(iMap, ELECTROSTATIC_PAIR).ne.0 ) then
303 < !            thisRCut = getElectrostaticCutOff(map_i,map_j)
303 >             ! thisRCut = getElectrostaticCutOff(map_i,map_j)
304               if (thisRcut > actualCutoff) actualCutoff = thisRcut
305            endif
306            
307            if ( iand(iMap, STICKY_PAIR).ne.0 ) then
308 < !             thisRCut = getStickyCutOff(map_i,map_j)
308 >             ! thisRCut = getStickyCutOff(map_i,map_j)
309                if (thisRcut > actualCutoff) actualCutoff = thisRcut
310             endif
311            
312             if ( iand(iMap, STICKYPOWER_PAIR).ne.0 ) then
313 < !              thisRCut = getStickyPowerCutOff(map_i,map_j)
313 >              ! thisRCut = getStickyPowerCutOff(map_i,map_j)
314                if (thisRcut > actualCutoff) actualCutoff = thisRcut
315             endif
316            
317             if ( iand(iMap, GAYBERNE_PAIR).ne.0 ) then
318 < !              thisRCut = getGayberneCutOff(map_i,map_j)
318 >              ! thisRCut = getGayberneCutOff(map_i,map_j)
319                if (thisRcut > actualCutoff) actualCutoff = thisRcut
320             endif
321            
# Line 316 | Line 338 | contains
338   !              thisRCut = getShapeLJCutOff(map_i,map_j)
339                if (thisRcut > actualCutoff) actualCutoff = thisRcut
340             endif
341 <           InteractionMap(map_i, map_j)%rList = actualCutoff
342 <           InteractionMap(map_i, map_j)%rListSq = actualCutoff * actualCutoff
341 >           InteractionMap(map_i, map_j)%rCut = actualCutoff
342 >           InteractionMap(map_i, map_j)%rCutSq = actualCutoff * actualCutoff
343 >           InteractionMap(map_i, map_j)%rListSq = (actualCutoff + skinThickness)**2
344 >
345 >           InteractionMap(map_j, map_i)%rCut = InteractionMap(map_i, map_j)%rCut
346 >           InteractionMap(map_j, map_i)%rCutSq = InteractionMap(map_i, map_j)%rCutSq
347 >           InteractionMap(map_j, map_i)%rListSq = InteractionMap(map_i, map_j)%rListSq
348          end do
349       end do
350 <          haveRlist = .true.
350 >     ! now the groups
351 >
352 >
353 >
354 >     haveRlist = .true.
355    end subroutine createRcuts
356  
357  
# Line 367 | Line 398 | contains
398      error = 0
399  
400      if (.not. haveInteractionMap) then
401 <
402 <       myStatus = 0
372 <
401 >      
402 >       myStatus = 0      
403         call createInteractionMap(myStatus)
404 <
404 >      
405         if (myStatus .ne. 0) then
406            write(default_error, *) 'createInteractionMap failed in doForces!'
407            error = -1
# Line 723 | Line 753 | contains
753   #endif
754         outer: do i = istart, iend
755  
756 + #ifdef IS_MPI
757 +             me_i = atid_row(i)
758 + #else
759 +             me_i = atid(i)
760 + #endif
761 +
762            if (update_nlist) point(i) = nlist + 1
763  
764            n_in_i = groupStartRow(i+1) - groupStartRow(i)
# Line 750 | Line 786 | contains
786               endif
787  
788   #ifdef IS_MPI
789 +             me_j = atid_col(j)
790               call get_interatomic_vector(q_group_Row(:,i), &
791                    q_group_Col(:,j), d_grp, rgrpsq)
792   #else
793 +             me_j = atid(j)
794               call get_interatomic_vector(q_group(:,i), &
795                    q_group(:,j), d_grp, rgrpsq)
796   #endif

Diff Legend

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