strbnd_geomgrad Subroutine

public subroutine strbnd_geomgrad(bds, grad)

Uses

  • proc~~strbnd_geomgrad~~UsesGraph proc~strbnd_geomgrad strbnd_geomgrad module~mod_jacobian_mat mod_jacobian_mat proc~strbnd_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

Arguments

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


Calls

proc~~strbnd_geomgrad~~CallsGraph proc~strbnd_geomgrad strbnd_geomgrad proc~rij_jacobian Rij_jacobian proc~strbnd_geomgrad->proc~rij_jacobian proc~simple_angle_jacobian simple_angle_jacobian proc~strbnd_geomgrad->proc~simple_angle_jacobian

Called by

proc~~strbnd_geomgrad~~CalledByGraph proc~strbnd_geomgrad strbnd_geomgrad proc~ommp_strbnd_geomgrad ommp_strbnd_geomgrad proc~ommp_strbnd_geomgrad->proc~strbnd_geomgrad proc~ommp_full_bnd_geomgrad ommp_full_bnd_geomgrad proc~ommp_full_bnd_geomgrad->proc~strbnd_geomgrad proc~c_ommp_strbnd_geomgrad C_ommp_strbnd_geomgrad proc~c_ommp_strbnd_geomgrad->proc~ommp_strbnd_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 strbnd_geomgrad(bds, grad)
        use mod_jacobian_mat, only: Rij_jacobian, simple_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

        integer(ip) :: i, ia, ib, ic
        real(rp) :: d_l1, d_l2, d_thet, l1, l2, thet, g1, g2, g3
        real(rp), dimension(3) :: a, b, c, &
                                  J1_a, J1_b, &
                                  J2_b, J2_c, &
                                  J3_a, J3_b, J3_c
        logical :: sk_a, sk_b, sk_c
        
        if(.not. bds%use_strbnd) return

        !$omp parallel do default(shared) schedule(dynamic) &
        !$omp private(i,ia,ib,ic,sk_a,sk_b,sk_c,a,b,c,l1,l2,d_l1,d_l2,thet,d_thet) &
        !$omp private(J1_a,J1_b,J2_b,J2_c,J3_a,J3_b,J3_c,g1,g2,g3)
        do i=1, bds%nstrbnd
            ia = bds%strbndat(1,i)
            ib = bds%strbndat(2,i)
            ic = bds%strbndat(3,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)
                if(sk_a .and. sk_b .and. sk_c) cycle
            else
                sk_a = .false.
                sk_b = .false.
                sk_c = .false.
            end if

            a = bds%top%cmm(:, ia)
            b = bds%top%cmm(:, ib)
            c = bds%top%cmm(:, ic)

            call Rij_jacobian(a, b, l1, J1_a, J1_b)
            call Rij_jacobian(b, c, l2, J2_b, J2_c)
            call simple_angle_jacobian(a, b, c, thet, J3_a, J3_b, J3_c)
            
            d_l1 = l1 - bds%strbndl10(i)
            d_l2 = l2 - bds%strbndl20(i)
            d_thet = thet - bds%strbndthet0(i) 
           
            g1 = bds%strbndk1(i) * d_thet
            g2 = bds%strbndk2(i) * d_thet
            g3 = bds%strbndk1(i) * d_l1 + bds%strbndk2(i) * d_l2

            if(.not. sk_a) then
                !$omp atomic update
                grad(1,ia) = grad(1,ia) + J1_a(1) * g1 + J3_a(1) * g3
                !$omp atomic update
                grad(2,ia) = grad(2,ia) + J1_a(2) * g1 + J3_a(2) * g3
                !$omp atomic update
                grad(3,ia) = grad(3,ia) + J1_a(3) * g1 + J3_a(3) * g3
            end if

            if(.not. sk_b) then
                !$omp atomic update
                grad(1,ib) = grad(1,ib) + J1_b(1) * g1 + J2_b(1) * g2 + J3_b(1) * g3
                !$omp atomic update
                grad(2,ib) = grad(2,ib) + J1_b(2) * g1 + J2_b(2) * g2 + J3_b(2) * g3
                !$omp atomic update
                grad(3,ib) = grad(3,ib) + J1_b(3) * g1 + J2_b(3) * g2 + J3_b(3) * g3
            end if

            if(.not. sk_c) then
                !$omp atomic update
                grad(1,ic) = grad(1,ic) + J2_c(1) * g2 + J3_c(1) * g3
                !$omp atomic update
                grad(2,ic) = grad(2,ic) + J2_c(2) * g2 + J3_c(2) * g3
                !$omp atomic update
                grad(3,ic) = grad(3,ic) + J2_c(3) * g2 + J3_c(3) * g3
            end if
        end do

    end subroutine strbnd_geomgrad