elec_prop_M2M Subroutine

private subroutine elec_prop_M2M(eel, do_V, do_E, do_Egrd, do_EHes)

Computes the electric potential, field and field gradients of static multipoles at all sites (polarizable sites are a subset of static ones)

Arguments

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

Electrostatics data structure
logical, intent(in) :: do_V

Flags to enable/disable the calculation of different components
logical, intent(in) :: do_E

Flags to enable/disable the calculation of different components
logical, intent(in) :: do_Egrd

Flags to enable/disable the calculation of different components
logical, intent(in) :: do_EHes

Flags to enable/disable the calculation of different components

Calls

proc~~elec_prop_m2m~~CallsGraph proc~elec_prop_m2m elec_prop_M2M proc~preapare_fmm_static preapare_fmm_static proc~elec_prop_m2m->proc~preapare_fmm_static proc~cart_propfar_at_ipart cart_propfar_at_ipart proc~elec_prop_m2m->proc~cart_propfar_at_ipart proc~q_elec_prop q_elec_prop proc~elec_prop_m2m->proc~q_elec_prop proc~coulomb_kernel coulomb_kernel proc~elec_prop_m2m->proc~coulomb_kernel proc~mu_elec_prop mu_elec_prop proc~elec_prop_m2m->proc~mu_elec_prop proc~quad_elec_prop quad_elec_prop proc~elec_prop_m2m->proc~quad_elec_prop interface~mallocate mallocate proc~preapare_fmm_static->interface~mallocate proc~fmm_solve_for_multipoles fmm_solve_for_multipoles proc~preapare_fmm_static->proc~fmm_solve_for_multipoles interface~mfree mfree proc~preapare_fmm_static->interface~mfree proc~fmm_l2l fmm_l2l proc~cart_propfar_at_ipart->proc~fmm_l2l proc~ntot_sph_harm ntot_sph_harm proc~cart_propfar_at_ipart->proc~ntot_sph_harm proc~fatal_error fatal_error proc~coulomb_kernel->proc~fatal_error proc~i_alloc1 i_alloc1 interface~mallocate->proc~i_alloc1 proc~i_alloc3 i_alloc3 interface~mallocate->proc~i_alloc3 proc~l_alloc2 l_alloc2 interface~mallocate->proc~l_alloc2 proc~r_alloc1 r_alloc1 interface~mallocate->proc~r_alloc1 proc~r_alloc3 r_alloc3 interface~mallocate->proc~r_alloc3 proc~i_alloc2 i_alloc2 interface~mallocate->proc~i_alloc2 proc~r_alloc2 r_alloc2 interface~mallocate->proc~r_alloc2 proc~l_alloc1 l_alloc1 interface~mallocate->proc~l_alloc1 proc~make_vfact make_vfact proc~fmm_l2l->proc~make_vfact proc~fmm_l2l_rotation_work fmm_l2l_rotation_work proc~fmm_l2l->proc~fmm_l2l_rotation_work proc~fmm_solve_for_multipoles->proc~fatal_error proc~tree_p2m tree_p2m proc~fmm_solve_for_multipoles->proc~tree_p2m proc~tree_m2m tree_m2m proc~fmm_solve_for_multipoles->proc~tree_m2m proc~tree_l2l tree_l2l proc~fmm_solve_for_multipoles->proc~tree_l2l proc~tree_m2l tree_m2l proc~fmm_solve_for_multipoles->proc~tree_m2l proc~time_pull time_pull proc~fmm_solve_for_multipoles->proc~time_pull proc~time_push time_push proc~fmm_solve_for_multipoles->proc~time_push proc~i_free2 i_free2 interface~mfree->proc~i_free2 proc~r_free3 r_free3 interface~mfree->proc~r_free3 proc~r_free2 r_free2 interface~mfree->proc~r_free2 proc~l_free1 l_free1 interface~mfree->proc~l_free1 proc~l_free2 l_free2 interface~mfree->proc~l_free2 proc~i_free3 i_free3 interface~mfree->proc~i_free3 proc~i_free1 i_free1 interface~mfree->proc~i_free1 proc~r_free1 r_free1 interface~mfree->proc~r_free1 proc~ommp_message ommp_message proc~fatal_error->proc~ommp_message proc~close_output close_output proc~fatal_error->proc~close_output proc~chk_free chk_free proc~i_free2->proc~chk_free proc~chk_alloc chk_alloc proc~i_alloc1->proc~chk_alloc proc~memory_init memory_init proc~i_alloc1->proc~memory_init proc~i_alloc3->proc~chk_alloc proc~i_alloc3->proc~memory_init proc~l_alloc2->proc~chk_alloc proc~l_alloc2->proc~memory_init proc~carttosph carttosph proc~fmm_l2l_rotation_work->proc~carttosph proc~fmm_sph_rotate_oxz_work fmm_sph_rotate_oxz_work proc~fmm_l2l_rotation_work->proc~fmm_sph_rotate_oxz_work proc~trgev trgev proc~fmm_l2l_rotation_work->proc~trgev proc~fmm_l2l_ztranslate_work fmm_l2l_ztranslate_work proc~fmm_l2l_rotation_work->proc~fmm_l2l_ztranslate_work proc~fmm_sph_rotate_oz_adj_work fmm_sph_rotate_oz_adj_work proc~fmm_l2l_rotation_work->proc~fmm_sph_rotate_oz_adj_work proc~fmm_sph_rotate_oz_work fmm_sph_rotate_oz_work proc~fmm_l2l_rotation_work->proc~fmm_sph_rotate_oz_work proc~tree_p2m->proc~ntot_sph_harm proc~fmm_m2m fmm_m2m proc~tree_p2m->proc~fmm_m2m proc~fmm_error fmm_error proc~tree_p2m->proc~fmm_error proc~r_free3->proc~chk_free proc~r_free2->proc~chk_free proc~l_free1->proc~chk_free proc~l_free2->proc~chk_free proc~close_output->proc~ommp_message proc~i_free3->proc~chk_free proc~tree_m2m->proc~ntot_sph_harm proc~tree_m2m->proc~fmm_m2m proc~i_free1->proc~chk_free proc~tree_l2l->proc~fmm_l2l proc~tree_l2l->proc~ntot_sph_harm proc~tree_m2l->proc~ntot_sph_harm proc~fmm_m2l fmm_m2l proc~tree_m2l->proc~fmm_m2l proc~time_pull->proc~fatal_error proc~time_pull->proc~ommp_message proc~mem_stat mem_stat proc~time_pull->proc~mem_stat proc~r_free1->proc~chk_free proc~r_alloc1->proc~chk_alloc proc~r_alloc1->proc~memory_init proc~r_alloc3->proc~chk_alloc proc~r_alloc3->proc~memory_init proc~i_alloc2->proc~chk_alloc proc~i_alloc2->proc~memory_init proc~r_alloc2->proc~chk_alloc proc~r_alloc2->proc~memory_init proc~l_alloc1->proc~chk_alloc proc~l_alloc1->proc~memory_init proc~time_push->proc~fatal_error proc~time_push->proc~mem_stat proc~mem_stat->proc~memory_init proc~fmm_m2m_rotation_work fmm_m2m_rotation_work proc~fmm_m2m->proc~fmm_m2m_rotation_work proc~chk_free->proc~fatal_error proc~chk_alloc->proc~fatal_error proc~fmm_m2l_rotation_work fmm_m2l_rotation_work proc~fmm_m2l->proc~fmm_m2l_rotation_work proc~fmm_m2m_rotation_work->proc~carttosph proc~fmm_m2m_rotation_work->proc~fmm_sph_rotate_oxz_work proc~fmm_m2m_rotation_work->proc~trgev proc~fmm_m2m_rotation_work->proc~fmm_sph_rotate_oz_adj_work proc~fmm_m2m_rotation_work->proc~fmm_sph_rotate_oz_work proc~fmm_m2m_ztranslate_work fmm_m2m_ztranslate_work proc~fmm_m2m_rotation_work->proc~fmm_m2m_ztranslate_work proc~fmm_m2l_rotation_work->proc~carttosph proc~fmm_m2l_rotation_work->proc~fmm_sph_rotate_oxz_work proc~fmm_m2l_rotation_work->proc~trgev proc~fmm_m2l_rotation_work->proc~fmm_sph_rotate_oz_adj_work proc~fmm_m2l_rotation_work->proc~fmm_sph_rotate_oz_work proc~fmm_m2l_ztranslate_work fmm_m2l_ztranslate_work proc~fmm_m2l_rotation_work->proc~fmm_m2l_ztranslate_work

