Compute some derived quantities from the input that
are used during the calculation. The upstream code have
to provide cmm, q, pol, adjacency matrix and in
the case of AMOEBA also multipoles rotation information, and
polarization group information.
This routine
* compute connectivity lists from connected atoms
* invert polar_mm list creating mm_polar
* populate cpol list of coordinates
* compute factors for thole damping
* scales by 1/3 AMOEBA quadrupoles (?)
* Build list for polarization groups, compute groups connectivity
* performs multipoles rotation
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(ommp_system), | intent(inout) | :: | sys_obj |
The system object to bi initialized |
subroutine mmpol_prepare(sys_obj)
!! Compute some derived quantities from the input that
!! are used during the calculation. The upstream code have
!! to provide cmm, q, pol, adjacency matrix and in
!! the case of AMOEBA also multipoles rotation information, and
!! polarization group information.
!! This routine
!! * compute connectivity lists from connected atoms
!! * invert polar_mm list creating mm_polar
!! * populate cpol list of coordinates
!! * compute factors for thole damping
!! * scales by 1/3 AMOEBA quadrupoles (?)
!! * Build list for polarization groups, compute groups connectivity
!! * performs multipoles rotation
use mod_adjacency_mat, only: build_conn_upto_n, matcpy, reverse_grp_tab
use mod_io, only: ommp_message
use mod_profiling, only: time_push, time_pull
use mod_constants, only: OMMP_VERBOSE_DEBUG
use mod_electrostatics, only: thole_init, remove_null_pol, &
make_screening_lists
implicit none
type(ommp_system), intent(inout) :: sys_obj
!! The system object to bi initialized
integer(ip) :: i
type(yale_sparse) :: adj, pg_adj
call time_push()
call ommp_message("Building connectivity lists", OMMP_VERBOSE_DEBUG)
! compute connectivity lists from connected atoms
if(size(sys_obj%top%conn) < 4) then
call matcpy(sys_obj%top%conn(1), adj)
deallocate(sys_obj%top%conn)
call build_conn_upto_n(adj, 4, sys_obj%top%conn, .false.)
end if
call remove_null_pol(sys_obj%eel)
! invert mm_polar list creating mm_polar
sys_obj%eel%mm_polar(:) = 0
if(sys_obj%eel%pol_atoms > 0) then
call ommp_message("Creating MM->polar and polar->MM lists", &
OMMP_VERBOSE_DEBUG)
!$omp parallel do default(shared) private(i)
do i = 1, sys_obj%eel%pol_atoms
sys_obj%eel%mm_polar(sys_obj%eel%polar_mm(i)) = i
end do
call ommp_message("Populating coordinates of polarizable atoms", &
OMMP_VERBOSE_DEBUG)
! populate cpol list of coordinates
!$omp parallel do default(shared) private(i)
do i = 1, sys_obj%eel%pol_atoms
sys_obj%eel%cpol(:,i) = sys_obj%top%cmm(:, sys_obj%eel%polar_mm(i))
end do
call ommp_message("Setting Thole factors", OMMP_VERBOSE_DEBUG)
! compute factors for thole damping
call thole_init(sys_obj%eel)
else
sys_obj%eel%thole = 0.0
end if
if(sys_obj%amoeba) then
! Copy multipoles from q to q0
sys_obj%eel%q0 = sys_obj%eel%q
! scales by 1/3 AMOEBA quadrupoles (?)
! Mysterious division of multipoles by three
! FL told me that it was done like that in
! Tinker
sys_obj%eel%q0(5:10,:) = sys_obj%eel%q0(5:10,:) / 3.0_rp
! polarization groups connectivity list
call reverse_grp_tab(sys_obj%eel%mmat_polgrp, &
sys_obj%eel%polgrp_mmat)
call build_pg_adjacency_matrix(sys_obj%eel, pg_adj)
call build_conn_upto_n(pg_adj, 3, sys_obj%eel%pg_conn, .true.)
! performs multipoles rotation
call rotate_multipoles(sys_obj%eel)
end if
call ommp_message("Building screening lists", OMMP_VERBOSE_DEBUG)
call time_push()
call make_screening_lists(sys_obj%eel)
call time_pull("Preparing screening lists")
call ommp_message("MMPol initialization (mmpol_prepare) completed.", OMMP_VERBOSE_DEBUG)
call time_pull('MMPol object initialization (mmpol_prepare)')
end subroutine mmpol_prepare