| 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.25 2005-07-29 17:34:06 tim Exp $, $Date: 2005-07-29 17:34:06 $, $Name: not supported by cvs2svn $, $Revision: 1.25 $ |
| 49 |
|
|
| 50 |
|
|
| 51 |
|
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 |
| 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 |
| 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 |
| 157 |
|
real(kind=dp) :: myRcut |
| 158 |
< |
! Test Types |
| 158 |
> |
!! Test Types |
| 159 |
|
logical :: i_is_LJ |
| 160 |
|
logical :: i_is_Elect |
| 161 |
|
logical :: i_is_Sticky |
| 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!") |
| 247 |
|
end do |
| 248 |
|
|
| 249 |
|
end do |
| 250 |
+ |
|
| 251 |
+ |
haveInteractionMap = .true. |
| 252 |
|
end subroutine createInteractionMap |
| 253 |
|
|
| 254 |
+ |
! 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) |
| 256 |
+ |
real(kind=dp), intent(in), optional :: defaultRList |
| 257 |
+ |
integer :: iMap |
| 258 |
+ |
integer :: map_i,map_j |
| 259 |
+ |
real(kind=dp) :: thisRCut = 0.0_dp |
| 260 |
+ |
real(kind=dp) :: actualCutoff = 0.0_dp |
| 261 |
+ |
integer, intent(out) :: stat |
| 262 |
+ |
integer :: nAtypes |
| 263 |
+ |
integer :: myStatus |
| 264 |
+ |
|
| 265 |
+ |
stat = 0 |
| 266 |
+ |
if (.not. haveInteractionMap) then |
| 267 |
+ |
|
| 268 |
+ |
call createInteractionMap(myStatus) |
| 269 |
+ |
|
| 270 |
+ |
if (myStatus .ne. 0) then |
| 271 |
+ |
write(default_error, *) 'createInteractionMap failed in doForces!' |
| 272 |
+ |
stat = -1 |
| 273 |
+ |
return |
| 274 |
+ |
endif |
| 275 |
+ |
endif |
| 276 |
+ |
|
| 277 |
+ |
|
| 278 |
+ |
nAtypes = getSize(atypes) |
| 279 |
+ |
!! If we pass a default rcut, set all atypes to that cutoff distance |
| 280 |
+ |
if(present(defaultRList)) then |
| 281 |
+ |
InteractionMap(:,:)%rList = defaultRList |
| 282 |
+ |
InteractionMap(:,:)%rListSq = defaultRList*defaultRList |
| 283 |
+ |
haveRlist = .true. |
| 284 |
+ |
return |
| 285 |
+ |
end if |
| 286 |
+ |
|
| 287 |
+ |
do map_i = 1,nAtypes |
| 288 |
+ |
do map_j = map_i,nAtypes |
| 289 |
+ |
iMap = InteractionMap(map_i, map_j)%InteractionHash |
| 290 |
+ |
|
| 291 |
+ |
if ( iand(iMap, LJ_PAIR).ne.0 ) then |
| 292 |
+ |
! thisRCut = getLJCutOff(map_i,map_j) |
| 293 |
+ |
if (thisRcut > actualCutoff) actualCutoff = thisRcut |
| 294 |
+ |
endif |
| 295 |
+ |
|
| 296 |
+ |
if ( iand(iMap, ELECTROSTATIC_PAIR).ne.0 ) then |
| 297 |
+ |
! thisRCut = getElectrostaticCutOff(map_i,map_j) |
| 298 |
+ |
if (thisRcut > actualCutoff) actualCutoff = thisRcut |
| 299 |
+ |
endif |
| 300 |
+ |
|
| 301 |
+ |
if ( iand(iMap, STICKY_PAIR).ne.0 ) then |
| 302 |
+ |
! thisRCut = getStickyCutOff(map_i,map_j) |
| 303 |
+ |
if (thisRcut > actualCutoff) actualCutoff = thisRcut |
| 304 |
+ |
endif |
| 305 |
+ |
|
| 306 |
+ |
if ( iand(iMap, STICKYPOWER_PAIR).ne.0 ) then |
| 307 |
+ |
! thisRCut = getStickyPowerCutOff(map_i,map_j) |
| 308 |
+ |
if (thisRcut > actualCutoff) actualCutoff = thisRcut |
| 309 |
+ |
endif |
| 310 |
+ |
|
| 311 |
+ |
if ( iand(iMap, GAYBERNE_PAIR).ne.0 ) then |
| 312 |
+ |
! thisRCut = getGayberneCutOff(map_i,map_j) |
| 313 |
+ |
if (thisRcut > actualCutoff) actualCutoff = thisRcut |
| 314 |
+ |
endif |
| 315 |
+ |
|
| 316 |
+ |
if ( iand(iMap, GAYBERNE_LJ).ne.0 ) then |
| 317 |
+ |
! thisRCut = getGaybrneLJCutOff(map_i,map_j) |
| 318 |
+ |
if (thisRcut > actualCutoff) actualCutoff = thisRcut |
| 319 |
+ |
endif |
| 320 |
+ |
|
| 321 |
+ |
if ( iand(iMap, EAM_PAIR).ne.0 ) then |
| 322 |
+ |
! thisRCut = getEAMCutOff(map_i,map_j) |
| 323 |
+ |
if (thisRcut > actualCutoff) actualCutoff = thisRcut |
| 324 |
+ |
endif |
| 325 |
+ |
|
| 326 |
+ |
if ( iand(iMap, SHAPE_PAIR).ne.0 ) then |
| 327 |
+ |
! thisRCut = getShapeCutOff(map_i,map_j) |
| 328 |
+ |
if (thisRcut > actualCutoff) actualCutoff = thisRcut |
| 329 |
+ |
endif |
| 330 |
+ |
|
| 331 |
+ |
if ( iand(iMap, SHAPE_LJ).ne.0 ) then |
| 332 |
+ |
! thisRCut = getShapeLJCutOff(map_i,map_j) |
| 333 |
+ |
if (thisRcut > actualCutoff) actualCutoff = thisRcut |
| 334 |
+ |
endif |
| 335 |
+ |
InteractionMap(map_i, map_j)%rList = actualCutoff |
| 336 |
+ |
InteractionMap(map_i, map_j)%rListSq = actualCutoff * actualCutoff |
| 337 |
+ |
end do |
| 338 |
+ |
end do |
| 339 |
+ |
haveRlist = .true. |
| 340 |
+ |
end subroutine createRcuts |
| 341 |
|
|
| 342 |
|
|
| 343 |
|
!!! THIS GOES AWAY FOR SIZE DEPENDENT CUTOFF |
| 647 |
|
integer :: localError |
| 648 |
|
integer :: propPack_i, propPack_j |
| 649 |
|
integer :: loopStart, loopEnd, loop |
| 650 |
< |
|
| 650 |
> |
integer :: iMap |
| 651 |
|
real(kind=dp) :: listSkin = 1.0 |
| 652 |
|
|
| 653 |
|
!! initialize local variables |
| 739 |
|
#endif |
| 740 |
|
outer: do i = istart, iend |
| 741 |
|
|
| 742 |
+ |
#ifdef IS_MPI |
| 743 |
+ |
me_i = atid_row(i) |
| 744 |
+ |
#else |
| 745 |
+ |
me_i = atid(i) |
| 746 |
+ |
#endif |
| 747 |
+ |
|
| 748 |
|
if (update_nlist) point(i) = nlist + 1 |
| 749 |
|
|
| 750 |
|
n_in_i = groupStartRow(i+1) - groupStartRow(i) |
| 772 |
|
endif |
| 773 |
|
|
| 774 |
|
#ifdef IS_MPI |
| 775 |
+ |
me_j = atid_col(j) |
| 776 |
|
call get_interatomic_vector(q_group_Row(:,i), & |
| 777 |
|
q_group_Col(:,j), d_grp, rgrpsq) |
| 778 |
|
#else |
| 779 |
+ |
me_j = atid(j) |
| 780 |
|
call get_interatomic_vector(q_group(:,i), & |
| 781 |
|
q_group(:,j), d_grp, rgrpsq) |
| 782 |
|
#endif |
| 783 |
|
|
| 784 |
< |
if (rgrpsq < rlistsq) then |
| 784 |
> |
if (rgrpsq < InteractionMap(me_i,me_j)%rListsq) then |
| 785 |
|
if (update_nlist) then |
| 786 |
|
nlist = nlist + 1 |
| 787 |
|
|
| 1117 |
|
endif |
| 1118 |
|
|
| 1119 |
|
if ( iand(iMap, GAYBERNE_LJ).ne.0 ) then |
| 1120 |
< |
call do_gblj_pair(i, j, d, r, rijsq, sw, vpair, fpair, & |
| 1121 |
< |
pot, A, f, t, do_pot) |
| 1120 |
> |
! call do_gblj_pair(i, j, d, r, rijsq, sw, vpair, fpair, & |
| 1121 |
> |
! pot, A, f, t, do_pot) |
| 1122 |
|
endif |
| 1123 |
|
|
| 1124 |
|
if ( iand(iMap, EAM_PAIR).ne.0 ) then |
