imptorsion_geomgrad Subroutine

public subroutine imptorsion_geomgrad(bds, grad)

Uses

  • proc~~imptorsion_geomgrad~~UsesGraph proc~imptorsion_geomgrad imptorsion_geomgrad module~mod_jacobian_mat mod_jacobian_mat proc~imptorsion_geomgrad->module~mod_jacobian_mat module~mod_memory mod_memory module~mod_jacobian_mat->module~mod_memory module~mod_io mod_io module~mod_memory->module~mod_io iso_c_binding iso_c_binding module~mod_memory->iso_c_binding module~mod_constants mod_constants module~mod_memory->module~mod_constants module~mod_io->module~mod_constants module~mod_constants->iso_c_binding

Compute torsion potential

Arguments

Type IntentOptional 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


Calls

proc~~imptorsion_geomgrad~~CallsGraph proc~imptorsion_geomgrad imptorsion_geomgrad proc~torsion_angle_jacobian torsion_angle_jacobian proc~imptorsion_geomgrad->proc~torsion_angle_jacobian proc~cross_product cross_product proc~torsion_angle_jacobian->proc~cross_product proc~versor_der versor_der proc~torsion_angle_jacobian->proc~versor_der proc~vec_skw vec_skw proc~torsion_angle_jacobian->proc~vec_skw

Called by

proc~~imptorsion_geomgrad~~CalledByGraph proc~imptorsion_geomgrad imptorsion_geomgrad proc~ommp_imptorsion_geomgrad ommp_imptorsion_geomgrad proc~ommp_imptorsion_geomgrad->proc~imptorsion_geomgrad proc~ommp_full_bnd_geomgrad ommp_full_bnd_geomgrad proc~ommp_full_bnd_geomgrad->proc~imptorsion_geomgrad proc~c_ommp_imptorsion_geomgrad C_ommp_imptorsion_geomgrad proc~c_ommp_imptorsion_geomgrad->proc~ommp_imptorsion_geomgrad proc~c_ommp_full_bnd_geomgrad C_ommp_full_bnd_geomgrad proc~c_ommp_full_bnd_geomgrad->proc~ommp_full_bnd_geomgrad proc~ommp_full_geomgrad ommp_full_geomgrad proc~ommp_full_geomgrad->proc~ommp_full_bnd_geomgrad proc~c_ommp_full_geomgrad C_ommp_full_geomgrad proc~c_ommp_full_geomgrad->proc~ommp_full_geomgrad proc~ommptest_fakeqm_internal_geomgrad ommptest_fakeqm_internal_geomgrad proc~ommptest_fakeqm_internal_geomgrad->proc~ommp_full_geomgrad proc~ommptest_totalqmmm_geomgrad ommptest_totalqmmm_geomgrad proc~ommptest_totalqmmm_geomgrad->proc~ommp_full_geomgrad proc~ommptest_fakeqm_linkatom_geomgrad ommptest_fakeqm_linkatom_geomgrad proc~ommptest_fakeqm_linkatom_geomgrad->proc~ommp_full_geomgrad

Contents

Source Code


Source Code

    subroutine imptorsion_geomgrad(bds, grad)
        !! Compute torsion potential
        use mod_jacobian_mat, only: torsion_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)
        !! 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)
        integer(ip) :: i, j, ia, ib, ic, id
        logical :: sk_a, sk_b, sk_c, sk_d

        if (.not. bds%use_imptorsion) return

        !$omp parallel do default(shared) &
        !$omp private(i, ia, ib, ic, id, sk_a, sk_b, sk_c, sk_d, j, thet, J_a, J_b, J_c, J_d, g)
        do i = 1, bds%nimptorsion
            ! Atoms that define the dihedral angle
            ia = bds%imptorsionat(1, i)
            ib = bds%imptorsionat(2, i)
            ic = bds%imptorsionat(3, i)
            id = bds%imptorsionat(4, 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 torsion_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)

            do j = 1, 3
                if (bds%imptorsn(j, i) < 1) exit
                g = -real(bds%imptorsn(j, i)) * sin(real(bds%imptorsn(j, i)) * thet &
                                                    - bds%imptorsphase(j, i)) &
                                                * bds%imptorsamp(j, i)
                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 do
    end subroutine imptorsion_geomgrad