49 |
|
use status |
50 |
|
use atype_module |
51 |
|
use vector_class |
52 |
+ |
use fForceOptions |
53 |
|
#ifdef IS_MPI |
54 |
|
use mpiSimulation |
55 |
|
#endif |
66 |
|
logical, save :: haveRcut = .false. |
67 |
|
logical, save :: haveMixingMap = .false. |
68 |
|
logical, save :: useGeometricDistanceMixing = .false. |
69 |
+ |
logical, save :: cleanArrays = .true. |
70 |
+ |
logical, save :: arraysAllocated = .false. |
71 |
|
|
72 |
|
|
70 |
– |
|
71 |
– |
|
73 |
|
character(len = statusMsgSize) :: errMesg |
74 |
|
integer :: sc_err |
75 |
|
|
76 |
|
character(len = 200) :: errMsg |
77 |
|
character(len=*), parameter :: RoutineName = "Sutton-Chen MODULE" |
78 |
|
!! Logical that determines if eam arrays should be zeroed |
78 |
– |
logical :: cleanme = .true. |
79 |
|
logical :: nmflag = .false. |
80 |
|
|
81 |
|
|
140 |
|
public :: do_SC_pair |
141 |
|
public :: newSCtype |
142 |
|
public :: calc_SC_prepair_rho |
143 |
+ |
public :: calc_SC_preforce_Frho |
144 |
|
public :: clean_SC |
145 |
|
public :: destroySCtypes |
146 |
|
public :: getSCCut |
147 |
|
! public :: setSCDefaultCutoff |
148 |
|
! public :: setSCUniformCutoff |
149 |
< |
public :: useGeometricMixing |
149 |
> |
|
150 |
|
|
151 |
|
contains |
152 |
|
|
176 |
|
|
177 |
|
! check to see if this is the first time into |
178 |
|
if (.not.associated(SCList%SCTypes)) then |
179 |
< |
call getMatchingElementList(atypes, "is_SuttonChen", .true., nSCtypes, MatchList) |
179 |
> |
call getMatchingElementList(atypes, "is_SC", .true., nSCtypes, MatchList) |
180 |
|
SCList%nSCtypes = nSCtypes |
181 |
|
allocate(SCList%SCTypes(nSCTypes)) |
182 |
|
nAtypes = getSize(atypes) |
221 |
|
real(kind=dp) :: cutValue |
222 |
|
|
223 |
|
scID = SCList%atidToSCType(atomID) |
224 |
< |
cutValue = SCList%SCTypes(scID)%sc_rcut |
224 |
> |
cutValue = 2.0_dp * SCList%SCTypes(scID)%alpha |
225 |
|
end function getSCCut |
226 |
|
|
227 |
|
|
243 |
|
if (.not. allocated(MixingMap)) then |
244 |
|
allocate(MixingMap(nSCtypes, nSCtypes)) |
245 |
|
endif |
246 |
< |
|
246 |
> |
useGeometricDistanceMixing = usesGeometricDistanceMixing() |
247 |
|
do i = 1, nSCtypes |
248 |
|
|
249 |
|
e1 = SCList%SCtypes(i)%epsilon |
290 |
|
|
291 |
|
!! routine checks to see if array is allocated, deallocates array if allocated |
292 |
|
!! and then creates the array to the required size |
293 |
< |
subroutine allocateSC(status) |
294 |
< |
integer, intent(out) :: status |
293 |
> |
subroutine allocateSC() |
294 |
> |
integer :: status |
295 |
|
|
296 |
|
#ifdef IS_MPI |
297 |
|
integer :: nAtomsInRow |
299 |
|
#endif |
300 |
|
integer :: alloc_stat |
301 |
|
|
302 |
< |
|
302 |
> |
|
303 |
|
status = 0 |
304 |
|
#ifdef IS_MPI |
305 |
|
nAtomsInRow = getNatomsInRow(plan_atom_row) |
312 |
|
allocate(frho(nlocal),stat=alloc_stat) |
313 |
|
if (alloc_stat /= 0) then |
314 |
|
status = -1 |
314 |
– |
return |
315 |
|
end if |
316 |
|
|
317 |
|
if (allocated(rho)) deallocate(rho) |
318 |
|
allocate(rho(nlocal),stat=alloc_stat) |
319 |
|
if (alloc_stat /= 0) then |
320 |
|
status = -1 |
321 |
– |
return |
321 |
|
end if |
322 |
|
|
323 |
|
if (allocated(dfrhodrho)) deallocate(dfrhodrho) |
324 |
|
allocate(dfrhodrho(nlocal),stat=alloc_stat) |
325 |
|
if (alloc_stat /= 0) then |
326 |
|
status = -1 |
328 |
– |
return |
327 |
|
end if |
328 |
|
|
329 |
|
#ifdef IS_MPI |
332 |
|
allocate(rho_tmp(nlocal),stat=alloc_stat) |
333 |
|
if (alloc_stat /= 0) then |
334 |
|
status = -1 |
337 |
– |
return |
335 |
|
end if |
336 |
|
|
337 |
|
|
339 |
|
allocate(frho_row(nAtomsInRow),stat=alloc_stat) |
340 |
|
if (alloc_stat /= 0) then |
341 |
|
status = -1 |
345 |
– |
return |
342 |
|
end if |
343 |
|
if (allocated(rho_row)) deallocate(rho_row) |
344 |
|
allocate(rho_row(nAtomsInRow),stat=alloc_stat) |
345 |
|
if (alloc_stat /= 0) then |
346 |
|
status = -1 |
351 |
– |
return |
347 |
|
end if |
348 |
|
if (allocated(dfrhodrho_row)) deallocate(dfrhodrho_row) |
349 |
|
allocate(dfrhodrho_row(nAtomsInRow),stat=alloc_stat) |
350 |
|
if (alloc_stat /= 0) then |
351 |
|
status = -1 |
357 |
– |
return |
352 |
|
end if |
353 |
|
|
354 |
|
|
358 |
|
allocate(frho_col(nAtomsInCol),stat=alloc_stat) |
359 |
|
if (alloc_stat /= 0) then |
360 |
|
status = -1 |
367 |
– |
return |
361 |
|
end if |
362 |
|
if (allocated(rho_col)) deallocate(rho_col) |
363 |
|
allocate(rho_col(nAtomsInCol),stat=alloc_stat) |
364 |
|
if (alloc_stat /= 0) then |
365 |
|
status = -1 |
373 |
– |
return |
366 |
|
end if |
367 |
|
if (allocated(dfrhodrho_col)) deallocate(dfrhodrho_col) |
368 |
|
allocate(dfrhodrho_col(nAtomsInCol),stat=alloc_stat) |
369 |
|
if (alloc_stat /= 0) then |
370 |
|
status = -1 |
379 |
– |
return |
371 |
|
end if |
372 |
|
|
373 |
|
#endif |
374 |
< |
|
374 |
> |
if (status == -1) then |
375 |
> |
call handleError("SuttonChen:allocateSC","Error in allocating SC arrays") |
376 |
> |
end if |
377 |
> |
arraysAllocated = .true. |
378 |
|
end subroutine allocateSC |
379 |
|
|
380 |
|
!! C sets rcut to be the largest cutoff of any atype |
389 |
|
|
390 |
|
end subroutine setCutoffSC |
391 |
|
|
392 |
< |
subroutine useGeometricMixing() |
399 |
< |
useGeometricDistanceMixing = .true. |
400 |
< |
haveMixingMap = .false. |
401 |
< |
return |
402 |
< |
end subroutine useGeometricMixing |
403 |
< |
|
404 |
< |
|
405 |
< |
|
406 |
< |
|
407 |
< |
|
408 |
< |
|
409 |
< |
|
410 |
< |
|
411 |
< |
|
392 |
> |
!! This array allocates module arrays if needed and builds mixing map. |
393 |
|
subroutine clean_SC() |
394 |
< |
|
394 |
> |
if (.not.arraysAllocated) call allocateSC() |
395 |
> |
if (.not.haveMixingMap) call createMixingMap() |
396 |
|
! clean non-IS_MPI first |
397 |
|
frho = 0.0_dp |
398 |
|
rho = 0.0_dp |
412 |
|
|
413 |
|
|
414 |
|
!! Calculates rho_r |
415 |
< |
subroutine calc_sc_prepair_rho(atom1,atom2,d,r,rijsq) |
415 |
> |
subroutine calc_sc_prepair_rho(atom1,atom2,d,r,rijsq, rcut) |
416 |
|
integer :: atom1,atom2 |
417 |
|
real(kind = dp), dimension(3) :: d |
418 |
|
real(kind = dp), intent(inout) :: r |
419 |
< |
real(kind = dp), intent(inout) :: rijsq |
419 |
> |
real(kind = dp), intent(inout) :: rijsq, rcut |
420 |
|
! value of electron density rho do to atom i at atom j |
421 |
|
real(kind = dp) :: rho_i_at_j |
422 |
|
! value of electron density rho do to atom j at atom i |
433 |
|
|
434 |
|
! check to see if we need to be cleaned at the start of a force loop |
435 |
|
|
436 |
< |
|
437 |
< |
|
436 |
> |
if (cleanArrays) call clean_SC() |
437 |
> |
cleanArrays = .false. |
438 |
|
|
439 |
|
#ifdef IS_MPI |
440 |
|
Atid1 = Atid_row(Atom1) |
476 |
|
integer :: myid |
477 |
|
|
478 |
|
|
497 |
– |
cleanme = .true. |
479 |
|
!! Scatter the electron density from pre-pair calculation back to local atoms |
480 |
|
#ifdef IS_MPI |
481 |
|
call scatter(rho_row,rho,plan_atom_row,sc_err) |
530 |
|
|
531 |
|
|
532 |
|
!! Does Sutton-Chen pairwise Force calculation. |
533 |
< |
subroutine do_sc_pair(atom1, atom2, d, rij, r2, sw, vpair, fpair, & |
533 |
> |
subroutine do_sc_pair(atom1, atom2, d, rij, r2, rcut, sw, vpair, fpair, & |
534 |
|
pot, f, do_pot) |
535 |
|
!Arguments |
536 |
|
integer, intent(in) :: atom1, atom2 |
537 |
< |
real( kind = dp ), intent(in) :: rij, r2 |
537 |
> |
real( kind = dp ), intent(in) :: rij, r2, rcut |
538 |
|
real( kind = dp ) :: pot, sw, vpair |
539 |
|
real( kind = dp ), dimension(3,nLocal) :: f |
540 |
|
real( kind = dp ), intent(in), dimension(3) :: d |
559 |
|
!Local Variables |
560 |
|
|
561 |
|
! write(*,*) "Frho: ", Frho(atom1) |
562 |
+ |
|
563 |
+ |
cleanArrays = .true. |
564 |
|
|
582 |
– |
|
565 |
|
dvpdr = 0.0E0_DP |
566 |
|
|
567 |
|
|