potential_D2E Subroutine

public subroutine potential_D2E(eel, cpt, V, amoeba_P_insted_of_D_)

This subroutine computes the potential generated by the induced point dipoles to a set of arbitrary coordinates, without applying any screening rules. Note: for AMOEBA D dipoles should be used.

Arguments

Type IntentOptional Attributes Name
type(ommp_electrostatics_type), intent(in) :: eel

Electrostatics data structure

real(kind=rp), intent(in) :: cpt(:,:)

Coordinates at which the electric field is requested

real(kind=rp), intent(inout) :: V(:)

Electric field (results will be added)

logical, intent(in), optional :: amoeba_P_insted_of_D_

For AMOEBA FF, if true the potential of P dipoles is computed, otherwise potential of D dipoles is computed


Calls

proc~~potential_d2e~~CallsGraph proc~potential_d2e potential_D2E proc~fatal_error fatal_error proc~potential_d2e->proc~fatal_error proc~mu_elec_prop mu_elec_prop proc~potential_d2e->proc~mu_elec_prop proc~coulomb_kernel coulomb_kernel proc~potential_d2e->proc~coulomb_kernel proc~ommp_message ommp_message proc~fatal_error->proc~ommp_message proc~close_output close_output proc~fatal_error->proc~close_output proc~coulomb_kernel->proc~fatal_error proc~close_output->proc~ommp_message

Called by

proc~~potential_d2e~~CalledByGraph proc~potential_d2e potential_D2E proc~electrostatic_for_ene electrostatic_for_ene proc~electrostatic_for_ene->proc~potential_d2e proc~ommp_potential_mmpol2ext ommp_potential_mmpol2ext proc~ommp_potential_mmpol2ext->proc~potential_d2e proc~ommp_potential_pol2ext ommp_potential_pol2ext proc~ommp_potential_pol2ext->proc~potential_d2e proc~c_ommp_prepare_qm_ele_ene C_ommp_prepare_qm_ele_ene proc~c_ommp_prepare_qm_ele_ene->proc~electrostatic_for_ene proc~c_ommp_potential_pol2ext C_ommp_potential_pol2ext proc~c_ommp_potential_pol2ext->proc~ommp_potential_pol2ext proc~c_ommp_potential_mmpol2ext C_ommp_potential_mmpol2ext proc~c_ommp_potential_mmpol2ext->proc~ommp_potential_mmpol2ext

Contents

Source Code


Source Code

    subroutine potential_D2E(eel, cpt, V, amoeba_P_insted_of_D_)
        !! This subroutine computes the potential generated by the induced 
        !! point dipoles to a set of arbitrary coordinates, without applying
        !! any screening rules. Note: for AMOEBA D dipoles should be used. 
        
        implicit none

        type(ommp_electrostatics_type), intent(in) :: eel
        !! Electrostatics data structure
        real(rp), intent(inout) :: V(:)
        !! Electric field (results will be added)
        real(rp), intent(in) :: cpt(:,:)
        !! Coordinates at which the electric field is requested
        logical, optional, intent(in) :: amoeba_P_insted_of_D_
        !! For AMOEBA FF, if true the potential of P dipoles
        !! is computed, otherwise potential of D dipoles is computed

        integer(ip) :: i, j, n_cpt
        logical :: amoeba_P_insted_of_D
        real(rp) :: kernel(5), dr(3), tmpV, tmpE(3), tmpEgr(6), tmpHE(10)

        if(eel%pol_atoms < 1) return

        if(present(amoeba_P_insted_of_D_)) then
            amoeba_P_insted_of_D = amoeba_P_insted_of_D_
        else
            amoeba_P_insted_of_D = .false.
        end if

        if(.not. eel%ipd_done) call fatal_error("IPD should be computed before&
                                                & D2E potential.")
        n_cpt = size(cpt, 2)

        if(eel%amoeba) then
            if(.not. amoeba_P_insted_of_D) then
                !$omp parallel do default(shared) schedule(dynamic) collapse(2) &
                !$omp private(i,j,dr,kernel,tmpV,tmpE,tmpEgr,tmpHE) reduction(+:V)
                do i=1, eel%pol_atoms
                    do j=1, n_cpt
                        dr = cpt(:,j) - eel%cpol(:,i)
                        call coulomb_kernel(dr, 2, kernel(1:2))
                        tmpV = 0.0_rp
                        
                        call mu_elec_prop(eel%ipd(:,i,_amoeba_D_), &
                                        dr, kernel, .true., tmpV, &
                                        .false., tmpE, .false., tmpEgr, & 
                                        .false., tmpHE)

                        V(j) = V(j) + tmpV
                    end do
                end do
            else
                !$omp parallel do default(shared) schedule(dynamic) collapse(2) &
                !$omp private(i,j,dr,kernel,tmpV,tmpE,tmpEgr,tmpHE) reduction(+:V)
                do i=1, eel%pol_atoms
                    do j=1, n_cpt
                        dr = cpt(:,j) - eel%cpol(:,i)
                        call coulomb_kernel(dr, 2, kernel(1:2))
                        tmpV = 0.0_rp
                        
                        call mu_elec_prop(eel%ipd(:,i,_amoeba_P_), &
                                        dr, kernel, .true., tmpV, &
                                        .false., tmpE, .false., tmpEgr, & 
                                        .false., tmpHE)

                        V(j) = V(j) + tmpV
                    end do
                end do
            end if
        else
            !$omp parallel do default(shared) schedule(dynamic) collapse(2) &
            !$omp private(i,j,dr,kernel,tmpV,tmpE,tmpEgr,tmpHE) reduction(+:V)
            do i=1, eel%pol_atoms
                ! loop on sources
                do j=1, n_cpt
                    dr = cpt(:,j) - eel%cpol(:,i)
                    call coulomb_kernel(dr, 2, kernel(1:2))
                    tmpV = 0.0_rp
                    
                    call mu_elec_prop(eel%ipd(:,i,1), &
                                      dr, kernel, .true., tmpV, &
                                      .false., tmpE, .false., tmpEgr, & 
                                      .false., tmpHE)
                    
                    V(j) = V(j) + tmpV
                end do
            end do
        end if
    end subroutine potential_D2E