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.
Type | Intent | Optional | 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 |
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