Prints a complete summary of all the quantities stored in the MMPol module
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(ommp_system), | intent(in) | :: | sys_obj | |||
character(len=*), | intent(in), | optional | :: | of_name |
subroutine mmpol_ommp_print_summary(sys_obj, of_name)
!! Prints a complete summary of all the quantities stored
!! in the MMPol module
use mod_memory, only: mallocate, mfree
use mod_io, only: iof_mmpol, print_matrix, print_int_vec
use mod_utils, only: sort_ivec
implicit none
type(ommp_system), intent(in) :: sys_obj
character(len=*), intent(in), optional :: of_name
integer(ip) :: of_unit
integer(ip) :: i, j, grp, igrp, lst(1000), ilst
real(rp), allocatable :: polar(:) ! Polarizabilities of all atoms
integer(ip), allocatable :: tmp(:)
character(len=OMMP_STR_CHAR_MAX) :: str
if(present(of_name)) then
of_unit = 101
open(unit=of_unit, &
file=of_name(1:len(trim(of_name))), &
action='write')
else
of_unit = iof_mmpol
end if
call mallocate('mmpol_ommp_print_summary [polar]', &
sys_obj%top%mm_atoms, polar)
polar = 0.0_rp
do i=1, sys_obj%eel%pol_atoms
polar(sys_obj%eel%polar_mm(i)) = sys_obj%eel%pol(i)
end do
write(of_unit, '(A, 4F8.4)') 'mscale: ', sys_obj%eel%mscale
if(sys_obj%eel%pol_atoms > 0) then
write(of_unit, '(A, 4F8.4)') 'pscale: ', sys_obj%eel%pscale
if(sys_obj%amoeba) write(of_unit, '(A, 4F8.4)') 'pscale (intra): ', &
sys_obj%eel%pscale_intra
write(of_unit, '(A, 4F8.4)') 'dscale: ', sys_obj%eel%dscale
write(of_unit, '(A, 4F8.4)') 'uscale: ', sys_obj%eel%uscale
end if
call print_matrix(.true., 'coordinates:', sys_obj%top%cmm, of_unit)
if(sys_obj%amoeba) then
call print_matrix(.true., 'multipoles - non rotated:', &
sys_obj%eel%q0, of_unit)
end if
call print_matrix(.true., 'multipoles :', sys_obj%eel%q, of_unit)
call print_matrix(.true., 'coordinates of polarizable atoms:', &
sys_obj%eel%cpol, of_unit)
call print_matrix(.false., 'polarizabilities:', polar, of_unit)
call print_matrix(.false., 'thole factors:', sys_obj%eel%thole, of_unit)
call print_int_vec('mm_polar list:', sys_obj%eel%mm_polar, &
of_unit)
call print_int_vec('polar_mm list:', sys_obj%eel%polar_mm, &
of_unit)
! write the connectivity information for each atom:
1000 format(t3,'connectivity information for the ',i8,'-th atom:')
do i = 1, sys_obj%top%mm_atoms
write(of_unit, 1000) i
do j=1, 4
if(j == 4 .and. .not. sys_obj%amoeba) cycle
write(str, "('1-', I1, ' neighors:')") j+1
call sort_ivec(sys_obj%top%conn(j)%ci(&
sys_obj%top%conn(j)%ri(i):&
sys_obj%top%conn(j)%ri(i+1)-1), tmp)
call print_int_vec(trim(str), tmp, of_unit)
end do
if(sys_obj%amoeba) then
do j=1, 4
ilst = 1
do igrp = &
sys_obj%eel%pg_conn(j)%ri(sys_obj%eel%mmat_polgrp(i)), &
sys_obj%eel%pg_conn(j)%ri(sys_obj%eel%mmat_polgrp(i)+1)-1
grp = sys_obj%eel%pg_conn(j)%ci(igrp)
lst(ilst:ilst+sys_obj%eel%polgrp_mmat%ri(grp+1) - &
sys_obj%eel%polgrp_mmat%ri(grp)-1) = &
sys_obj%eel%polgrp_mmat%ci(&
sys_obj%eel%polgrp_mmat%ri(grp):&
sys_obj%eel%polgrp_mmat%ri(grp+1)-1)
ilst = ilst+&
sys_obj%eel%polgrp_mmat%ri(grp+1)- &
sys_obj%eel%polgrp_mmat%ri(grp)
end do
write(str, "('1-', I1, ' polarization neighors:')") j
if(ilst == 1) ilst = 0
call sort_ivec(lst(1:ilst-1), tmp)
! needed to addres the empty array case
call print_int_vec(trim(str), tmp, of_unit)
end do
end if
end do
if(present(of_name)) close(of_unit)
if(allocated(tmp)) &
call mfree('mmpol_ommp_print_summary [tmp]', tmp)
call mfree('mmpol_ommp_print_summary [polar]', polar)
end subroutine mmpol_ommp_print_summary