q_elec_prop Subroutine

public subroutine q_elec_prop(q, dr, kernel, do_V, V, do_E, E, do_grdE, grdE, do_HE, HE)

TODO Computes the electric potential of a charge at position from the charge itself. Pre-computed kernel should be provided as input. The result is added to .

Arguments

Type IntentOptional Attributes Name
real(kind=rp), intent(in) :: q

Charge

real(kind=rp), intent(in) :: dr(3)

Distance vector

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

Array of coulomb kernel (either damped or undamped)

logical, intent(in) :: do_V

Flags to enable/disable calculation of different electrostatic properties

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

Electric potential

logical, intent(in) :: do_E

Flags to enable/disable calculation of different electrostatic properties

real(kind=rp), intent(inout) :: E(3)

Electric potential

logical, intent(in) :: do_grdE

Flags to enable/disable calculation of different electrostatic properties

real(kind=rp), intent(inout) :: grdE(6)

Electric potential

logical, intent(in) :: do_HE

Flags to enable/disable calculation of different electrostatic properties

real(kind=rp), intent(inout) :: HE(10)

Electric potential


Called by

proc~~q_elec_prop~~CalledByGraph proc~q_elec_prop q_elec_prop proc~elec_prop_m2m elec_prop_M2M proc~elec_prop_m2m->proc~q_elec_prop proc~electrostatic_for_ene electrostatic_for_ene proc~electrostatic_for_ene->proc~q_elec_prop proc~potential_m2e potential_M2E proc~electrostatic_for_ene->proc~potential_m2e proc~field_m2e field_M2E proc~field_m2e->proc~q_elec_prop proc~potential_m2e->proc~q_elec_prop proc~elec_prop_m2d elec_prop_M2D proc~elec_prop_m2d->proc~q_elec_prop proc~electrostatic_for_grad electrostatic_for_grad proc~electrostatic_for_grad->proc~q_elec_prop proc~electrostatic_for_grad->proc~field_m2e proc~prepare_fixedelec prepare_fixedelec proc~prepare_fixedelec->proc~elec_prop_m2m 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~ommp_field_mm2ext ommp_field_mm2ext proc~ommp_field_mm2ext->proc~field_m2e proc~c_ommp_field_mm2ext C_ommp_field_mm2ext proc~c_ommp_field_mm2ext->proc~field_m2e proc~ommp_field_mmpol2ext ommp_field_mmpol2ext proc~ommp_field_mmpol2ext->proc~field_m2e proc~ommp_potential_mmpol2ext ommp_potential_mmpol2ext proc~ommp_potential_mmpol2ext->proc~potential_m2e proc~ommp_potential_mm2ext ommp_potential_mm2ext proc~ommp_potential_mm2ext->proc~potential_m2e proc~prepare_polelec prepare_polelec proc~prepare_polelec->proc~elec_prop_m2d proc~c_ommp_prepare_qm_ele_grd C_ommp_prepare_qm_ele_grd proc~c_ommp_prepare_qm_ele_grd->proc~electrostatic_for_grad proc~fixedelec_geomgrad fixedelec_geomgrad proc~fixedelec_geomgrad->proc~prepare_fixedelec proc~c_ommp_field_mmpol2ext C_ommp_field_mmpol2ext proc~c_ommp_field_mmpol2ext->proc~ommp_field_mmpol2ext proc~energy_mm_mm energy_MM_MM proc~energy_mm_mm->proc~prepare_fixedelec proc~c_ommp_potential_mm2ext C_ommp_potential_mm2ext proc~c_ommp_potential_mm2ext->proc~ommp_potential_mm2ext proc~c_ommp_potential_mmpol2ext C_ommp_potential_mmpol2ext proc~c_ommp_potential_mmpol2ext->proc~ommp_potential_mmpol2ext proc~polelec_geomgrad polelec_geomgrad proc~polelec_geomgrad->proc~prepare_polelec proc~ommp_get_polelec_energy ommp_get_polelec_energy proc~ommp_get_polelec_energy->proc~prepare_polelec proc~energy_mm_pol energy_MM_pol proc~ommp_get_polelec_energy->proc~energy_mm_pol proc~ommp_set_external_field ommp_set_external_field proc~ommp_set_external_field->proc~prepare_polelec proc~energy_mm_pol->proc~prepare_polelec proc~ommp_fixedelec_geomgrad ommp_fixedelec_geomgrad proc~ommp_fixedelec_geomgrad->proc~fixedelec_geomgrad proc~ommp_polelec_geomgrad ommp_polelec_geomgrad proc~ommp_polelec_geomgrad->proc~polelec_geomgrad proc~ommp_get_fixedelec_energy ommp_get_fixedelec_energy proc~ommp_get_fixedelec_energy->proc~energy_mm_mm proc~ommp_get_full_ele_energy ommp_get_full_ele_energy proc~ommp_get_full_ele_energy->proc~ommp_get_polelec_energy proc~ommp_get_full_ele_energy->proc~ommp_get_fixedelec_energy proc~ommp_full_geomgrad ommp_full_geomgrad proc~ommp_full_geomgrad->proc~fixedelec_geomgrad proc~ommp_full_geomgrad->proc~polelec_geomgrad proc~c_ommp_get_polelec_energy C_ommp_get_polelec_energy proc~c_ommp_get_polelec_energy->proc~ommp_get_polelec_energy proc~ommp_set_external_field_nomm ommp_set_external_field_nomm proc~ommp_set_external_field_nomm->proc~ommp_set_external_field proc~c_ommp_set_external_field_nomm C_ommp_set_external_field_nomm proc~c_ommp_set_external_field_nomm->proc~ommp_set_external_field proc~c_ommp_set_external_field C_ommp_set_external_field proc~c_ommp_set_external_field->proc~ommp_set_external_field proc~c_ommp_fixedelec_geomgrad C_ommp_fixedelec_geomgrad proc~c_ommp_fixedelec_geomgrad->proc~ommp_fixedelec_geomgrad proc~c_ommp_polelec_geomgrad C_ommp_polelec_geomgrad proc~c_ommp_polelec_geomgrad->proc~ommp_polelec_geomgrad proc~c_ommp_get_full_ele_energy C_ommp_get_full_ele_energy proc~c_ommp_get_full_ele_energy->proc~ommp_get_full_ele_energy proc~c_ommp_get_fixedelec_energy C_ommp_get_fixedelec_energy proc~c_ommp_get_fixedelec_energy->proc~ommp_get_fixedelec_energy proc~ommp_get_full_energy ommp_get_full_energy proc~ommp_get_full_energy->proc~ommp_get_full_ele_energy proc~c_ommp_full_geomgrad C_ommp_full_geomgrad proc~c_ommp_full_geomgrad->proc~ommp_full_geomgrad proc~c_ommp_get_full_energy C_ommp_get_full_energy proc~c_ommp_get_full_energy->proc~ommp_get_full_energy

