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

Comparing trunk/OOPSE-2.0/src/UseTheForce/doForces.F90 (file contents):
Revision 2261 by gezelter, Tue Jun 28 13:58:45 2005 UTC vs.
Revision 2266 by chuckv, Thu Jul 28 22:12:45 2005 UTC

# Line 45 | Line 45
45  
46   !! @author Charles F. Vardeman II
47   !! @author Matthew Meineke
48 < !! @version $Id: doForces.F90,v 1.22 2005-06-28 13:58:45 gezelter Exp $, $Date: 2005-06-28 13:58:45 $, $Name: not supported by cvs2svn $, $Revision: 1.22 $
48 > !! @version $Id: doForces.F90,v 1.24 2005-07-28 22:12:45 chuckv Exp $, $Date: 2005-07-28 22:12:45 $, $Name: not supported by cvs2svn $, $Revision: 1.24 $
49  
50  
51   module doForces
# Line 82 | Line 82 | module doForces
82    logical, save :: haveNeighborList = .false.
83    logical, save :: haveSIMvariables = .false.
84    logical, save :: haveSaneForceField = .false.
85 +  logical, save :: haveInteractionMap = .false.
86  
87    logical, save :: FF_uses_DirectionalAtoms
88    logical, save :: FF_uses_LennardJones
# Line 121 | Line 122 | module doForces
122    public :: init_FF
123    public :: do_force_loop
124   !  public :: setRlistDF
125 +  !public :: addInteraction
126 +  !public :: setInteractionHash
127 +  !public :: getInteractionHash
128 +  public :: createInteractionMap
129 +  public :: createRcuts
130  
131   #ifdef PROFILE
132    public :: getforcetime
# Line 131 | Line 137 | module doForces
137  
138    type, public :: Interaction
139       integer :: InteractionHash
140 <     real(kind=dp) :: rCut
140 >     real(kind=dp) :: rList = 0.0_dp
141 >     real(kind=dp) :: rListSq = 0.0_dp
142    end type Interaction
143    
144    type(Interaction), dimension(:,:),allocatable :: InteractionMap
145    
146 <  !public :: addInteraction
140 <  !public :: setInteractionHash
141 <  !public :: getInteractionHash
142 <  public :: createInteractionMap
146 >
147    
148   contains
149  
150  
151    subroutine createInteractionMap(status)
152      integer :: nAtypes
153 <    integer :: status
153 >    integer, intent(out) :: status
154      integer :: i
155      integer :: j
156      integer :: ihash
# Line 167 | Line 171 | contains
171      logical :: j_is_EAM
172      logical :: j_is_Shape
173      
174 +    status = 0
175      
176      if (.not. associated(atypes)) then
177         call handleError("atype", "atypes was not present before call of createDefaultInteractionMap!")
# Line 244 | Line 249 | contains
249      end do
250    end subroutine createInteractionMap
251  
252 + ! 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.
253 +  subroutine createRcuts(defaultRList,stat)
254 +    real(kind=dp), intent(in), optional :: defaultRList
255 +    integer :: iMap
256 +    integer :: map_i,map_j
257 +    real(kind=dp) :: thisRCut = 0.0_dp
258 +    real(kind=dp) :: actualCutoff = 0.0_dp
259 +    integer, intent(out) :: stat
260 +    integer :: nAtypes
261 +    integer :: myStatus
262 +
263 +    stat = 0
264 +    if (.not. haveInteractionMap) then
265 +
266 +       call createInteractionMap(myStatus)
267 +
268 +       if (myStatus .ne. 0) then
269 +          write(default_error, *) 'createInteractionMap failed in doForces!'
270 +          stat = -1
271 +          return
272 +       endif
273 +    endif
274 +
275 +
276 +    nAtypes = getSize(atypes)
277 + ! If we pass a default rcut, set all atypes to that cutoff distance
278 +    if(present(defaultRList)) then
279 +       InteractionMap(:,:)%rList = defaultRList
280 +       InteractionMap(:,:)%rListSq = defaultRList*defaultRList
281 +       haveRlist = .true.
282 +       return
283 +    end if
284 +
285 +    do map_i = 1,nAtypes
286 +       do map_j = map_i,nAtypes
287 +          iMap = InteractionMap(map_i, map_j)%InteractionHash
288 +          
289 +          if ( iand(iMap, LJ_PAIR).ne.0 ) then
290 + !            thisRCut = getLJCutOff(map_i,map_j)
291 +             if (thisRcut > actualCutoff) actualCutoff = thisRcut
292 +          endif
293 +          
294 +          if ( iand(iMap, ELECTROSTATIC_PAIR).ne.0 ) then
295 + !            thisRCut = getElectrostaticCutOff(map_i,map_j)
296 +             if (thisRcut > actualCutoff) actualCutoff = thisRcut
297 +          endif
298 +          
299 +          if ( iand(iMap, STICKY_PAIR).ne.0 ) then
300 + !             thisRCut = getStickyCutOff(map_i,map_j)
301 +              if (thisRcut > actualCutoff) actualCutoff = thisRcut
302 +           endif
303 +          
304 +           if ( iand(iMap, STICKYPOWER_PAIR).ne.0 ) then
305 + !              thisRCut = getStickyPowerCutOff(map_i,map_j)
306 +              if (thisRcut > actualCutoff) actualCutoff = thisRcut
307 +           endif
308 +          
309 +           if ( iand(iMap, GAYBERNE_PAIR).ne.0 ) then
310 + !              thisRCut = getGayberneCutOff(map_i,map_j)
311 +              if (thisRcut > actualCutoff) actualCutoff = thisRcut
312 +           endif
313 +          
314 +           if ( iand(iMap, GAYBERNE_LJ).ne.0 ) then
315 + !              thisRCut = getGaybrneLJCutOff(map_i,map_j)
316 +              if (thisRcut > actualCutoff) actualCutoff = thisRcut
317 +           endif
318 +          
319 +           if ( iand(iMap, EAM_PAIR).ne.0 ) then      
320 + !              thisRCut = getEAMCutOff(map_i,map_j)
321 +              if (thisRcut > actualCutoff) actualCutoff = thisRcut
322 +           endif
323 +          
324 +           if ( iand(iMap, SHAPE_PAIR).ne.0 ) then      
325 + !              thisRCut = getShapeCutOff(map_i,map_j)
326 +              if (thisRcut > actualCutoff) actualCutoff = thisRcut
327 +           endif
328 +          
329 +           if ( iand(iMap, SHAPE_LJ).ne.0 ) then      
330 + !              thisRCut = getShapeLJCutOff(map_i,map_j)
331 +              if (thisRcut > actualCutoff) actualCutoff = thisRcut
332 +           endif
333 +           InteractionMap(map_i, map_j)%rList = actualCutoff
334 +           InteractionMap(map_i, map_j)%rListSq = actualCutoff * actualCutoff
335 +        end do
336 +     end do
337 +          haveRlist = .true.
338 +  end subroutine createRcuts
339  
340  
341   !!! THIS GOES AWAY FOR SIZE DEPENDENT CUTOFF
# Line 553 | Line 645 | contains
645      integer :: localError
646      integer :: propPack_i, propPack_j
647      integer :: loopStart, loopEnd, loop
648 <
648 >    integer :: iMap
649      real(kind=dp) :: listSkin = 1.0  
650  
651      !! initialize local variables  
# Line 651 | Line 743 | contains
743  
744            if (update_nlist) then
745   #ifdef IS_MPI
746 +             me_i = atid_row(i)
747               jstart = 1
748               jend = nGroupsInCol
749   #else
750 +             me_i = atid(i)
751               jstart = i+1
752               jend = nGroups
753   #endif
# Line 672 | Line 766 | contains
766               endif
767  
768   #ifdef IS_MPI
769 +             me_j = atid_col(j)
770               call get_interatomic_vector(q_group_Row(:,i), &
771                    q_group_Col(:,j), d_grp, rgrpsq)
772   #else
773 +             me_j = atid(j)
774               call get_interatomic_vector(q_group(:,i), &
775                    q_group(:,j), d_grp, rgrpsq)
776   #endif
777  
778 <             if (rgrpsq < rlistsq) then
778 >             if (rgrpsq < InteractionMap(me_i,me_j)%rListsq) then
779                  if (update_nlist) then
780                     nlist = nlist + 1
781  
# Line 1015 | Line 1111 | contains
1111      endif
1112      
1113      if ( iand(iMap, GAYBERNE_LJ).ne.0 ) then
1114 <       call do_gblj_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1115 <            pot, A, f, t, do_pot)
1114 > !      call do_gblj_pair(i, j, d, r, rijsq, sw, vpair, fpair, &
1115 > !           pot, A, f, t, do_pot)
1116      endif
1117  
1118      if ( iand(iMap, EAM_PAIR).ne.0 ) then      

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines