Takes in input a QM Helper object, with initialized atom types, and using a parameter file, it generates a OMMP System object that corresponds to the QM system. It is used for internal testing pourpose but other creative things are always possible.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(ommp_qm_helper), | intent(in) | :: | qmh | |||
character(len=*), | intent(in) | :: | prm_file | |||
type(ommp_system), | intent(inout), | pointer | :: | sys |
subroutine ommp_system_from_qm_helper(qmh, prm_file, sys)
!! Takes in input a QM Helper object, with initialized
!! atom types, and using a parameter file, it generates a
!! OMMP System object that corresponds to the QM system.
!! It is used for internal testing pourpose but other
!! creative things are always possible.
use mod_mmpol, only: mmpol_prepare, mmpol_init, mmpol_init_nonbonded, &
mmpol_init_bonded
use mod_topology, only: check_conn_matrix
use mod_adjacency_mat, only: build_conn_upto_n
use mod_prm, only: check_keyword, assign_pol, assign_mpoles, &
assign_vdw, assign_bond, assign_angle, assign_urey, &
assign_strbnd, assign_opb, assign_pitors, &
assign_torsion, assign_tortors, assign_angtor, &
assign_strtor, assign_imptorsion, get_prm_ff_type
use mod_constants, only: OMMP_STR_CHAR_MAX
use mod_io, only: fatal_error, large_file_read
use mod_utils, only: str_to_lower, str_uncomment
implicit none
type(ommp_system), intent(inout), pointer :: sys
type(ommp_qm_helper), intent(in) :: qmh
character(len=*), intent(in) :: prm_file
integer(ommp_integer) :: i, ist
character(len=OMMP_STR_CHAR_MAX), allocatable :: prm_buf(:)
if(.not. (qmh%qm_top%attype_initialized .and. &
qmh%qm_top%atz_initialized)) &
call ommp_fatal("In order to convert QM Helper to &
&OMMP System atom types and atomic &
&number should be set.")
allocate(sys)
! Load prm file in RAM
call large_file_read(prm_file, prm_buf)
! Remove comments from prm file
!$omp parallel do
do i=1, size(prm_buf)
prm_buf(i) = str_to_lower(prm_buf(i))
prm_buf(i) = str_uncomment(prm_buf(i), '!')
end do
call mmpol_init(sys, get_prm_ff_type(prm_buf), &
qmh%qm_top%mm_atoms, qmh%qm_top%mm_atoms)
do i=1, sys%top%mm_atoms
sys%eel%polar_mm(i) = i
end do
! Copy the topology from QM system!
sys%top%cmm = qmh%qm_top%cmm
sys%top%use_frozen = qmh%qm_top%use_frozen
sys%top%frozen = qmh%qm_top%frozen
sys%top%atz_initialized = qmh%qm_top%atz_initialized
sys%top%atz = qmh%qm_top%atz
sys%top%attype_initialized = qmh%qm_top%attype_initialized
sys%top%attype = qmh%qm_top%attype
! Create connectivities from adjacency matrix
call build_conn_upto_n(qmh%qm_top%conn(1), 4, sys%top%conn, .false.)
! Now assign parameters
if( .not. check_keyword(prm_buf)) then
call ommp_fatal("PRM file cannot be completely understood")
end if
call ommp_message("QMH->SYS Assigning electrostatic parameters", OMMP_VERBOSE_DEBUG)
call assign_pol(sys%eel, prm_buf)
call assign_mpoles(sys%eel, prm_buf)
call ommp_message("QMH->SYS Assigning non-bonded parameters", OMMP_VERBOSE_DEBUG)
call mmpol_init_nonbonded(sys)
call assign_vdw(sys%vdw, sys%top, prm_buf)
call ommp_message("QMH->SYS Assigning bonded parameters", OMMP_VERBOSE_DEBUG)
call mmpol_init_bonded(sys)
call check_conn_matrix(sys%top, 4)
call assign_bond(sys%bds, prm_buf)
call assign_angle(sys%bds, prm_buf)
call assign_urey(sys%bds, prm_buf)
call assign_strbnd(sys%bds, prm_buf)
call assign_opb(sys%bds, prm_buf)
call assign_pitors(sys%bds, prm_buf)
call assign_torsion(sys%bds, prm_buf)
call assign_imptorsion(sys%bds, prm_buf)
call assign_tortors(sys%bds, prm_buf)
call assign_angtor(sys%bds, prm_buf)
call assign_strtor(sys%bds, prm_buf)
deallocate(prm_buf)
call mmpol_prepare(sys)
call ommp_message('QMH->SYS Completed', OMMP_VERBOSE_DEBUG)
end subroutine ommp_system_from_qm_helper