Contents

Source Code


Source Code

    subroutine q_elec_prop(q, dr, kernel, &
                           do_V, V, do_E, E, do_grdE, grdE, do_HE, HE)
        !! TODO
        !! Computes the electric potential of a charge \(q\) at position
        !! \(\mathbf{dr}\) from the charge itself. Pre-computed kernel should
        !! be provided as input. 
        !! The result is added to \(V\).
        !! $$ V_q = q \frac{f(r)}{||\mathbf{dr}||} $$

        implicit none

        real(rp), intent(in) :: q
        !! Charge
        real(rp), intent(in) :: dr(3)
        !! Distance vector
        real(rp), intent(in) :: kernel(:)
        !! Array of coulomb kernel (either damped or undamped)
        logical, intent(in) :: do_V, do_E, do_grdE, do_HE
        !! Flags to enable/disable calculation of different electrostatic 
        !! properties
        real(rp), intent(inout) :: V, E(3), grdE(6), HE(10)
        !! Electric potential
        
        if(do_V) then
            V =  V + kernel(1) * q
        end if
        
        if(do_E) then
            E = E + q * kernel(2) * dr
        end if

        if(do_grdE) then
            ! xx
            grdE(1) = grdE(1) + q * 3.0_rp * kernel(3) * dr(1) * dr(1) - q * kernel(2)
            ! xy
            grdE(2) = grdE(2) + q * 3.0_rp * kernel(3) * dr(1) * dr(2)
            ! yy
            grdE(3) = grdE(3) + q * 3.0_rp * kernel(3) * dr(2) * dr(2) - q * kernel(2)
            ! xz
            grdE(4) = grdE(4) + q * 3.0_rp * kernel(3) * dr(1) * dr(3)
            ! yz
            grdE(5) = grdE(5) + q * 3.0_rp * kernel(3) * dr(2) * dr(3)
            ! zz
            grdE(6) = grdE(6) + q * 3.0_rp * kernel(3) * dr(3) * dr(3) - q * kernel(2)
        end if

        if(do_HE) then
            ! xxx
            HE(1) = HE(1) + (15.0_rp * kernel(4) * dr(1) * dr(1) * dr(1) - &
                             9.0_rp * kernel(3) * dr(1)) * q
            ! xxy
            HE(2) = HE(2) + (15.0_rp * kernel(4) * dr(1) * dr(1) * dr(2) - &
                             3.0_rp * kernel(3) * dr(2)) * q
            ! xxz
            HE(3) = HE(3) + (15.0_rp * kernel(4) * dr(1) * dr(1) * dr(3) - &
                             3.0_rp * kernel(3) * dr(3)) * q
            ! xyy
            HE(4) = HE(4) + (15.0_rp * kernel(4) * dr(2) * dr(2) * dr(1) - & 
                             3.0_rp * kernel(3) * dr(1)) * q 
            ! xyz
            HE(5) = HE(5) + 15.0_rp * kernel(4) * dr(1) * dr(2) * dr(3) * q
            
            ! xzz
            HE(6) = HE(6) + (15.0_rp * kernel(4) * dr(3) * dr(3) * dr(1) - &
                             3.0_rp * kernel(3) * dr(1)) * q
            ! yyy
            HE(7) = HE(7) + (15.0_rp * kernel(4) * dr(2) * dr(2) * dr(2) - & 
                             9.0_rp * kernel(3) * dr(2)) * q 
            ! yyz
            HE(8) = HE(8) + (15.0_rp * kernel(4) * dr(2) * dr(2) * dr(3) - &
                             3.0_rp * kernel(3) * dr(3)) * q
            ! yzz
            HE(9) = HE(9) + (15.0_rp * kernel(4) * dr(3) * dr(3) * dr(2) - &
                             3.0_rp * kernel(3) * dr(2)) * q
            ! zzz
            HE(10) = HE(10) + (15.0_rp * kernel(4) * dr(3) * dr(3) * dr(3) - &
                               9.0_rp * kernel(3) * dr(3)) * q
        end if

    end subroutine q_elec_prop