Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(ommp_bonded_type), | intent(in) | :: | bds | |||
real(kind=rp), | intent(inout) | :: | grad(3,bds%top%mm_atoms) |
Gradients of bond stretching terms of potential energy |
subroutine opb_geomgrad(bds, grad)
use mod_jacobian_mat, only: opb_angle_jacobian
implicit none
type(ommp_bonded_type), intent(in) :: bds
! Bonded potential data structure
real(rp), intent(inout) :: grad(3,bds%top%mm_atoms)
!! Gradients of bond stretching terms of potential energy
real(rp) :: thet, g, J_a(3), J_b(3), J_c(3), J_d(3)
integer(ip) :: i, ia, ib, ic, id
logical :: sk_a, sk_b, sk_c, sk_d
if(.not. bds%use_opb) return
!$omp parallel do default(shared) schedule(dynamic)&
!$omp private(i,ia,ib,ic,id,sk_a,sk_b,sk_c,sk_d,thet,J_a,J_b,J_c,J_d,g)
do i=1, bds%nopb
ia = bds%opbat(2,i)
ib = bds%opbat(4,i)
ic = bds%opbat(3,i)
id = bds%opbat(1,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)
if(sk_a .and. sk_b .and. sk_c .and. sk_d) cycle
else
sk_a = .false.
sk_b = .false.
sk_c = .false.
sk_d = .false.
end if
call opb_angle_jacobian(bds%top%cmm(:,ia), &
bds%top%cmm(:,ib), &
bds%top%cmm(:,ic), &
bds%top%cmm(:,id), &
thet, J_a, J_b, J_c, J_d)
g = bds%kopb(i) * thet * (2.0 + 3.0*bds%opb_cubic*thet &
+ 4.0*bds%opb_quartic*thet**2 + 5.0*bds%opb_pentic*thet**3 &
+ 6.0*bds%opb_sextic*thet**4)
if(.not. sk_a) then
!$omp atomic update
grad(1,ia) = grad(1,ia) + J_a(1) * g
!$omp atomic update
grad(2,ia) = grad(2,ia) + J_a(2) * g
!$omp atomic update
grad(3,ia) = grad(3,ia) + J_a(3) * g
end if
if(.not. sk_b) then
!$omp atomic update
grad(1,ib) = grad(1,ib) + J_b(1) * g
!$omp atomic update
grad(2,ib) = grad(2,ib) + J_b(2) * g
!$omp atomic update
grad(3,ib) = grad(3,ib) + J_b(3) * g
end if
if(.not. sk_c) then
!$omp atomic update
grad(1,ic) = grad(1,ic) + J_c(1) * g
!$omp atomic update
grad(2,ic) = grad(2,ic) + J_c(2) * g
!$omp atomic update
grad(3,ic) = grad(3,ic) + J_c(3) * g
end if
if(.not. sk_d) then
!$omp atomic update
grad(1,id) = grad(1,id) + J_d(1) * g
!$omp atomic update
grad(2,id) = grad(2,id) + J_d(2) * g
!$omp atomic update
grad(3,id) = grad(3,id) + J_d(3) * g
end if
end do
end subroutine opb_geomgrad