45 |
|
|
46 |
|
!! @author Charles F. Vardeman II |
47 |
|
!! @author Matthew Meineke |
48 |
< |
!! @version $Id: doForces.F90,v 1.25 2005-07-29 17:34:06 tim Exp $, $Date: 2005-07-29 17:34:06 $, $Name: not supported by cvs2svn $, $Revision: 1.25 $ |
48 |
> |
!! @version $Id: doForces.F90,v 1.27 2005-08-09 19:40:56 chuckv Exp $, $Date: 2005-08-09 19:40:56 $, $Name: not supported by cvs2svn $, $Revision: 1.27 $ |
49 |
|
|
50 |
|
|
51 |
|
module doForces |
116 |
|
logical, save :: SIM_uses_PBC |
117 |
|
logical, save :: SIM_uses_molecular_cutoffs |
118 |
|
|
119 |
– |
!!!GO AWAY--------- |
120 |
– |
!!!!!real(kind=dp), save :: rlist, rlistsq |
119 |
|
|
120 |
|
public :: init_FF |
121 |
|
public :: do_force_loop |
124 |
|
!public :: setInteractionHash |
125 |
|
!public :: getInteractionHash |
126 |
|
public :: createInteractionMap |
127 |
< |
public :: createRcuts |
127 |
> |
public :: createGroupCutoffs |
128 |
|
|
129 |
|
#ifdef PROFILE |
130 |
|
public :: getforcetime |
132 |
|
real :: forceTimeInitial, forceTimeFinal |
133 |
|
integer :: nLoops |
134 |
|
#endif |
137 |
– |
|
138 |
– |
type, public :: Interaction |
139 |
– |
integer :: InteractionHash |
140 |
– |
real(kind=dp) :: rList = 0.0_dp |
141 |
– |
real(kind=dp) :: rListSq = 0.0_dp |
142 |
– |
end type Interaction |
135 |
|
|
136 |
< |
type(Interaction), dimension(:,:),allocatable :: InteractionMap |
137 |
< |
|
136 |
> |
!! Variables for cutoff mapping and interaction mapping |
137 |
> |
! Bit hash to determine pair-pair interactions. |
138 |
> |
integer, dimension(:,:),allocatable :: InteractionHash |
139 |
> |
!! Cuttoffs in OOPSE are handled on a Group-Group pair basis. |
140 |
> |
! Largest cutoff for atypes for all potentials |
141 |
> |
real(kind=dp), dimension(:), allocatable :: atypeMaxCuttoff |
142 |
> |
! Largest cutoff for groups |
143 |
> |
real(kind=dp), dimension(:), allocatable :: groupMaxCutoff |
144 |
> |
! Group to Gtype transformation Map |
145 |
> |
integer,dimension(:), allocatable :: groupToGtype |
146 |
> |
! Group Type Max Cutoff |
147 |
> |
real(kind=dp), dimension(:), allocatable :: gtypeMaxCutoff |
148 |
> |
! GroupType definition |
149 |
> |
type ::gtype |
150 |
> |
real(kind=dp) :: rcut ! Group Cutoff |
151 |
> |
real(kind=dp) :: rcutsq ! Group Cutoff Squared |
152 |
> |
real(kind=dp) :: rlistsq ! List cutoff Squared |
153 |
> |
end type gtype |
154 |
|
|
155 |
+ |
type(gtype), dimension(:,:), allocatable :: gtypeCutoffMap |
156 |
|
|
157 |
|
contains |
158 |
|
|
180 |
|
logical :: j_is_EAM |
181 |
|
logical :: j_is_Shape |
182 |
|
|
183 |
< |
status = 0 |
184 |
< |
|
183 |
> |
status = 0 |
184 |
> |
|
185 |
|
if (.not. associated(atypes)) then |
186 |
< |
call handleError("atype", "atypes was not present before call of createDefaultInteractionMap!") |
186 |
> |
call handleError("atype", "atypes was not present before call of createDefaultInteractionHash!") |
187 |
|
status = -1 |
188 |
|
return |
189 |
|
endif |
195 |
|
return |
196 |
|
end if |
197 |
|
|
198 |
< |
if (.not. allocated(InteractionMap)) then |
199 |
< |
allocate(InteractionMap(nAtypes,nAtypes)) |
198 |
> |
if (.not. allocated(InteractionHash)) then |
199 |
> |
allocate(InteractionHash(nAtypes,nAtypes)) |
200 |
|
endif |
201 |
|
|
202 |
|
do i = 1, nAtypes |
250 |
|
if (i_is_LJ .and. j_is_Shape) iHash = ior(iHash, SHAPE_LJ) |
251 |
|
|
252 |
|
|
253 |
< |
InteractionMap(i,j)%InteractionHash = iHash |
254 |
< |
InteractionMap(j,i)%InteractionHash = iHash |
253 |
> |
InteractionHash(i,j) = iHash |
254 |
> |
InteractionHash(j,i) = iHash |
255 |
|
|
256 |
|
end do |
257 |
|
|
260 |
|
haveInteractionMap = .true. |
261 |
|
end subroutine createInteractionMap |
262 |
|
|
263 |
< |
! 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. |
255 |
< |
subroutine createRcuts(defaultRList,stat) |
263 |
> |
subroutine createGroupCutoffs(skinThickness,defaultrList,stat) |
264 |
|
real(kind=dp), intent(in), optional :: defaultRList |
265 |
+ |
real(kind-dp), intent(in), :: skinThickenss |
266 |
+ |
! Query each potential and return the cutoff for that potential. We |
267 |
+ |
! build the neighbor list based on the largest cutoff value for that |
268 |
+ |
! atype. Each potential can decide whether to calculate the force for |
269 |
+ |
! that atype based upon it's own cutoff. |
270 |
+ |
|
271 |
+ |
|
272 |
+ |
real(kind=dp), intent(in), optional :: defaultRCut, defaultSkinThickness |
273 |
+ |
|
274 |
|
integer :: iMap |
275 |
|
integer :: map_i,map_j |
276 |
|
real(kind=dp) :: thisRCut = 0.0_dp |
291 |
|
endif |
292 |
|
endif |
293 |
|
|
277 |
– |
|
294 |
|
nAtypes = getSize(atypes) |
295 |
|
!! If we pass a default rcut, set all atypes to that cutoff distance |
296 |
|
if(present(defaultRList)) then |
297 |
< |
InteractionMap(:,:)%rList = defaultRList |
298 |
< |
InteractionMap(:,:)%rListSq = defaultRList*defaultRList |
297 |
> |
InteractionMap(:,:)%rCut = defaultRCut |
298 |
> |
InteractionMap(:,:)%rCutSq = defaultRCut*defaultRCut |
299 |
> |
InteractionMap(:,:)%rListSq = (defaultRCut+defaultSkinThickness)**2 |
300 |
|
haveRlist = .true. |
301 |
|
return |
302 |
|
end if |
349 |
|
! thisRCut = getShapeLJCutOff(map_i,map_j) |
350 |
|
if (thisRcut > actualCutoff) actualCutoff = thisRcut |
351 |
|
endif |
352 |
< |
InteractionMap(map_i, map_j)%rList = actualCutoff |
353 |
< |
InteractionMap(map_i, map_j)%rListSq = actualCutoff * actualCutoff |
352 |
> |
InteractionMap(map_i, map_j)%rCut = actualCutoff |
353 |
> |
InteractionMap(map_i, map_j)%rCutSq = actualCutoff * actualCutoff |
354 |
> |
InteractionMap(map_i, map_j)%rListSq = (actualCutoff + skinThickness)**2 |
355 |
> |
|
356 |
> |
InteractionMap(map_j, map_i)%rCut = InteractionMap(map_i, map_j)%rCut |
357 |
> |
InteractionMap(map_j, map_i)%rCutSq = InteractionMap(map_i, map_j)%rCutSq |
358 |
> |
InteractionMap(map_j, map_i)%rListSq = InteractionMap(map_i, map_j)%rListSq |
359 |
|
end do |
360 |
|
end do |
361 |
< |
haveRlist = .true. |
340 |
< |
end subroutine createRcuts |
361 |
> |
! now the groups |
362 |
|
|
363 |
|
|
343 |
– |
!!! THIS GOES AWAY FOR SIZE DEPENDENT CUTOFF |
344 |
– |
!!$ subroutine setRlistDF( this_rlist ) |
345 |
– |
!!$ |
346 |
– |
!!$ real(kind=dp) :: this_rlist |
347 |
– |
!!$ |
348 |
– |
!!$ rlist = this_rlist |
349 |
– |
!!$ rlistsq = rlist * rlist |
350 |
– |
!!$ |
351 |
– |
!!$ haveRlist = .true. |
352 |
– |
!!$ |
353 |
– |
!!$ end subroutine setRlistDF |
364 |
|
|
365 |
+ |
haveRlist = .true. |
366 |
+ |
end subroutine createGroupCutoffs |
367 |
|
|
368 |
|
subroutine setSimVariables() |
369 |
|
SIM_uses_DirectionalAtoms = SimUsesDirectionalAtoms() |
395 |
|
error = 0 |
396 |
|
|
397 |
|
if (.not. haveInteractionMap) then |
398 |
< |
|
399 |
< |
myStatus = 0 |
388 |
< |
|
398 |
> |
|
399 |
> |
myStatus = 0 |
400 |
|
call createInteractionMap(myStatus) |
401 |
< |
|
401 |
> |
|
402 |
|
if (myStatus .ne. 0) then |
403 |
|
write(default_error, *) 'createInteractionMap failed in doForces!' |
404 |
|
error = -1 |
1010 |
|
#else |
1011 |
|
me_i = atid(i) |
1012 |
|
#endif |
1013 |
< |
iMap = InteractionMap(me_i, me_j)%InteractionHash |
1013 |
> |
iMap = InteractionHash(me_i,me_j) |
1014 |
|
|
1015 |
|
if ( iand(iMap, ELECTROSTATIC_PAIR).ne.0 ) then |
1016 |
|
|