22 neighbor_kind_pairs_type
25 fist_nonbond_env_type,&
68 #include "./base/base_uses.f90"
75 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'manybody_potential'
94 particle_set, cell, pot_manybody, para_env, mm_section)
96 TYPE(fist_nonbond_env_type),
POINTER :: fist_nonbond_env
97 TYPE(atomic_kind_type),
POINTER :: atomic_kind_set(:)
98 TYPE(distribution_1d_type),
POINTER :: local_particles
99 TYPE(particle_type),
POINTER :: particle_set(:)
100 TYPE(cell_type),
POINTER :: cell
101 REAL(
dp),
INTENT(INOUT) :: pot_manybody
102 TYPE(mp_para_env_type),
OPTIONAL,
POINTER :: para_env
103 TYPE(section_vals_type),
POINTER :: mm_section
105 CHARACTER(LEN=*),
PARAMETER :: routinen =
'energy_manybody'
107 INTEGER :: atom_a, atom_b, handle, i, iend, ifirst, igrp, ikind, ilast, ilist, indexa, &
108 ipair, iparticle, iparticle_local, istart, iunique, jkind, junique, mpair, nkinds, &
109 nloc_size, npairs, nparticle, nparticle_local, nr_h3o, nr_o, nr_oh, nunique
110 INTEGER,
DIMENSION(:),
POINTER :: glob_loc_list_a, unique_list_a, work_list
111 INTEGER,
DIMENSION(:, :),
POINTER :: glob_loc_list,
list, sort_list
112 LOGICAL :: any_allegro, any_deepmd, any_gal, &
113 any_gal21, any_nequip, any_quip, &
114 any_siepmann, any_tersoff
115 REAL(kind=
dp) :: drij, embed, pot_allegro, pot_deepmd, &
116 pot_loc, pot_nequip, pot_quip, qr, &
118 REAL(kind=
dp),
DIMENSION(3) :: cell_v, cvi
119 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: glob_cell_v
120 REAL(kind=
dp),
POINTER :: fembed(:)
121 TYPE(allegro_pot_type),
POINTER :: allegro
122 TYPE(eam_pot_type),
POINTER :: eam
123 TYPE(eam_type),
DIMENSION(:),
POINTER :: eam_data
124 TYPE(fist_neighbor_type),
POINTER :: nonbonded
125 TYPE(gal21_pot_type),
POINTER :: gal21
126 TYPE(gal_pot_type),
POINTER :: gal
127 TYPE(neighbor_kind_pairs_type),
POINTER :: neighbor_kind_pair
128 TYPE(nequip_pot_type),
POINTER :: nequip
129 TYPE(pair_potential_pp_type),
POINTER :: potparm
130 TYPE(pair_potential_single_type),
POINTER :: pot
131 TYPE(pos_type),
DIMENSION(:),
POINTER :: r_last_update_pbc
132 TYPE(siepmann_pot_type),
POINTER :: siepmann
133 TYPE(tersoff_pot_type),
POINTER :: tersoff
135 NULLIFY (eam, siepmann, tersoff, gal, gal21)
136 any_tersoff = .false.
137 any_siepmann = .false.
142 any_allegro = .false.
144 CALL timeset(routinen, handle)
146 potparm=potparm, eam_data=eam_data)
148 DO ikind = 1,
SIZE(atomic_kind_set)
149 pot => potparm%pot(ikind, ikind)%pot
150 DO i = 1,
SIZE(pot%type)
151 IF (pot%type(i) /=
ea_type) cycle
152 eam => pot%set(i)%eam
153 nparticle =
SIZE(particle_set)
154 ALLOCATE (fembed(nparticle))
156 cpassert(
ASSOCIATED(eam_data))
158 nparticle_local = local_particles%n_el(ikind)
159 DO iparticle_local = 1, nparticle_local
160 iparticle = local_particles%list(ikind)%array(iparticle_local)
161 indexa = int(eam_data(iparticle)%rho/eam%drhoar) + 1
162 IF (indexa > eam%npoints - 1) indexa = eam%npoints - 1
163 qr = eam_data(iparticle)%rho - eam%rhoval(indexa)
165 embed = eam%frho(indexa) + qr*eam%frhop(indexa)
166 fembed(iparticle) = eam%frhop(indexa) + qr*(eam%frhop(indexa + 1) - eam%frhop(indexa))/eam%drhoar
168 pot_manybody = pot_manybody + embed
171 CALL para_env%sum(fembed)
172 DO iparticle = 1, nparticle
173 IF (particle_set(iparticle)%atomic_kind%kind_number == ikind)
THEN
174 eam_data(iparticle)%f_embed = fembed(iparticle)
182 DO ikind = 1,
SIZE(atomic_kind_set)
183 DO jkind = ikind,
SIZE(atomic_kind_set)
184 pot => potparm%pot(ikind, jkind)%pot
185 any_tersoff = any_tersoff .OR. any(pot%type ==
tersoff_type)
186 any_quip = any_quip .OR. any(pot%type ==
quip_type)
187 any_nequip = any_nequip .OR. any(pot%type ==
nequip_type)
188 any_allegro = any_allegro .OR. any(pot%type ==
allegro_type)
189 any_deepmd = any_deepmd .OR. any(pot%type ==
deepmd_type)
190 any_siepmann = any_siepmann .OR. any(pot%type ==
siepmann_type)
191 any_gal = any_gal .OR. any(pot%type ==
gal_type)
192 any_gal21 = any_gal21 .OR. any(pot%type ==
gal21_type)
199 fist_nonbond_env, pot_quip, para_env)
200 pot_manybody = pot_manybody + pot_quip
204 NULLIFY (glob_loc_list, glob_cell_v, glob_loc_list_a)
205 CALL setup_nequip_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
207 nequip, glob_loc_list_a, r_last_update_pbc, pot_nequip, &
208 fist_nonbond_env, para_env)
209 pot_manybody = pot_manybody + pot_nequip
213 IF (any_allegro)
THEN
214 NULLIFY (glob_loc_list, glob_cell_v, glob_loc_list_a, unique_list_a)
218 allegro, glob_loc_list_a, r_last_update_pbc, pot_allegro, &
219 fist_nonbond_env, unique_list_a)
220 pot_manybody = pot_manybody + pot_allegro
226 fist_nonbond_env, pot_deepmd, para_env)
227 pot_manybody = pot_manybody + pot_deepmd
231 IF (any_tersoff)
THEN
232 NULLIFY (glob_loc_list, glob_cell_v, glob_loc_list_a)
233 CALL setup_tersoff_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
234 DO ilist = 1, nonbonded%nlists
235 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
236 npairs = neighbor_kind_pair%npairs
237 IF (npairs == 0) cycle
238 kind_group_loop:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
239 istart = neighbor_kind_pair%grp_kind_start(igrp)
240 iend = neighbor_kind_pair%grp_kind_end(igrp)
241 ikind = neighbor_kind_pair%ij_kind(1, igrp)
242 jkind = neighbor_kind_pair%ij_kind(2, igrp)
243 list => neighbor_kind_pair%list
244 cvi = neighbor_kind_pair%cell_vector
245 pot => potparm%pot(ikind, jkind)%pot
246 DO i = 1,
SIZE(pot%type)
248 rab2_max = pot%set(i)%tersoff%rcutsq
249 cell_v = matmul(cell%hmat, cvi)
250 pot => potparm%pot(ikind, jkind)%pot
251 tersoff => pot%set(i)%tersoff
252 npairs = iend - istart + 1
253 IF (npairs /= 0)
THEN
254 ALLOCATE (sort_list(2, npairs), work_list(npairs))
255 sort_list =
list(:, istart:iend)
258 CALL sort(sort_list(1, :), npairs, work_list)
260 work_list(ipair) = sort_list(2, work_list(ipair))
262 sort_list(2, :) = work_list
265 DO ipair = 1, npairs - 1
266 IF (sort_list(1, ipair + 1) /= sort_list(1, ipair)) nunique = nunique + 1
269 junique = sort_list(1, ipair)
271 DO iunique = 1, nunique
273 IF (glob_loc_list_a(ifirst) > atom_a) cycle
274 DO mpair = ifirst,
SIZE(glob_loc_list_a)
275 IF (glob_loc_list_a(mpair) == atom_a)
EXIT
278 DO mpair = ifirst,
SIZE(glob_loc_list_a)
279 IF (glob_loc_list_a(mpair) /= atom_a)
EXIT
283 IF (ifirst /= 0) nloc_size = ilast - ifirst + 1
284 DO WHILE (ipair <= npairs)
285 IF (sort_list(1, ipair) /= junique)
EXIT
286 atom_b = sort_list(2, ipair)
289 rij(:) = r_last_update_pbc(atom_b)%r(:) - r_last_update_pbc(atom_a)%r(:) + cell_v
290 drij = dot_product(rij, rij)
292 IF (drij > rab2_max) cycle
294 CALL tersoff_energy(pot_loc, tersoff, r_last_update_pbc, atom_a, atom_b, nloc_size, &
295 glob_loc_list(:, ifirst:ilast), glob_cell_v(:, ifirst:ilast), cell_v, drij)
296 pot_manybody = pot_manybody + 0.5_dp*pot_loc
299 IF (ipair <= npairs) junique = sort_list(1, ipair)
301 DEALLOCATE (sort_list, work_list)
304 END DO kind_group_loop
310 IF (any_siepmann)
THEN
311 NULLIFY (glob_loc_list, glob_cell_v, glob_loc_list_a)
316 DO ilist = 1, nonbonded%nlists
317 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
318 npairs = neighbor_kind_pair%npairs
319 IF (npairs == 0) cycle
320 kind_group_loop_2:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
321 istart = neighbor_kind_pair%grp_kind_start(igrp)
322 iend = neighbor_kind_pair%grp_kind_end(igrp)
323 ikind = neighbor_kind_pair%ij_kind(1, igrp)
324 jkind = neighbor_kind_pair%ij_kind(2, igrp)
325 list => neighbor_kind_pair%list
326 cvi = neighbor_kind_pair%cell_vector
327 pot => potparm%pot(ikind, jkind)%pot
328 DO i = 1,
SIZE(pot%type)
330 rab2_max = pot%set(i)%siepmann%rcutsq
331 cell_v = matmul(cell%hmat, cvi)
332 pot => potparm%pot(ikind, jkind)%pot
333 siepmann => pot%set(i)%siepmann
334 npairs = iend - istart + 1
335 IF (npairs /= 0)
THEN
336 ALLOCATE (sort_list(2, npairs), work_list(npairs))
337 sort_list =
list(:, istart:iend)
340 CALL sort(sort_list(1, :), npairs, work_list)
342 work_list(ipair) = sort_list(2, work_list(ipair))
344 sort_list(2, :) = work_list
347 DO ipair = 1, npairs - 1
348 IF (sort_list(1, ipair + 1) /= sort_list(1, ipair)) nunique = nunique + 1
351 junique = sort_list(1, ipair)
353 DO iunique = 1, nunique
355 IF (glob_loc_list_a(ifirst) > atom_a) cycle
356 DO mpair = ifirst,
SIZE(glob_loc_list_a)
357 IF (glob_loc_list_a(mpair) == atom_a)
EXIT
360 DO mpair = ifirst,
SIZE(glob_loc_list_a)
361 IF (glob_loc_list_a(mpair) /= atom_a)
EXIT
365 IF (ifirst /= 0) nloc_size = ilast - ifirst + 1
366 DO WHILE (ipair <= npairs)
367 IF (sort_list(1, ipair) /= junique)
EXIT
368 atom_b = sort_list(2, ipair)
371 rij(:) = r_last_update_pbc(atom_b)%r(:) - r_last_update_pbc(atom_a)%r(:) + cell_v
372 drij = dot_product(rij, rij)
374 IF (drij > rab2_max) cycle
376 CALL siepmann_energy(pot_loc, siepmann, r_last_update_pbc, atom_a, atom_b, nloc_size, &
377 glob_loc_list(:, ifirst:ilast), cell_v, cell, drij, &
378 particle_set, nr_oh, nr_h3o, nr_o)
379 pot_manybody = pot_manybody + pot_loc
382 IF (ipair <= npairs) junique = sort_list(1, ipair)
384 DEALLOCATE (sort_list, work_list)
387 END DO kind_group_loop_2
391 print_h3o=.false., print_o=.false.)
393 print_h3o=.true., print_o=.false.)
395 print_h3o=.false., print_o=.true.)
400 NULLIFY (glob_loc_list, glob_cell_v, glob_loc_list_a)
401 CALL setup_gal_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
402 DO ilist = 1, nonbonded%nlists
403 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
404 npairs = neighbor_kind_pair%npairs
405 IF (npairs == 0) cycle
406 kind_group_loop_3:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
407 istart = neighbor_kind_pair%grp_kind_start(igrp)
408 iend = neighbor_kind_pair%grp_kind_end(igrp)
409 ikind = neighbor_kind_pair%ij_kind(1, igrp)
410 jkind = neighbor_kind_pair%ij_kind(2, igrp)
411 list => neighbor_kind_pair%list
412 cvi = neighbor_kind_pair%cell_vector
413 pot => potparm%pot(ikind, jkind)%pot
414 DO i = 1,
SIZE(pot%type)
416 rab2_max = pot%set(i)%gal%rcutsq
417 cell_v = matmul(cell%hmat, cvi)
418 pot => potparm%pot(ikind, jkind)%pot
419 gal => pot%set(i)%gal
420 npairs = iend - istart + 1
421 IF (npairs /= 0)
THEN
422 ALLOCATE (sort_list(2, npairs), work_list(npairs))
423 sort_list =
list(:, istart:iend)
426 CALL sort(sort_list(1, :), npairs, work_list)
428 work_list(ipair) = sort_list(2, work_list(ipair))
430 sort_list(2, :) = work_list
433 DO ipair = 1, npairs - 1
434 IF (sort_list(1, ipair + 1) /= sort_list(1, ipair)) nunique = nunique + 1
437 junique = sort_list(1, ipair)
439 DO iunique = 1, nunique
441 IF (glob_loc_list_a(ifirst) > atom_a) cycle
442 DO mpair = ifirst,
SIZE(glob_loc_list_a)
443 IF (glob_loc_list_a(mpair) == atom_a)
EXIT
446 DO mpair = ifirst,
SIZE(glob_loc_list_a)
447 IF (glob_loc_list_a(mpair) /= atom_a)
EXIT
451 IF (ifirst /= 0) nloc_size = ilast - ifirst + 1
452 DO WHILE (ipair <= npairs)
453 IF (sort_list(1, ipair) /= junique)
EXIT
454 atom_b = sort_list(2, ipair)
457 rij(:) = r_last_update_pbc(atom_b)%r(:) - r_last_update_pbc(atom_a)%r(:) + cell_v
458 drij = dot_product(rij, rij)
460 IF (drij > rab2_max) cycle
462 CALL gal_energy(pot_loc, gal, r_last_update_pbc, atom_a, atom_b, &
463 cell, particle_set, mm_section)
465 pot_manybody = pot_manybody + pot_loc
468 IF (ipair <= npairs) junique = sort_list(1, ipair)
470 DEALLOCATE (sort_list, work_list)
473 END DO kind_group_loop_3
480 NULLIFY (glob_loc_list, glob_cell_v, glob_loc_list_a)
481 CALL setup_gal21_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
482 DO ilist = 1, nonbonded%nlists
483 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
484 npairs = neighbor_kind_pair%npairs
485 IF (npairs == 0) cycle
486 kind_group_loop_5:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
487 istart = neighbor_kind_pair%grp_kind_start(igrp)
488 iend = neighbor_kind_pair%grp_kind_end(igrp)
489 ikind = neighbor_kind_pair%ij_kind(1, igrp)
490 jkind = neighbor_kind_pair%ij_kind(2, igrp)
491 list => neighbor_kind_pair%list
492 cvi = neighbor_kind_pair%cell_vector
493 pot => potparm%pot(ikind, jkind)%pot
494 DO i = 1,
SIZE(pot%type)
496 rab2_max = pot%set(i)%gal21%rcutsq
497 cell_v = matmul(cell%hmat, cvi)
498 pot => potparm%pot(ikind, jkind)%pot
499 gal21 => pot%set(i)%gal21
500 npairs = iend - istart + 1
501 IF (npairs /= 0)
THEN
502 ALLOCATE (sort_list(2, npairs), work_list(npairs))
503 sort_list =
list(:, istart:iend)
506 CALL sort(sort_list(1, :), npairs, work_list)
508 work_list(ipair) = sort_list(2, work_list(ipair))
510 sort_list(2, :) = work_list
513 DO ipair = 1, npairs - 1
514 IF (sort_list(1, ipair + 1) /= sort_list(1, ipair)) nunique = nunique + 1
517 junique = sort_list(1, ipair)
519 DO iunique = 1, nunique
521 IF (glob_loc_list_a(ifirst) > atom_a) cycle
522 DO mpair = ifirst,
SIZE(glob_loc_list_a)
523 IF (glob_loc_list_a(mpair) == atom_a)
EXIT
526 DO mpair = ifirst,
SIZE(glob_loc_list_a)
527 IF (glob_loc_list_a(mpair) /= atom_a)
EXIT
531 IF (ifirst /= 0) nloc_size = ilast - ifirst + 1
532 DO WHILE (ipair <= npairs)
533 IF (sort_list(1, ipair) /= junique)
EXIT
534 atom_b = sort_list(2, ipair)
537 rij(:) = r_last_update_pbc(atom_b)%r(:) - r_last_update_pbc(atom_a)%r(:) + cell_v
538 drij = dot_product(rij, rij)
540 IF (drij > rab2_max) cycle
542 CALL gal21_energy(pot_loc, gal21, r_last_update_pbc, atom_a, atom_b, &
543 cell, particle_set, mm_section)
545 pot_manybody = pot_manybody + pot_loc
548 IF (ipair <= npairs) junique = sort_list(1, ipair)
550 DEALLOCATE (sort_list, work_list)
553 END DO kind_group_loop_5
558 CALL timestop(handle)
574 f_nonbond, pv_nonbond, use_virial)
576 TYPE(fist_nonbond_env_type),
POINTER :: fist_nonbond_env
577 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
578 TYPE(cell_type),
POINTER :: cell
579 REAL(kind=
dp),
DIMENSION(:, :),
INTENT(INOUT) :: f_nonbond, pv_nonbond
580 LOGICAL,
INTENT(IN) :: use_virial
582 CHARACTER(LEN=*),
PARAMETER :: routinen =
'force_nonbond_manybody'
584 INTEGER :: atom_a, atom_b, handle, i, i_a, i_b, iend, ifirst, igrp, ikind, ilast, ilist, &
585 ipair, istart, iunique, jkind, junique, kind_a, kind_b, mpair, nkinds, nloc_size, npairs, &
587 INTEGER,
ALLOCATABLE,
DIMENSION(:, :) :: eam_kinds_index
588 INTEGER,
DIMENSION(:),
POINTER :: glob_loc_list_a, work_list
589 INTEGER,
DIMENSION(:, :),
POINTER :: glob_loc_list,
list, sort_list
590 LOGICAL :: any_allegro, any_deepmd, any_gal, &
591 any_gal21, any_nequip, any_quip, &
592 any_siepmann, any_tersoff
593 REAL(kind=
dp) :: f_eam,
fac, fr(3), ptens11, ptens12, ptens13, ptens21, ptens22, ptens23, &
594 ptens31, ptens32, ptens33, rab(3), rab2, rab2_max, rtmp(3)
595 REAL(kind=
dp),
DIMENSION(3) :: cell_v, cvi
596 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: glob_cell_v
597 TYPE(eam_pot_type),
POINTER :: eam_a, eam_b
598 TYPE(eam_type),
DIMENSION(:),
POINTER :: eam_data
599 TYPE(fist_neighbor_type),
POINTER :: nonbonded
600 TYPE(gal21_pot_type),
POINTER :: gal21
601 TYPE(gal_pot_type),
POINTER :: gal
602 TYPE(neighbor_kind_pairs_type),
POINTER :: neighbor_kind_pair
603 TYPE(pair_potential_pp_type),
POINTER :: potparm
604 TYPE(pair_potential_single_type),
POINTER :: pot
605 TYPE(pos_type),
DIMENSION(:),
POINTER :: r_last_update_pbc
606 TYPE(siepmann_pot_type),
POINTER :: siepmann
607 TYPE(tersoff_pot_type),
POINTER :: tersoff
609 any_tersoff = .false.
612 any_allegro = .false.
613 any_siepmann = .false.
617 CALL timeset(routinen, handle)
618 NULLIFY (eam_a, eam_b, tersoff, siepmann, gal, gal21)
621 natom_types=nkinds, eam_data=eam_data, r_last_update_pbc=r_last_update_pbc)
625 ptens11 = 0.0_dp; ptens12 = 0.0_dp; ptens13 = 0.0_dp
626 ptens21 = 0.0_dp; ptens22 = 0.0_dp; ptens23 = 0.0_dp
627 ptens31 = 0.0_dp; ptens32 = 0.0_dp; ptens33 = 0.0_dp
630 nkinds =
SIZE(potparm%pot, 1)
631 ALLOCATE (eam_kinds_index(nkinds, nkinds))
634 DO jkind = ikind, nkinds
635 DO i = 1,
SIZE(potparm%pot(ikind, jkind)%pot%type)
636 IF (potparm%pot(ikind, jkind)%pot%type(i) ==
ea_type)
THEN
638 cpassert(eam_kinds_index(ikind, jkind) == -1)
639 cpassert(eam_kinds_index(jkind, ikind) == -1)
640 eam_kinds_index(ikind, jkind) = i
641 eam_kinds_index(jkind, ikind) = i
647 DO jkind = ikind, nkinds
648 any_deepmd = any_deepmd .OR. any(potparm%pot(ikind, jkind)%pot%type ==
deepmd_type)
660 DO ilist = 1, nonbonded%nlists
661 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
662 npairs = neighbor_kind_pair%npairs
663 IF (npairs == 0) cycle
664 kind_group_loop1:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
665 istart = neighbor_kind_pair%grp_kind_start(igrp)
666 iend = neighbor_kind_pair%grp_kind_end(igrp)
667 ikind = neighbor_kind_pair%ij_kind(1, igrp)
668 jkind = neighbor_kind_pair%ij_kind(2, igrp)
669 list => neighbor_kind_pair%list
670 cvi = neighbor_kind_pair%cell_vector
671 pot => potparm%pot(ikind, jkind)%pot
673 rab2_max = pot%rcutsq
674 cell_v = matmul(cell%hmat, cvi)
675 any_tersoff = any_tersoff .OR. any(pot%type ==
tersoff_type)
676 any_siepmann = any_siepmann .OR. any(pot%type ==
siepmann_type)
677 any_deepmd = any_deepmd .OR. any(pot%type ==
deepmd_type)
678 any_gal = any_gal .OR. any(pot%type ==
gal_type)
679 any_gal21 = any_gal21 .OR. any(pot%type ==
gal21_type)
680 any_nequip = any_nequip .OR. any(pot%type ==
nequip_type)
681 any_allegro = any_allegro .OR. any(pot%type ==
allegro_type)
682 i = eam_kinds_index(ikind, jkind)
685 cpassert(
ASSOCIATED(eam_data))
686 DO ipair = istart, iend
687 atom_a =
list(1, ipair)
688 atom_b =
list(2, ipair)
690 IF (atom_a == atom_b)
fac = 0.5_dp
691 kind_a = particle_set(atom_a)%atomic_kind%kind_number
692 kind_b = particle_set(atom_b)%atomic_kind%kind_number
693 i_a = eam_kinds_index(kind_a, kind_a)
694 i_b = eam_kinds_index(kind_b, kind_b)
695 eam_a => potparm%pot(kind_a, kind_a)%pot%set(i_a)%eam
696 eam_b => potparm%pot(kind_b, kind_b)%pot%set(i_b)%eam
700 rab = r_last_update_pbc(atom_b)%r - r_last_update_pbc(atom_a)%r
702 rab2 = rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3)
703 IF (rab2 <= rab2_max)
THEN
704 CALL get_force_eam(rab2, eam_a, eam_b, eam_data, atom_a, atom_b, f_eam)
707 fr(1) = -f_eam*rab(1)
708 fr(2) = -f_eam*rab(2)
709 fr(3) = -f_eam*rab(3)
710 f_nonbond(1, atom_a) = f_nonbond(1, atom_a) - fr(1)
711 f_nonbond(2, atom_a) = f_nonbond(2, atom_a) - fr(2)
712 f_nonbond(3, atom_a) = f_nonbond(3, atom_a) - fr(3)
714 f_nonbond(1, atom_b) = f_nonbond(1, atom_b) + fr(1)
715 f_nonbond(2, atom_b) = f_nonbond(2, atom_b) + fr(2)
716 f_nonbond(3, atom_b) = f_nonbond(3, atom_b) + fr(3)
718 ptens11 = ptens11 + rab(1)*fr(1)
719 ptens21 = ptens21 + rab(2)*fr(1)
720 ptens31 = ptens31 + rab(3)*fr(1)
721 ptens12 = ptens12 + rab(1)*fr(2)
722 ptens22 = ptens22 + rab(2)*fr(2)
723 ptens32 = ptens32 + rab(3)*fr(2)
724 ptens13 = ptens13 + rab(1)*fr(3)
725 ptens23 = ptens23 + rab(2)*fr(3)
726 ptens33 = ptens33 + rab(3)*fr(3)
730 END DO kind_group_loop1
732 DEALLOCATE (eam_kinds_index)
735 IF (any_tersoff)
THEN
736 NULLIFY (glob_loc_list, glob_cell_v, glob_loc_list_a)
737 CALL setup_tersoff_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
738 DO ilist = 1, nonbonded%nlists
739 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
740 npairs = neighbor_kind_pair%npairs
741 IF (npairs == 0) cycle
742 kind_group_loop2:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
743 istart = neighbor_kind_pair%grp_kind_start(igrp)
744 iend = neighbor_kind_pair%grp_kind_end(igrp)
745 ikind = neighbor_kind_pair%ij_kind(1, igrp)
746 jkind = neighbor_kind_pair%ij_kind(2, igrp)
747 list => neighbor_kind_pair%list
748 cvi = neighbor_kind_pair%cell_vector
749 pot => potparm%pot(ikind, jkind)%pot
752 rab2_max = pot%rcutsq
753 cell_v = matmul(cell%hmat, cvi)
754 DO i = 1,
SIZE(pot%type)
757 npairs = iend - istart + 1
758 tersoff => pot%set(i)%tersoff
759 ALLOCATE (sort_list(2, npairs), work_list(npairs))
760 sort_list =
list(:, istart:iend)
763 CALL sort(sort_list(1, :), npairs, work_list)
765 work_list(ipair) = sort_list(2, work_list(ipair))
767 sort_list(2, :) = work_list
770 DO ipair = 1, npairs - 1
771 IF (sort_list(1, ipair + 1) /= sort_list(1, ipair)) nunique = nunique + 1
774 junique = sort_list(1, ipair)
776 DO iunique = 1, nunique
778 IF (glob_loc_list_a(ifirst) > atom_a) cycle
779 DO mpair = ifirst,
SIZE(glob_loc_list_a)
780 IF (glob_loc_list_a(mpair) == atom_a)
EXIT
783 DO mpair = ifirst,
SIZE(glob_loc_list_a)
784 IF (glob_loc_list_a(mpair) /= atom_a)
EXIT
788 IF (ifirst /= 0) nloc_size = ilast - ifirst + 1
789 DO WHILE (ipair <= npairs)
790 IF (sort_list(1, ipair) /= junique)
EXIT
791 atom_b = sort_list(2, ipair)
793 rtmp = r_last_update_pbc(atom_b)%r(:) - r_last_update_pbc(atom_a)%r(:) + cell_v
795 IF (dot_product(rtmp, rtmp) <= tersoff%rcutsq)
THEN
797 nloc_size, glob_loc_list(:, ifirst:ilast), glob_cell_v(:, ifirst:ilast), &
798 atom_a, atom_b, f_nonbond, pv_nonbond, use_virial, tersoff%rcutsq)
802 IF (ipair <= npairs) junique = sort_list(1, ipair)
804 DEALLOCATE (sort_list, work_list)
807 END DO kind_group_loop2
812 IF (any_siepmann)
THEN
813 NULLIFY (glob_loc_list, glob_cell_v, glob_loc_list_a)
815 DO ilist = 1, nonbonded%nlists
816 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
817 npairs = neighbor_kind_pair%npairs
818 IF (npairs == 0) cycle
819 kind_group_loop3:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
820 istart = neighbor_kind_pair%grp_kind_start(igrp)
821 iend = neighbor_kind_pair%grp_kind_end(igrp)
822 ikind = neighbor_kind_pair%ij_kind(1, igrp)
823 jkind = neighbor_kind_pair%ij_kind(2, igrp)
824 list => neighbor_kind_pair%list
825 cvi = neighbor_kind_pair%cell_vector
826 pot => potparm%pot(ikind, jkind)%pot
829 rab2_max = pot%rcutsq
830 cell_v = matmul(cell%hmat, cvi)
831 DO i = 1,
SIZE(pot%type)
834 npairs = iend - istart + 1
835 siepmann => pot%set(i)%siepmann
836 ALLOCATE (sort_list(2, npairs), work_list(npairs))
837 sort_list =
list(:, istart:iend)
840 CALL sort(sort_list(1, :), npairs, work_list)
842 work_list(ipair) = sort_list(2, work_list(ipair))
844 sort_list(2, :) = work_list
847 DO ipair = 1, npairs - 1
848 IF (sort_list(1, ipair + 1) /= sort_list(1, ipair)) nunique = nunique + 1
851 junique = sort_list(1, ipair)
853 DO iunique = 1, nunique
855 IF (glob_loc_list_a(ifirst) > atom_a) cycle
856 DO mpair = ifirst,
SIZE(glob_loc_list_a)
857 IF (glob_loc_list_a(mpair) == atom_a)
EXIT
860 DO mpair = ifirst,
SIZE(glob_loc_list_a)
861 IF (glob_loc_list_a(mpair) /= atom_a)
EXIT
865 IF (ifirst /= 0) nloc_size = ilast - ifirst + 1
866 DO WHILE (ipair <= npairs)
867 IF (sort_list(1, ipair) /= junique)
EXIT
868 atom_b = sort_list(2, ipair)
870 rtmp = r_last_update_pbc(atom_b)%r(:) - r_last_update_pbc(atom_a)%r(:) + cell_v
872 IF (dot_product(rtmp, rtmp) <= siepmann%rcutsq)
THEN
874 atom_a, atom_b, f_nonbond, use_virial, siepmann%rcutsq, &
877 nloc_size, glob_loc_list(:, ifirst:ilast), &
878 atom_a, atom_b, f_nonbond, use_virial, siepmann%rcutsq, &
883 IF (ipair <= npairs) junique = sort_list(1, ipair)
885 DEALLOCATE (sort_list, work_list)
888 END DO kind_group_loop3
895 NULLIFY (glob_loc_list, glob_cell_v, glob_loc_list_a)
896 CALL setup_gal_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
897 DO ilist = 1, nonbonded%nlists
898 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
899 npairs = neighbor_kind_pair%npairs
900 IF (npairs == 0) cycle
901 kind_group_loop4:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
902 istart = neighbor_kind_pair%grp_kind_start(igrp)
903 iend = neighbor_kind_pair%grp_kind_end(igrp)
904 ikind = neighbor_kind_pair%ij_kind(1, igrp)
905 jkind = neighbor_kind_pair%ij_kind(2, igrp)
906 list => neighbor_kind_pair%list
907 cvi = neighbor_kind_pair%cell_vector
908 pot => potparm%pot(ikind, jkind)%pot
911 rab2_max = pot%rcutsq
912 cell_v = matmul(cell%hmat, cvi)
913 DO i = 1,
SIZE(pot%type)
916 npairs = iend - istart + 1
917 gal => pot%set(i)%gal
918 ALLOCATE (sort_list(2, npairs), work_list(npairs))
919 sort_list =
list(:, istart:iend)
922 CALL sort(sort_list(1, :), npairs, work_list)
924 work_list(ipair) = sort_list(2, work_list(ipair))
926 sort_list(2, :) = work_list
929 DO ipair = 1, npairs - 1
930 IF (sort_list(1, ipair + 1) /= sort_list(1, ipair)) nunique = nunique + 1
933 junique = sort_list(1, ipair)
935 DO iunique = 1, nunique
937 IF (glob_loc_list_a(ifirst) > atom_a) cycle
938 DO mpair = ifirst,
SIZE(glob_loc_list_a)
939 IF (glob_loc_list_a(mpair) == atom_a)
EXIT
942 DO mpair = ifirst,
SIZE(glob_loc_list_a)
943 IF (glob_loc_list_a(mpair) /= atom_a)
EXIT
947 IF (ifirst /= 0) nloc_size = ilast - ifirst + 1
948 DO WHILE (ipair <= npairs)
949 IF (sort_list(1, ipair) /= junique)
EXIT
950 atom_b = sort_list(2, ipair)
952 rtmp = r_last_update_pbc(atom_b)%r(:) - r_last_update_pbc(atom_a)%r(:) + cell_v
954 IF (dot_product(rtmp, rtmp) <= gal%rcutsq)
THEN
956 atom_a, atom_b, f_nonbond, use_virial, &
961 IF (ipair <= npairs) junique = sort_list(1, ipair)
963 DEALLOCATE (sort_list, work_list)
966 END DO kind_group_loop4
973 NULLIFY (glob_loc_list, glob_cell_v, glob_loc_list_a)
974 CALL setup_gal21_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
975 DO ilist = 1, nonbonded%nlists
976 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
977 npairs = neighbor_kind_pair%npairs
978 IF (npairs == 0) cycle
979 kind_group_loop6:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
980 istart = neighbor_kind_pair%grp_kind_start(igrp)
981 iend = neighbor_kind_pair%grp_kind_end(igrp)
982 ikind = neighbor_kind_pair%ij_kind(1, igrp)
983 jkind = neighbor_kind_pair%ij_kind(2, igrp)
984 list => neighbor_kind_pair%list
985 cvi = neighbor_kind_pair%cell_vector
986 pot => potparm%pot(ikind, jkind)%pot
989 rab2_max = pot%rcutsq
990 cell_v = matmul(cell%hmat, cvi)
991 DO i = 1,
SIZE(pot%type)
994 npairs = iend - istart + 1
995 gal21 => pot%set(i)%gal21
996 ALLOCATE (sort_list(2, npairs), work_list(npairs))
997 sort_list =
list(:, istart:iend)
1000 CALL sort(sort_list(1, :), npairs, work_list)
1001 DO ipair = 1, npairs
1002 work_list(ipair) = sort_list(2, work_list(ipair))
1004 sort_list(2, :) = work_list
1007 DO ipair = 1, npairs - 1
1008 IF (sort_list(1, ipair + 1) /= sort_list(1, ipair)) nunique = nunique + 1
1011 junique = sort_list(1, ipair)
1013 DO iunique = 1, nunique
1015 IF (glob_loc_list_a(ifirst) > atom_a) cycle
1016 DO mpair = ifirst,
SIZE(glob_loc_list_a)
1017 IF (glob_loc_list_a(mpair) == atom_a)
EXIT
1020 DO mpair = ifirst,
SIZE(glob_loc_list_a)
1021 IF (glob_loc_list_a(mpair) /= atom_a)
EXIT
1025 IF (ifirst /= 0) nloc_size = ilast - ifirst + 1
1026 DO WHILE (ipair <= npairs)
1027 IF (sort_list(1, ipair) /= junique)
EXIT
1028 atom_b = sort_list(2, ipair)
1030 rtmp = r_last_update_pbc(atom_b)%r(:) - r_last_update_pbc(atom_a)%r(:) + cell_v
1032 IF (dot_product(rtmp, rtmp) <= gal21%rcutsq)
THEN
1034 atom_a, atom_b, f_nonbond, pv_nonbond, use_virial, &
1039 IF (ipair <= npairs) junique = sort_list(1, ipair)
1041 DEALLOCATE (sort_list, work_list)
1044 END DO kind_group_loop6
1050 IF (any_nequip)
THEN
1055 IF (any_allegro)
THEN
1059 IF (use_virial)
THEN
1060 pv_nonbond(1, 1) = pv_nonbond(1, 1) + ptens11
1061 pv_nonbond(1, 2) = pv_nonbond(1, 2) + ptens12
1062 pv_nonbond(1, 3) = pv_nonbond(1, 3) + ptens13
1063 pv_nonbond(2, 1) = pv_nonbond(2, 1) + ptens21
1064 pv_nonbond(2, 2) = pv_nonbond(2, 2) + ptens22
1065 pv_nonbond(2, 3) = pv_nonbond(2, 3) + ptens23
1066 pv_nonbond(3, 1) = pv_nonbond(3, 1) + ptens31
1067 pv_nonbond(3, 2) = pv_nonbond(3, 2) + ptens32
1068 pv_nonbond(3, 3) = pv_nonbond(3, 3) + ptens33
1070 CALL timestop(handle)
static GRID_HOST_DEVICE double fac(const int i)
Factorial function, e.g. fac(5) = 5! = 120.
Define the atomic kind types and their sub types.
Handles all functions related to the CELL.
stores a lists of integer that are local to a processor. The idea is that these integers represent ob...
Define the neighbor list data types and the corresponding functionality.
subroutine, public fist_nonbond_env_get(fist_nonbond_env, potparm14, potparm, nonbonded, rlist_cut, rlist_lowsq, aup, lup, ei_scale14, vdw_scale14, shift_cutoff, do_electrostatics, r_last_update, r_last_update_pbc, rshell_last_update_pbc, rcore_last_update_pbc, cell_last_update, num_update, last_update, counter, natom_types, long_range_correction, ij_kind_full_fac, eam_data, quip_data, nequip_data, allegro_data, deepmd_data, charges)
sets a fist_nonbond_env
Defines the basic variable types.
integer, parameter, public dp
An array-based list which grows on demand. When the internal array is full, a new array of twice the ...
subroutine, public setup_allegro_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, unique_list_a, cell)
...
subroutine, public allegro_add_force_virial(fist_nonbond_env, f_nonbond, pv_nonbond, use_virial)
...
subroutine, public destroy_allegro_arrays(glob_loc_list, glob_cell_v, glob_loc_list_a, unique_list_a)
...
subroutine, public allegro_energy_store_force_virial(nonbonded, particle_set, cell, atomic_kind_set, potparm, allegro, glob_loc_list_a, r_last_update_pbc, pot_allegro, fist_nonbond_env, unique_list_a)
...
subroutine, public deepmd_energy_store_force_virial(particle_set, cell, atomic_kind_set, potparm, fist_nonbond_env, pot_deepmd, para_env)
...
subroutine, public deepmd_add_force_virial(fist_nonbond_env, force, pv_nonbond, use_virial)
...
subroutine, public get_force_eam(rab2, eam_a, eam_b, eam_data, atom_a, atom_b, f_eam)
...
Implementation of the GAL21 potential.
subroutine, public destroy_gal21_arrays(glob_loc_list, glob_cell_v, glob_loc_list_a)
...
subroutine, public gal21_forces(gal21, r_last_update_pbc, iparticle, jparticle, f_nonbond, pv_nonbond, use_virial, cell, particle_set)
forces generated by the GAL2119 potential
subroutine, public setup_gal21_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
...
subroutine, public gal21_energy(pot_loc, gal21, r_last_update_pbc, iparticle, jparticle, cell, particle_set, mm_section)
Main part of the energy evaluation of GAL2119.
Implementation of the GAL19 potential.
subroutine, public destroy_gal_arrays(glob_loc_list, glob_cell_v, glob_loc_list_a)
...
subroutine, public setup_gal_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
...
subroutine, public gal_energy(pot_loc, gal, r_last_update_pbc, iparticle, jparticle, cell, particle_set, mm_section)
Main part of the energy evaluation of GAL19.
subroutine, public gal_forces(gal, r_last_update_pbc, iparticle, jparticle, f_nonbond, use_virial, cell, particle_set)
forces generated by the GAL19 potential
subroutine, public nequip_add_force_virial(fist_nonbond_env, f_nonbond, pv_nonbond, use_virial)
...
subroutine, public nequip_energy_store_force_virial(nonbonded, particle_set, cell, atomic_kind_set, potparm, nequip, glob_loc_list_a, r_last_update_pbc, pot_nequip, fist_nonbond_env, para_env)
...
subroutine, public destroy_nequip_arrays(glob_loc_list, glob_cell_v, glob_loc_list_a)
...
subroutine, public setup_nequip_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
...
subroutine, public energy_manybody(fist_nonbond_env, atomic_kind_set, local_particles, particle_set, cell, pot_manybody, para_env, mm_section)
computes the embedding contribution to the energy
subroutine, public force_nonbond_manybody(fist_nonbond_env, particle_set, cell, f_nonbond, pv_nonbond, use_virial)
...
subroutine, public quip_add_force_virial(fist_nonbond_env, force, virial)
...
subroutine, public quip_energy_store_force_virial(particle_set, cell, atomic_kind_set, potparm, fist_nonbond_env, pot_quip, para_env)
...
implementation of dipole and three-body part of Siepmann-Sprik potential dipole term: 3rd term in Eq....
subroutine, public siepmann_energy(pot_loc, siepmann, r_last_update_pbc, atom_a, atom_b, nloc_size, full_loc_list, cell_v, cell, drij, particle_set, nr_oh, nr_h3o, nr_o)
energy of two-body dipole term and three-body term
subroutine, public setup_siepmann_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
...
subroutine, public siepmann_forces_v2(siepmann, r_last_update_pbc, cell_v, cell, iparticle, jparticle, f_nonbond, use_virial, rcutsq, particle_set)
forces generated by the dipole term
subroutine, public siepmann_forces_v3(siepmann, r_last_update_pbc, cell_v, n_loc_size, full_loc_list, iparticle, jparticle, f_nonbond, use_virial, rcutsq, cell, particle_set)
forces generated by the three-body term
subroutine, public print_nr_ions_siepmann(nr_ions, mm_section, para_env, print_oh, print_h3o, print_o)
prints the number of OH- ions or H3O+ ions near surface
subroutine, public destroy_siepmann_arrays(glob_loc_list, glob_cell_v, glob_loc_list_a)
...
subroutine, public setup_tersoff_arrays(nonbonded, potparm, glob_loc_list, glob_cell_v, glob_loc_list_a, cell)
...
subroutine, public tersoff_forces(tersoff, r_last_update_pbc, cell_v, n_loc_size, full_loc_list, loc_cell_v, iparticle, jparticle, f_nonbond, pv_nonbond, use_virial, rcutsq)
...
subroutine, public destroy_tersoff_arrays(glob_loc_list, glob_cell_v, glob_loc_list_a)
...
subroutine, public tersoff_energy(pot_loc, tersoff, r_last_update_pbc, atom_a, atom_b, nloc_size, full_loc_list, loc_cell_v, cell_v, drij)
...
Interface to the message passing library MPI.
integer, parameter, public allegro_type
integer, parameter, public gal_type
integer, parameter, public nequip_type
integer, parameter, public deepmd_type
integer, parameter, public quip_type
integer, parameter, public siepmann_type
integer, parameter, public gal21_type
integer, parameter, public ea_type
integer, parameter, public tersoff_type
Define the data structure for the particle information.
All kind of helpful little routines.