| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| type(ommp_bonded_type), | intent(in) | :: | bds | |||
| real(kind=rp), | intent(inout) | :: | grad(3,bds%top%mm_atoms) |
improper torsion potential, result will be added to V |
subroutine pitors_geomgrad(bds, grad)
use mod_jacobian_mat, only: pitors_angle_jacobian
use mod_constants, only : pi
implicit none
type(ommp_bonded_type), intent(in) :: bds
! Bonded potential data structure
real(rp), intent(inout) :: grad(3,bds%top%mm_atoms)
!! improper torsion potential, result will be added to V
real(rp) :: thet, g, J_a(3), J_b(3), J_c(3), J_d(3), J_e(3), J_f(3)
integer(ip) :: i, ia, ib, ic, id, ie, if_
logical :: sk_a, sk_b, sk_c, sk_d, sk_e, sk_f
if(.not. bds%use_pitors) return
!$omp parallel do default(shared) schedule(dynamic) &
!$omp private(i,ia,ib,ic,id,ie,if_,sk_a,sk_b,sk_c,sk_d,sk_e,sk_f) &
!$omp private(J_a,J_b,J_c,J_d,J_e,J_f,g,thet)
do i=1, bds%npitors
ia = bds%pitorsat(1,i)
ic = bds%pitorsat(2,i)
id = bds%pitorsat(3,i)
ib = bds%pitorsat(4,i)
ie = bds%pitorsat(5,i)
if_ = bds%pitorsat(6,i)
if(bds%top%use_frozen) then
sk_a = bds%top%frozen(ia)
sk_b = bds%top%frozen(ib)
sk_c = bds%top%frozen(ic)
sk_d = bds%top%frozen(id)
sk_e = bds%top%frozen(ie)
sk_f = bds%top%frozen(if_)
if(sk_a .and. sk_b .and. sk_c .and. sk_d .and. sk_e .and. sk_f) cycle
else
sk_a = .false.
sk_b = .false.
sk_c = .false.
sk_d = .false.
sk_e = .false.
sk_f = .false.
end if
call pitors_angle_jacobian(bds%top%cmm(:,ia), &
bds%top%cmm(:,ib), &
bds%top%cmm(:,ic), &
bds%top%cmm(:,id), &
bds%top%cmm(:,ie), &
bds%top%cmm(:,if_), &
thet, J_a, J_b, J_c, J_d, J_e, J_f)
g = -2.0 * bds%kpitors(i) * sin(2.0*thet-pi)
if(.not. sk_a) then
!$omp atomic update
grad(1,ia) = grad(1,ia) + g * J_a(1)
!$omp atomic update
grad(2,ia) = grad(2,ia) + g * J_a(2)
!$omp atomic update
grad(3,ia) = grad(3,ia) + g * J_a(3)
end if
if(.not. sk_b) then
!$omp atomic update
grad(1,ib) = grad(1,ib) + g * J_b(1)
!$omp atomic update
grad(2,ib) = grad(2,ib) + g * J_b(2)
!$omp atomic update
grad(3,ib) = grad(3,ib) + g * J_b(3)
end if
if(.not. sk_c) then
!$omp atomic update
grad(1,ic) = grad(1,ic) + g * J_c(1)
!$omp atomic update
grad(2,ic) = grad(2,ic) + g * J_c(2)
!$omp atomic update
grad(3,ic) = grad(3,ic) + g * J_c(3)
end if
if(.not. sk_d) then
!$omp atomic update
grad(1,id) = grad(1,id) + g * J_d(1)
!$omp atomic update
grad(2,id) = grad(2,id) + g * J_d(2)
!$omp atomic update
grad(3,id) = grad(3,id) + g * J_d(3)
end if
if(.not. sk_e) then
!$omp atomic update
grad(1,ie) = grad(1,ie) + g * J_e(1)
!$omp atomic update
grad(2,ie) = grad(2,ie) + g * J_e(2)
!$omp atomic update
grad(3,ie) = grad(3,ie) + g * J_e(3)
end if
if(.not. sk_f) then
!$omp atomic update
grad(1,if_) = grad(1,if_) + g * J_f(1)
!$omp atomic update
grad(2,if_) = grad(2,if_) + g * J_f(2)
!$omp atomic update
grad(3,if_) = grad(3,if_) + g * J_f(3)
end if
end do
end subroutine pitors_geomgrad