Called by

proc~~elec_prop_m2m~~CalledByGraph proc~elec_prop_m2m elec_prop_M2M proc~prepare_fixedelec prepare_fixedelec proc~prepare_fixedelec->proc~elec_prop_m2m proc~fixedelec_geomgrad fixedelec_geomgrad proc~fixedelec_geomgrad->proc~prepare_fixedelec proc~energy_mm_mm energy_MM_MM proc~energy_mm_mm->proc~prepare_fixedelec proc~ommp_fixedelec_geomgrad ommp_fixedelec_geomgrad proc~ommp_fixedelec_geomgrad->proc~fixedelec_geomgrad proc~ommp_get_fixedelec_energy ommp_get_fixedelec_energy proc~ommp_get_fixedelec_energy->proc~energy_mm_mm proc~ommp_full_geomgrad ommp_full_geomgrad proc~ommp_full_geomgrad->proc~fixedelec_geomgrad proc~c_ommp_fixedelec_geomgrad C_ommp_fixedelec_geomgrad proc~c_ommp_fixedelec_geomgrad->proc~ommp_fixedelec_geomgrad proc~ommp_get_full_ele_energy ommp_get_full_ele_energy proc~ommp_get_full_ele_energy->proc~ommp_get_fixedelec_energy proc~c_ommp_get_fixedelec_energy C_ommp_get_fixedelec_energy proc~c_ommp_get_fixedelec_energy->proc~ommp_get_fixedelec_energy proc~c_ommp_full_geomgrad C_ommp_full_geomgrad proc~c_ommp_full_geomgrad->proc~ommp_full_geomgrad proc~ommp_get_full_energy ommp_get_full_energy proc~ommp_get_full_energy->proc~ommp_get_full_ele_energy 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_full_energy C_ommp_get_full_energy proc~c_ommp_get_full_energy->proc~ommp_get_full_energy

