Compute torsion potential
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(ommp_bonded_type), | intent(in) | :: | bds | |||
real(kind=rp), | intent(inout) | :: | V |
torsion potential, result will be added to V |
subroutine torsion_potential(bds, V)
!! Compute torsion potential
use mod_constants, only: pi, eps_rp
implicit none
type(ommp_bonded_type), intent(in) :: bds
! Bonded potential data structure
real(rp), intent(inout) :: V
!! torsion potential, result will be added to V
real(rp) :: thet, costhet
integer(ip) :: i, j
if(.not. bds%use_torsion) return
!$omp parallel do default(shared) &
!$omp private(i,costhet,thet,j) reduction(+:V)
do i=1, bds%ntorsion
! Atoms that defines the dihedral angle
costhet = cos_torsion(bds%top, bds%torsionat(:,i))
if(costhet + 1.0 <= eps_rp) then
thet = pi
else if(abs(costhet - 1.0) <= eps_rp) then
thet = 0.0
else
thet = acos(costhet)
end if
do j=1, 6
if(bds%torsn(j,i) < 1) exit
V = V + bds%torsamp(j,i) * (1+cos(real(bds%torsn(j,i))*thet &
- bds%torsphase(j,i)))
end do
end do
end subroutine torsion_potential