17 neighbor_kind_pairs_type
22 pair_potential_single_type,&
26 #include "./base/base_uses.f90"
33 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'manybody_tersoff'
51 SUBROUTINE tersoff_energy(pot_loc, tersoff, r_last_update_pbc, atom_a, atom_b, nloc_size, &
52 full_loc_list, loc_cell_v, cell_v, drij)
54 REAL(kind=
dp),
INTENT(OUT) :: pot_loc
55 TYPE(tersoff_pot_type),
POINTER :: tersoff
56 TYPE(pos_type),
DIMENSION(:),
POINTER :: r_last_update_pbc
57 INTEGER,
INTENT(IN) :: atom_a, atom_b, nloc_size
58 INTEGER,
DIMENSION(2, 1:nloc_size) :: full_loc_list
59 REAL(kind=
dp),
DIMENSION(3, 1:nloc_size) :: loc_cell_v
60 REAL(kind=
dp),
DIMENSION(3) :: cell_v
63 REAL(kind=
dp) :: b_ij, f_a, f_c, f_r
65 b_ij = ter_b_ij(tersoff, r_last_update_pbc, atom_a, atom_b, nloc_size, &
66 full_loc_list, loc_cell_v, cell_v, tersoff%rcutsq)
67 f_c = ter_f_c(tersoff, drij)
68 f_a = ter_f_a(tersoff, drij)
69 f_r = ter_f_r(tersoff, drij)
70 pot_loc = f_c*(f_r + b_ij*f_a)
81 FUNCTION ter_f_c(tersoff, r)
82 TYPE(tersoff_pot_type),
POINTER :: tersoff
83 REAL(kind=
dp),
INTENT(IN) :: r
84 REAL(kind=
dp) :: ter_f_c
86 REAL(kind=
dp) :: bigd, bigr, rmd, rpd
90 rmd = tersoff%bigR - tersoff%bigD
91 rpd = tersoff%bigR + tersoff%bigD
93 IF (r < rmd) ter_f_c = 1.0_dp
94 IF (r > rpd) ter_f_c = 0.0_dp
95 IF ((r < rpd) .AND. (r > rmd))
THEN
96 ter_f_c = 0.5_dp*(1.0_dp - sin(0.5_dp*
pi*(r - bigr)/(bigd)))
107 FUNCTION ter_f_c_d(tersoff, r)
108 TYPE(tersoff_pot_type),
POINTER :: tersoff
109 REAL(kind=
dp),
INTENT(IN) :: r
110 REAL(kind=
dp) :: ter_f_c_d
112 REAL(kind=
dp) :: bigd, bigr, rmd, rpd
116 rmd = tersoff%bigR - tersoff%bigD
117 rpd = tersoff%bigR + tersoff%bigD
119 IF (r < rmd) ter_f_c_d = 0.0_dp
120 IF (r > rpd) ter_f_c_d = 0.0_dp
121 IF ((r < rpd) .AND. (r > rmd))
THEN
122 ter_f_c_d = (0.25_dp*
pi/bigd)*cos(0.5_dp*
pi*(r - bigr)/(bigd))/r
125 END FUNCTION ter_f_c_d
134 FUNCTION ter_f_r(tersoff, r)
135 TYPE(tersoff_pot_type),
POINTER :: tersoff
136 REAL(kind=
dp),
INTENT(IN) :: r
137 REAL(kind=
dp) :: ter_f_r
139 REAL(kind=
dp) :: a, lambda1
142 lambda1 = tersoff%lambda1
144 ter_f_r = a*exp(-lambda1*r)
155 FUNCTION ter_f_r_d(tersoff, r)
156 TYPE(tersoff_pot_type),
POINTER :: tersoff
157 REAL(kind=
dp),
INTENT(IN) :: r
158 REAL(kind=
dp) :: ter_f_r_d
160 REAL(kind=
dp) :: a, f_r, lambda1
163 lambda1 = tersoff%lambda1
164 f_r = a*exp(-lambda1*r)
166 ter_f_r_d = lambda1*f_r/r
168 END FUNCTION ter_f_r_d
177 FUNCTION ter_f_a(tersoff, r)
178 TYPE(tersoff_pot_type),
POINTER :: tersoff
179 REAL(kind=
dp),
INTENT(IN) :: r
180 REAL(kind=
dp) :: ter_f_a
182 REAL(kind=
dp) :: b, lambda2
185 lambda2 = tersoff%lambda2
187 ter_f_a = -b*exp(-lambda2*r)
198 FUNCTION ter_f_a_d(tersoff, r)
199 TYPE(tersoff_pot_type),
POINTER :: tersoff
200 REAL(kind=
dp),
INTENT(IN) :: r
201 REAL(kind=
dp) :: ter_f_a_d
203 REAL(kind=
dp) :: b, lambda2
206 lambda2 = tersoff%lambda2
208 ter_f_a_d = -b*lambda2*exp(-lambda2*r)/r
210 END FUNCTION ter_f_a_d
218 FUNCTION ter_a_ij(tersoff)
219 TYPE(tersoff_pot_type),
POINTER :: tersoff
220 REAL(kind=
dp) :: ter_a_ij
222 REAL(kind=
dp) :: alpha, n
225 alpha = tersoff%alpha
231 END FUNCTION ter_a_ij
247 FUNCTION ter_b_ij(tersoff, r_last_update_pbc, iparticle, jparticle, n_loc_size, &
248 full_loc_list, loc_cell_v, cell_v, rcutsq)
249 TYPE(tersoff_pot_type),
POINTER :: tersoff
250 TYPE(pos_type),
DIMENSION(:),
POINTER :: r_last_update_pbc
251 INTEGER,
INTENT(IN) :: iparticle, jparticle, n_loc_size
252 INTEGER,
DIMENSION(2, 1:n_loc_size) :: full_loc_list
253 REAL(kind=
dp),
DIMENSION(3, 1:n_loc_size) :: loc_cell_v
254 REAL(kind=
dp),
DIMENSION(3) :: cell_v
255 REAL(kind=
dp),
INTENT(IN) :: rcutsq
256 REAL(kind=
dp) :: ter_b_ij
258 REAL(kind=
dp) :: beta, n, zeta_ij
263 zeta_ij = ter_zeta_ij(tersoff, r_last_update_pbc, iparticle, jparticle, &
264 n_loc_size, full_loc_list, loc_cell_v, cell_v, rcutsq)
265 ter_b_ij = (1.0_dp + (beta*zeta_ij)**n)**(-0.5_dp/n)
267 END FUNCTION ter_b_ij
283 FUNCTION ter_b_ij_d(tersoff, r_last_update_pbc, iparticle, jparticle, n_loc_size, &
284 full_loc_list, loc_cell_v, cell_v, rcutsq)
285 TYPE(tersoff_pot_type),
POINTER :: tersoff
286 TYPE(pos_type),
DIMENSION(:),
POINTER :: r_last_update_pbc
287 INTEGER,
INTENT(IN) :: iparticle, jparticle, n_loc_size
288 INTEGER,
DIMENSION(2, 1:n_loc_size) :: full_loc_list
289 REAL(kind=
dp),
DIMENSION(3, 1:n_loc_size) :: loc_cell_v
290 REAL(kind=
dp),
DIMENSION(3) :: cell_v
291 REAL(kind=
dp),
INTENT(IN) :: rcutsq
292 REAL(kind=
dp) :: ter_b_ij_d
294 REAL(kind=
dp) :: beta, beta_n, n, zeta_ij, zeta_ij_n, &
300 zeta_ij = ter_zeta_ij(tersoff, r_last_update_pbc, iparticle, jparticle, n_loc_size, &
301 full_loc_list, loc_cell_v, cell_v, rcutsq)
303 IF (zeta_ij > 0.0_dp) zeta_ij_nm1 = zeta_ij**(n - 1.0_dp)
304 zeta_ij_n = zeta_ij**(n)
307 ter_b_ij_d = -0.5_dp*beta_n*zeta_ij_nm1* &
308 ((1.0_dp + beta_n*zeta_ij_n)**((-0.5_dp/n) - 1.0_dp))
310 END FUNCTION ter_b_ij_d
328 FUNCTION ter_zeta_ij(tersoff, r_last_update_pbc, iparticle, jparticle, n_loc_size, &
329 full_loc_list, loc_cell_v, cell_v, rcutsq)
330 TYPE(tersoff_pot_type),
POINTER :: tersoff
331 TYPE(pos_type),
DIMENSION(:),
POINTER :: r_last_update_pbc
332 INTEGER,
INTENT(IN) :: iparticle, jparticle, n_loc_size
333 INTEGER,
DIMENSION(2, 1:n_loc_size) :: full_loc_list
334 REAL(kind=
dp),
DIMENSION(3, 1:n_loc_size) :: loc_cell_v
335 REAL(kind=
dp),
DIMENSION(3) :: cell_v
336 REAL(kind=
dp),
INTENT(IN) :: rcutsq
337 REAL(kind=
dp) :: ter_zeta_ij
339 INTEGER :: ilist, kparticle
340 REAL(kind=
dp) :: cell_v_2(3), costheta, drij, drik, &
341 expterm, f_c, gterm, lambda3, n, &
342 rab2_max, rij(3), rik(3)
346 lambda3 = tersoff%lambda3
348 rij(:) = r_last_update_pbc(jparticle)%r(:) - r_last_update_pbc(iparticle)%r(:) + cell_v
349 drij = sqrt(dot_product(rij, rij))
351 DO ilist = 1, n_loc_size
352 kparticle = full_loc_list(2, ilist)
353 IF (kparticle == jparticle) cycle
354 cell_v_2 = loc_cell_v(:, ilist)
355 rik(:) = r_last_update_pbc(kparticle)%r(:) - r_last_update_pbc(iparticle)%r(:) + cell_v_2
356 drik = dot_product(rik, rik)
357 IF (drik > rab2_max) cycle
359 costheta = dot_product(rij, rik)/(drij*drik)
360 IF (costheta < -1.0_dp) costheta = -1.0_dp
361 IF (costheta > +1.0_dp) costheta = +1.0_dp
362 f_c = ter_f_c(tersoff, drik)
363 gterm = ter_g(tersoff, costheta)
364 expterm = exp((lambda3*(drij - drik))**3)
365 ter_zeta_ij = ter_zeta_ij + f_c*gterm*expterm
368 END FUNCTION ter_zeta_ij
389 SUBROUTINE ter_zeta_ij_d(tersoff, r_last_update_pbc, iparticle, jparticle, f_nonbond, pv_nonbond, prefactor, &
390 n_loc_size, full_loc_list, loc_cell_v, cell_v, rcutsq, use_virial)
391 TYPE(tersoff_pot_type),
POINTER :: tersoff
392 TYPE(pos_type),
DIMENSION(:),
POINTER :: r_last_update_pbc
393 INTEGER,
INTENT(IN) :: iparticle, jparticle
394 REAL(kind=
dp),
DIMENSION(:, :),
INTENT(INOUT) :: f_nonbond, pv_nonbond
395 REAL(kind=
dp),
INTENT(IN) :: prefactor
396 INTEGER,
INTENT(IN) :: n_loc_size
397 INTEGER,
DIMENSION(2, 1:n_loc_size) :: full_loc_list
398 REAL(kind=
dp),
DIMENSION(3, 1:n_loc_size) :: loc_cell_v
399 REAL(kind=
dp),
DIMENSION(3) :: cell_v
400 REAL(kind=
dp),
INTENT(IN) :: rcutsq
401 LOGICAL,
INTENT(IN) :: use_virial
403 INTEGER :: ilist, kparticle, nparticle
404 REAL(kind=
dp) :: costheta, drij, drik, expterm, &
405 expterm_d, f_c, f_c_d, gterm, gterm_d, &
407 REAL(kind=
dp),
DIMENSION(3) :: cell_v_2, dcosdri, dcosdrj, dcosdrk, &
408 dri, drj, drk, rij, rij_hat, rik, &
412 lambda3 = tersoff%lambda3
415 rij(:) = r_last_update_pbc(jparticle)%r(:) - r_last_update_pbc(iparticle)%r(:) + cell_v
416 drij = sqrt(dot_product(rij, rij))
417 rij_hat(:) = rij(:)/drij
419 nparticle =
SIZE(r_last_update_pbc)
420 DO ilist = 1, n_loc_size
421 kparticle = full_loc_list(2, ilist)
422 IF (kparticle == jparticle) cycle
423 cell_v_2 = loc_cell_v(:, ilist)
424 rik(:) = r_last_update_pbc(kparticle)%r(:) - r_last_update_pbc(iparticle)%r(:) + cell_v_2
425 drik = dot_product(rik, rik)
427 IF (drik > rab2_max) cycle
429 rik_hat(:) = rik(:)/drik
430 costheta = dot_product(rij, rik)/(drij*drik)
431 IF (costheta < -1.0_dp) costheta = -1.0_dp
432 IF (costheta > +1.0_dp) costheta = +1.0_dp
434 dcosdrj(:) = (1.0_dp/(drij))*(rik_hat(:) - costheta*rij_hat(:))
435 dcosdrk(:) = (1.0_dp/(drik))*(rij_hat(:) - costheta*rik_hat(:))
436 dcosdri(:) = -(dcosdrj(:) + dcosdrk(:))
438 f_c = ter_f_c(tersoff, drik)
439 f_c_d = ter_f_c_d(tersoff, drik)
440 gterm = ter_g(tersoff, costheta)
441 gterm_d = ter_g_d(tersoff, costheta)
442 expterm = exp((lambda3*(drij - drik))**3)
443 expterm_d = (3.0_dp)*(lambda3**3)*((drij - drik)**2)*expterm
445 dri = f_c_d*gterm*expterm*(rik) &
446 + f_c*gterm_d*expterm*(dcosdri) &
447 + f_c*gterm*expterm_d*(-rij_hat + rik_hat)
450 drj = f_c*gterm_d*expterm*(dcosdrj) &
451 + f_c*gterm*expterm_d*(rij_hat)
453 drk = f_c_d*gterm*expterm*(-rik) &
454 + f_c*gterm_d*expterm*(dcosdrk) &
455 + f_c*gterm*expterm_d*(-rik_hat)
457 f_nonbond(1, iparticle) = f_nonbond(1, iparticle) + prefactor*dri(1)
458 f_nonbond(2, iparticle) = f_nonbond(2, iparticle) + prefactor*dri(2)
459 f_nonbond(3, iparticle) = f_nonbond(3, iparticle) + prefactor*dri(3)
461 f_nonbond(1, jparticle) = f_nonbond(1, jparticle) + prefactor*drj(1)
462 f_nonbond(2, jparticle) = f_nonbond(2, jparticle) + prefactor*drj(2)
463 f_nonbond(3, jparticle) = f_nonbond(3, jparticle) + prefactor*drj(3)
465 f_nonbond(1, kparticle) = f_nonbond(1, kparticle) + prefactor*drk(1)
466 f_nonbond(2, kparticle) = f_nonbond(2, kparticle) + prefactor*drk(2)
467 f_nonbond(3, kparticle) = f_nonbond(3, kparticle) + prefactor*drk(3)
470 pv_nonbond(1, 1) = pv_nonbond(1, 1) + prefactor*(rij(1)*drj(1) + rik(1)*drk(1))
471 pv_nonbond(1, 2) = pv_nonbond(1, 2) + prefactor*(rij(1)*drj(2) + rik(1)*drk(2))
472 pv_nonbond(1, 3) = pv_nonbond(1, 3) + prefactor*(rij(1)*drj(3) + rik(1)*drk(3))
474 pv_nonbond(2, 1) = pv_nonbond(2, 1) + prefactor*(rij(2)*drj(1) + rik(2)*drk(1))
475 pv_nonbond(2, 2) = pv_nonbond(2, 2) + prefactor*(rij(2)*drj(2) + rik(2)*drk(2))
476 pv_nonbond(2, 3) = pv_nonbond(2, 3) + prefactor*(rij(2)*drj(3) + rik(2)*drk(3))
478 pv_nonbond(3, 1) = pv_nonbond(3, 1) + prefactor*(rij(3)*drj(1) + rik(3)*drk(1))
479 pv_nonbond(3, 2) = pv_nonbond(3, 2) + prefactor*(rij(3)*drj(2) + rik(3)*drk(2))
480 pv_nonbond(3, 3) = pv_nonbond(3, 3) + prefactor*(rij(3)*drj(3) + rik(3)*drk(3))
483 END SUBROUTINE ter_zeta_ij_d
492 FUNCTION ter_g(tersoff, costheta)
493 TYPE(tersoff_pot_type),
POINTER :: tersoff
494 REAL(kind=
dp),
INTENT(IN) :: costheta
495 REAL(kind=
dp) :: ter_g
497 REAL(kind=
dp) :: c, c2, d, d2, h
505 ter_g = 1.0_dp + (c2/d2) - (c2)/(d2 + (h - costheta)**2)
516 FUNCTION ter_g_d(tersoff, costheta)
517 TYPE(tersoff_pot_type),
POINTER :: tersoff
518 REAL(kind=
dp),
INTENT(IN) :: costheta
519 REAL(kind=
dp) :: ter_g_d
521 REAL(kind=
dp) :: c, c2, d, d2, h, hc, sintheta
530 sintheta = sqrt(1.0 - costheta**2)
534 ter_g_d = (-2.0_dp*c2*hc)/(d2 + hc**2)**2
556 full_loc_list, loc_cell_v, iparticle, jparticle, f_nonbond, pv_nonbond, &
558 TYPE(tersoff_pot_type),
POINTER :: tersoff
559 TYPE(pos_type),
DIMENSION(:),
POINTER :: r_last_update_pbc
560 REAL(kind=
dp),
DIMENSION(3) :: cell_v
561 INTEGER,
INTENT(IN) :: n_loc_size
562 INTEGER,
DIMENSION(2, 1:n_loc_size) :: full_loc_list
563 REAL(kind=
dp),
DIMENSION(3, 1:n_loc_size) :: loc_cell_v
564 INTEGER,
INTENT(IN) :: iparticle, jparticle
565 REAL(kind=
dp),
DIMENSION(:, :),
INTENT(INOUT) :: f_nonbond, pv_nonbond
566 LOGICAL,
INTENT(IN) :: use_virial
567 REAL(kind=
dp),
INTENT(IN) :: rcutsq
569 CHARACTER(LEN=*),
PARAMETER :: routinen =
'tersoff_forces'
572 REAL(kind=
dp) :: b_ij, b_ij_d, drij, f_a, f_a1, f_a2, &
573 f_a_d, f_c, f_c_d, f_r, f_r1, f_r2, &
574 f_r_d,
fac, prefactor, rij(3), &
577 CALL timeset(routinen, handle)
578 rij(:) = r_last_update_pbc(jparticle)%r(:) - r_last_update_pbc(iparticle)%r(:) + cell_v
579 drij = sqrt(dot_product(rij, rij))
580 rij_hat(:) = rij(:)/drij
583 b_ij = ter_b_ij(tersoff, r_last_update_pbc, iparticle, jparticle, n_loc_size, full_loc_list, loc_cell_v, cell_v, rcutsq)
584 b_ij_d = ter_b_ij_d(tersoff, r_last_update_pbc, iparticle, jparticle, n_loc_size, full_loc_list, loc_cell_v, cell_v, rcutsq)
585 f_a = ter_f_a(tersoff, drij)
586 f_a_d = ter_f_a_d(tersoff, drij)
587 f_c = ter_f_c(tersoff, drij)
588 f_c_d = ter_f_c_d(tersoff, drij)
589 f_r = ter_f_r(tersoff, drij)
590 f_r_d = ter_f_r_d(tersoff, drij)
595 f_nonbond(1, iparticle) = f_nonbond(1, iparticle) + f_r1*rij(1)
596 f_nonbond(2, iparticle) = f_nonbond(2, iparticle) + f_r1*rij(2)
597 f_nonbond(3, iparticle) = f_nonbond(3, iparticle) + f_r1*rij(3)
598 f_nonbond(1, jparticle) = f_nonbond(1, jparticle) - f_r1*rij(1)
599 f_nonbond(2, jparticle) = f_nonbond(2, jparticle) - f_r1*rij(2)
600 f_nonbond(3, jparticle) = f_nonbond(3, jparticle) - f_r1*rij(3)
603 pv_nonbond(1, 1) = pv_nonbond(1, 1) - f_r1*rij(1)*rij(1)
604 pv_nonbond(1, 2) = pv_nonbond(1, 2) - f_r1*rij(1)*rij(2)
605 pv_nonbond(1, 3) = pv_nonbond(1, 3) - f_r1*rij(1)*rij(3)
606 pv_nonbond(2, 1) = pv_nonbond(2, 1) - f_r1*rij(2)*rij(1)
607 pv_nonbond(2, 2) = pv_nonbond(2, 2) - f_r1*rij(2)*rij(2)
608 pv_nonbond(2, 3) = pv_nonbond(2, 3) - f_r1*rij(2)*rij(3)
609 pv_nonbond(3, 1) = pv_nonbond(3, 1) - f_r1*rij(3)*rij(1)
610 pv_nonbond(3, 2) = pv_nonbond(3, 2) - f_r1*rij(3)*rij(2)
611 pv_nonbond(3, 3) = pv_nonbond(3, 3) - f_r1*rij(3)*rij(3)
615 f_nonbond(1, iparticle) = f_nonbond(1, iparticle) + f_r2*rij(1)
616 f_nonbond(2, iparticle) = f_nonbond(2, iparticle) + f_r2*rij(2)
617 f_nonbond(3, iparticle) = f_nonbond(3, iparticle) + f_r2*rij(3)
618 f_nonbond(1, jparticle) = f_nonbond(1, jparticle) - f_r2*rij(1)
619 f_nonbond(2, jparticle) = f_nonbond(2, jparticle) - f_r2*rij(2)
620 f_nonbond(3, jparticle) = f_nonbond(3, jparticle) - f_r2*rij(3)
623 pv_nonbond(1, 1) = pv_nonbond(1, 1) - f_r2*rij(1)*rij(1)
624 pv_nonbond(1, 2) = pv_nonbond(1, 2) - f_r2*rij(1)*rij(2)
625 pv_nonbond(1, 3) = pv_nonbond(1, 3) - f_r2*rij(1)*rij(3)
626 pv_nonbond(2, 1) = pv_nonbond(2, 1) - f_r2*rij(2)*rij(1)
627 pv_nonbond(2, 2) = pv_nonbond(2, 2) - f_r2*rij(2)*rij(2)
628 pv_nonbond(2, 3) = pv_nonbond(2, 3) - f_r2*rij(2)*rij(3)
629 pv_nonbond(3, 1) = pv_nonbond(3, 1) - f_r2*rij(3)*rij(1)
630 pv_nonbond(3, 2) = pv_nonbond(3, 2) - f_r2*rij(3)*rij(2)
631 pv_nonbond(3, 3) = pv_nonbond(3, 3) - f_r2*rij(3)*rij(3)
635 f_a1 = f_c_d*b_ij*f_a*
fac
636 f_nonbond(1, iparticle) = f_nonbond(1, iparticle) + f_a1*rij(1)
637 f_nonbond(2, iparticle) = f_nonbond(2, iparticle) + f_a1*rij(2)
638 f_nonbond(3, iparticle) = f_nonbond(3, iparticle) + f_a1*rij(3)
639 f_nonbond(1, jparticle) = f_nonbond(1, jparticle) - f_a1*rij(1)
640 f_nonbond(2, jparticle) = f_nonbond(2, jparticle) - f_a1*rij(2)
641 f_nonbond(3, jparticle) = f_nonbond(3, jparticle) - f_a1*rij(3)
644 pv_nonbond(1, 1) = pv_nonbond(1, 1) - f_a1*rij(1)*rij(1)
645 pv_nonbond(1, 2) = pv_nonbond(1, 2) - f_a1*rij(1)*rij(2)
646 pv_nonbond(1, 3) = pv_nonbond(1, 3) - f_a1*rij(1)*rij(3)
647 pv_nonbond(2, 1) = pv_nonbond(2, 1) - f_a1*rij(2)*rij(1)
648 pv_nonbond(2, 2) = pv_nonbond(2, 2) - f_a1*rij(2)*rij(2)
649 pv_nonbond(2, 3) = pv_nonbond(2, 3) - f_a1*rij(2)*rij(3)
650 pv_nonbond(3, 1) = pv_nonbond(3, 1) - f_a1*rij(3)*rij(1)
651 pv_nonbond(3, 2) = pv_nonbond(3, 2) - f_a1*rij(3)*rij(2)
652 pv_nonbond(3, 3) = pv_nonbond(3, 3) - f_a1*rij(3)*rij(3)
656 f_a2 = f_c*b_ij*f_a_d*
fac
657 f_nonbond(1, iparticle) = f_nonbond(1, iparticle) + f_a2*rij(1)
658 f_nonbond(2, iparticle) = f_nonbond(2, iparticle) + f_a2*rij(2)
659 f_nonbond(3, iparticle) = f_nonbond(3, iparticle) + f_a2*rij(3)
660 f_nonbond(1, jparticle) = f_nonbond(1, jparticle) - f_a2*rij(1)
661 f_nonbond(2, jparticle) = f_nonbond(2, jparticle) - f_a2*rij(2)
662 f_nonbond(3, jparticle) = f_nonbond(3, jparticle) - f_a2*rij(3)
665 pv_nonbond(1, 1) = pv_nonbond(1, 1) - f_a2*rij(1)*rij(1)
666 pv_nonbond(1, 2) = pv_nonbond(1, 2) - f_a2*rij(1)*rij(2)
667 pv_nonbond(1, 3) = pv_nonbond(1, 3) - f_a2*rij(1)*rij(3)
668 pv_nonbond(2, 1) = pv_nonbond(2, 1) - f_a2*rij(2)*rij(1)
669 pv_nonbond(2, 2) = pv_nonbond(2, 2) - f_a2*rij(2)*rij(2)
670 pv_nonbond(2, 3) = pv_nonbond(2, 3) - f_a2*rij(2)*rij(3)
671 pv_nonbond(3, 1) = pv_nonbond(3, 1) - f_a2*rij(3)*rij(1)
672 pv_nonbond(3, 2) = pv_nonbond(3, 2) - f_a2*rij(3)*rij(2)
673 pv_nonbond(3, 3) = pv_nonbond(3, 3) - f_a2*rij(3)*rij(3)
677 prefactor = f_c*b_ij_d*f_a*
fac
678 CALL ter_zeta_ij_d(tersoff, r_last_update_pbc, iparticle, jparticle, f_nonbond, pv_nonbond, prefactor, &
679 n_loc_size, full_loc_list, loc_cell_v, cell_v, rcutsq, use_virial)
680 CALL timestop(handle)
696 TYPE(fist_neighbor_type),
POINTER :: nonbonded
697 TYPE(pair_potential_pp_type),
POINTER :: potparm
698 INTEGER,
DIMENSION(:, :),
POINTER :: glob_loc_list
699 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: glob_cell_v
700 INTEGER,
DIMENSION(:),
POINTER :: glob_loc_list_a
701 TYPE(cell_type),
POINTER :: cell
703 CHARACTER(LEN=*),
PARAMETER :: routinen =
'setup_tersoff_arrays'
705 INTEGER :: handle, i, iend, igrp, ikind, ilist, &
706 ipair, istart, jkind, nkinds, npairs, &
708 INTEGER,
DIMENSION(:),
POINTER :: work_list, work_list2
709 INTEGER,
DIMENSION(:, :),
POINTER ::
list
710 REAL(kind=
dp),
DIMENSION(3) :: cell_v, cvi
711 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: rwork_list
712 TYPE(neighbor_kind_pairs_type),
POINTER :: neighbor_kind_pair
713 TYPE(pair_potential_single_type),
POINTER :: pot
715 cpassert(.NOT.
ASSOCIATED(glob_loc_list))
716 cpassert(.NOT.
ASSOCIATED(glob_loc_list_a))
717 cpassert(.NOT.
ASSOCIATED(glob_cell_v))
718 CALL timeset(routinen, handle)
720 nkinds =
SIZE(potparm%pot, 1)
721 DO ilist = 1, nonbonded%nlists
722 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
723 npairs = neighbor_kind_pair%npairs
724 IF (npairs == 0) cycle
725 kind_group_loop1:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
726 istart = neighbor_kind_pair%grp_kind_start(igrp)
727 iend = neighbor_kind_pair%grp_kind_end(igrp)
728 ikind = neighbor_kind_pair%ij_kind(1, igrp)
729 jkind = neighbor_kind_pair%ij_kind(2, igrp)
730 pot => potparm%pot(ikind, jkind)%pot
731 npairs = iend - istart + 1
733 DO i = 1,
SIZE(pot%type)
734 IF (pot%type(i) ==
tersoff_type) npairs_tot = npairs_tot + npairs
736 END DO kind_group_loop1
738 ALLOCATE (work_list(npairs_tot))
739 ALLOCATE (work_list2(npairs_tot))
740 ALLOCATE (glob_loc_list(2, npairs_tot))
741 ALLOCATE (glob_cell_v(3, npairs_tot))
744 DO ilist = 1, nonbonded%nlists
745 neighbor_kind_pair => nonbonded%neighbor_kind_pairs(ilist)
746 npairs = neighbor_kind_pair%npairs
747 IF (npairs == 0) cycle
748 kind_group_loop2:
DO igrp = 1, neighbor_kind_pair%ngrp_kind
749 istart = neighbor_kind_pair%grp_kind_start(igrp)
750 iend = neighbor_kind_pair%grp_kind_end(igrp)
751 ikind = neighbor_kind_pair%ij_kind(1, igrp)
752 jkind = neighbor_kind_pair%ij_kind(2, igrp)
753 list => neighbor_kind_pair%list
754 cvi = neighbor_kind_pair%cell_vector
755 pot => potparm%pot(ikind, jkind)%pot
756 npairs = iend - istart + 1
758 cell_v = matmul(cell%hmat, cvi)
759 DO i = 1,
SIZE(pot%type)
763 glob_loc_list(:, npairs_tot + ipair) =
list(:, istart - 1 + ipair)
764 glob_cell_v(1:3, npairs_tot + ipair) = cell_v(1:3)
766 npairs_tot = npairs_tot + npairs
769 END DO kind_group_loop2
772 CALL sort(glob_loc_list(1, :), npairs_tot, work_list)
773 DO ipair = 1, npairs_tot
774 work_list2(ipair) = glob_loc_list(2, work_list(ipair))
776 glob_loc_list(2, :) = work_list2
777 DEALLOCATE (work_list2)
778 ALLOCATE (rwork_list(3, npairs_tot))
779 DO ipair = 1, npairs_tot
780 rwork_list(:, ipair) = glob_cell_v(:, work_list(ipair))
782 glob_cell_v = rwork_list
783 DEALLOCATE (rwork_list)
784 DEALLOCATE (work_list)
785 ALLOCATE (glob_loc_list_a(npairs_tot))
786 glob_loc_list_a = glob_loc_list(1, :)
787 CALL timestop(handle)
800 INTEGER,
DIMENSION(:, :),
POINTER :: glob_loc_list
801 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: glob_cell_v
802 INTEGER,
DIMENSION(:),
POINTER :: glob_loc_list_a
804 IF (
ASSOCIATED(glob_loc_list))
THEN
805 DEALLOCATE (glob_loc_list)
807 IF (
ASSOCIATED(glob_loc_list_a))
THEN
808 DEALLOCATE (glob_loc_list_a)
810 IF (
ASSOCIATED(glob_cell_v))
THEN
811 DEALLOCATE (glob_cell_v)
static GRID_HOST_DEVICE double fac(const int i)
Factorial function, e.g. fac(5) = 5! = 120.
Handles all functions related to the CELL.
Define the neighbor list data types and the corresponding functionality.
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_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)
...
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
integer, parameter, public tersoff_type
All kind of helpful little routines.