Contents

Source Code

elec_prop_M2M

Source Code

    subroutine elec_prop_M2M(eel, do_V, do_E, do_Egrd, do_EHes)
        !! Computes the electric potential, field and field gradients of 
        !! static multipoles at all sites (polarizable sites are a 
        !! subset of static ones)
        implicit none
        
        type(ommp_electrostatics_type), intent(inout) :: eel
        !! Electrostatics data structure
        logical, intent(in) :: do_V, do_E, do_Egrd, do_EHes
        !! Flags to enable/disable the calculation of different components


        real(rp) :: kernel(6), dr(3), tmpV, tmpE(3), tmpEgr(6), tmpHE(10), scalf
        integer(ip) :: i, j, idx, sidx, ikernel
        logical :: to_do, to_scale
        type(ommp_topology_type), pointer :: top

        top => eel%top

        if(do_EHes) then
            ikernel = 3 
        elseif(do_Egrd) then
            ikernel = 2
        elseif(do_E) then
            ikernel = 1
        elseif(do_V) then
            ikernel = 0
        else
            return
        end if
        if(eel%amoeba) ikernel = ikernel + 2

        if(eel%use_fmm) then
            call preapare_fmm_static(eel)

            ! TODO !
            !$omp parallel do default(shared)
            do i=1, top%mm_atoms 
                call cart_propfar_at_ipart(eel%fmm_static, i, &
                                        do_V, eel%V_M2M(i), &
                                        do_E, eel%E_M2M(:,i), &
                                        do_Egrd, eel%Egrd_M2M(:,i), &
                                        do_EHes, eel%EHes_M2M(:,i))
            end do

            if(allocated(eel%list_S_S_fmm_far)) then
                !$omp parallel do default(shared) schedule(dynamic) &
                !$omp private(i,j,idx,scalf,dr,kernel,tmpV,tmpE,tmpEgr,tmpHE)
                do i=1, top%mm_atoms
                    do idx=eel%list_S_S_fmm_far%ri(i), eel%list_S_S_fmm_far%ri(i+1)-1
                        j = eel%list_S_S_fmm_far%ci(idx)
                        scalf = eel%scalef_S_S_fmm_far(idx) - 1.0

                        dr = top%cmm(:,i) - top%cmm(:, j)
                        call coulomb_kernel(dr, ikernel, kernel)
                        
                        if(do_V) tmpV = 0.0_rp
                        if(do_E) tmpE = 0.0_rp
                        if(do_Egrd) tmpEgr = 0.0_rp
                        if(do_EHes) tmpHE = 0.0_rp

                        call q_elec_prop(eel%q(1,j), dr, kernel, &
                                            do_V, tmpV, & 
                                            do_E, tmpE, &
                                            do_Egrd, tmpEgr, &
                                            do_EHes, tmpHE)
                        if(eel%amoeba) then
                            call mu_elec_prop(eel%q(2:4,j), dr, kernel, &
                                                do_V, tmpV, & 
                                                do_E, tmpE, &
                                                do_Egrd, tmpEgr, &
                                                do_EHes, tmpHE)

                            call quad_elec_prop(eel%q(5:10,j), dr, kernel, &
                                                do_V, tmpV, & 
                                                do_E, tmpE, &
                                                do_Egrd, tmpEgr, &
                                                do_EHes, tmpHE)
                        end if

                        if(do_V) eel%V_M2M(i) = eel%V_M2M(i) + tmpV * scalf
                        if(do_E) eel%E_M2M(:,i) = eel%E_M2M(:,i) + tmpE * scalf
                        if(do_Egrd) eel%Egrd_M2M(:,i) = eel%Egrd_M2M(:,i) + tmpEgr * scalf
                        if(do_EHes) eel%EHes_M2M(:,i) = eel%EHes_M2M(:,i) + tmpHE * scalf
                    end do
                end do
            end if
            
            !$omp parallel do default(shared) schedule(dynamic) &
            !$omp private(i,j,idx,sidx,to_scale,to_do,scalf,dr,kernel,tmpV,tmpE,tmpEgr,tmpHE)
            do i=1, top%mm_atoms
                do idx=eel%fmm_near_field_list%ri(i), &
                       eel%fmm_near_field_list%ri(i+1)-1
                    j = eel%fmm_near_field_list%ci(idx)
                    
                    to_do = .true.
                    to_scale = .false.
                    scalf = 1.0

                    ! Check if the element should be scaled
                    do sidx=eel%list_S_S%ri(i), eel%list_S_S%ri(i+1)-1
                        if(eel%list_S_S%ci(sidx) == j) then
                            to_scale = .true.
                            exit
                        end if
                    end do

                    !If it should set the correct variables
                    if(to_scale) then
                        to_do = eel%todo_S_S(sidx)
                        scalf = eel%scalef_S_S(sidx)
                    end if

                    if(to_do) then
                        dr = top%cmm(:,i) - top%cmm(:, j)
                        call coulomb_kernel(dr, ikernel, kernel)
                        
                        if(do_V) tmpV = 0.0_rp
                        if(do_E) tmpE = 0.0_rp
                        if(do_Egrd) tmpEgr = 0.0_rp
                        if(do_EHes) tmpHE = 0.0_rp

                        call q_elec_prop(eel%q(1,j), dr, kernel, &
                                            do_V, tmpV, & 
                                            do_E, tmpE, &
                                            do_Egrd, tmpEgr, &
                                            do_EHes, tmpHE)
                        if(eel%amoeba) then
                            call mu_elec_prop(eel%q(2:4,j), dr, kernel, &
                                                do_V, tmpV, & 
                                                do_E, tmpE, &
                                                do_Egrd, tmpEgr, &
                                                do_EHes, tmpHE)

                            call quad_elec_prop(eel%q(5:10,j), dr, kernel, &
                                                do_V, tmpV, & 
                                                do_E, tmpE, &
                                                do_Egrd, tmpEgr, &
                                                do_EHes, tmpHE)
                        end if

                        if(do_V) eel%V_M2M(i) = eel%V_M2M(i) + tmpV * scalf
                        if(do_E) eel%E_M2M(:,i) = eel%E_M2M(:,i) + tmpE * scalf
                        if(do_Egrd) eel%Egrd_M2M(:,i) = eel%Egrd_M2M(:,i) + tmpEgr * scalf
                        if(do_EHes) eel%EHes_M2M(:,i) = eel%EHes_M2M(:,i) + tmpHE * scalf
                    end if
                end do
            end do
        else
        if(eel%amoeba) then
            !$omp parallel do default(shared) schedule(dynamic) &
            !$omp private(i,j,idx,to_do,to_scale,scalf,dr,kernel,tmpV,tmpE,tmpEgr,tmpHE)
            do j=1, top%mm_atoms
                ! loop on sources
                do i=1, top%mm_atoms
                    if(j == i) cycle
                    !loop on target

                    to_do = .true.
                    to_scale = .false.
                    scalf = 1.0

                    ! Check if the element should be scaled
                    do idx=eel%list_S_S%ri(i), eel%list_S_S%ri(i+1)-1
                        if(eel%list_S_S%ci(idx) == j) then
                            to_scale = .true.
                            exit
                        end if
                    end do

                    !If it should set the correct variables
                    if(to_scale) then
                        to_do = eel%todo_S_S(idx)
                        scalf = eel%scalef_S_S(idx)
                    end if

                    if(to_do) then
                        dr = top%cmm(:,j) - top%cmm(:, i)
                        call coulomb_kernel(dr, ikernel, kernel)
                        
                        if(do_V) tmpV = 0.0_rp
                        if(do_E) tmpE = 0.0_rp
                        if(do_Egrd) tmpEgr = 0.0_rp
                        if(do_EHes) tmpHE = 0.0_rp

                        call q_elec_prop(eel%q(1,i), dr, kernel, &
                                         do_V, tmpV, & 
                                         do_E, tmpE, &
                                         do_Egrd, tmpEgr, &
                                         do_EHes, tmpHE)

                        call mu_elec_prop(eel%q(2:4,i), dr, kernel, &
                                          do_V, tmpV, & 
                                          do_E, tmpE, &
                                          do_Egrd, tmpEgr, &
                                          do_EHes, tmpHE)

                        call quad_elec_prop(eel%q(5:10,i), dr, kernel, &
                                            do_V, tmpV, & 
                                            do_E, tmpE, &
                                            do_Egrd, tmpEgr, &
                                            do_EHes, tmpHE)

                        if(to_scale) then
                            if(do_V) eel%V_M2M(j) = eel%V_M2M(j) + tmpV * scalf
                            if(do_E) eel%E_M2M(:,j) = eel%E_M2M(:,j) + tmpE * scalf
                            if(do_Egrd) eel%Egrd_M2M(:,j) = eel%Egrd_M2M(:,j) + tmpEgr * scalf
                            if(do_EHes) eel%EHes_M2M(:,j) = eel%EHes_M2M(:,j) + tmpHE * scalf
                        else
                            if(do_V) eel%V_M2M(j) = eel%V_M2M(j) + tmpV
                            if(do_E) eel%E_M2M(:,j) = eel%E_M2M(:,j) + tmpE
                            if(do_Egrd) eel%Egrd_M2M(:,j) = eel%Egrd_M2M(:,j) + tmpEgr
                            if(do_EHes) eel%EHes_M2M(:,j) = eel%EHes_M2M(:,j) + tmpHE
                        end if
                    end if
                end do
            end do
        else
            !$omp parallel do default(shared) schedule(dynamic) &
            !$omp private(i,j,idx,to_do,to_scale,scalf,dr,kernel,tmpV,tmpE,tmpEgr,tmpHE) 
            do j=1, top%mm_atoms
                ! loop on sources
                do i=1, top%mm_atoms
                    if(j == i) cycle
                    !loop on target
                    to_do = .true.
                    to_scale = .false.
                    scalf = 1.0

                    ! Check if the element should be scaled
                    do idx=eel%list_S_S%ri(i), eel%list_S_S%ri(i+1)-1
                        if(eel%list_S_S%ci(idx) == j) then
                            to_scale = .true.
                            exit
                        end if
                    end do

                    !If it should set the correct variables
                    if(to_scale) then
                        to_do = eel%todo_S_S(idx)
                        scalf = eel%scalef_S_S(idx)
                    end if
                    
                    if(to_do) then
                        dr = top%cmm(:,j) - top%cmm(:, i)
                        call coulomb_kernel(dr, ikernel, kernel)
                        
                        if(do_V) tmpV = 0.0_rp
                        if(do_E) tmpE = 0.0_rp
                        if(do_Egrd) tmpEgr = 0.0_rp
                        if(do_EHes) tmpHE = 0.0_rp

                        call q_elec_prop(eel%q(1,i), dr, kernel, &
                                         do_V, tmpV, & 
                                         do_E, tmpE, &
                                         do_Egrd, tmpEgr, &
                                         do_EHes, tmpHE)

                        if(to_scale) then
                            if(do_V) eel%V_M2M(j) = eel%V_M2M(j) + tmpV * scalf
                            if(do_E) eel%E_M2M(:,j) = eel%E_M2M(:,j) + tmpE * scalf
                            if(do_Egrd) eel%Egrd_M2M(:,j) = eel%Egrd_M2M(:,j) + tmpEgr * scalf
                            if(do_EHes) eel%EHes_M2M(:,j) = eel%EHes_M2M(:,j) + tmpHE * scalf
                        else
                            if(do_V) eel%V_M2M(j) = eel%V_M2M(j) + tmpV
                            if(do_E) eel%E_M2M(:,j) = eel%E_M2M(:,j) + tmpE
                            if(do_Egrd) eel%Egrd_M2M(:,j) = eel%Egrd_M2M(:,j) + tmpEgr
                            if(do_EHes) eel%EHes_M2M(:,j) = eel%EHes_M2M(:,j) + tmpHE
                        end if
                    end if
                end do
            end do
        end if
        end if
    end subroutine elec_prop_M2M