mmpol_init_from_mmp Subroutine

public subroutine mmpol_init_from_mmp(input_file, sys_obj)

Uses

  • proc~~mmpol_init_from_mmp~~UsesGraph proc~mmpol_init_from_mmp mmpol_init_from_mmp module~mod_electrostatics mod_electrostatics proc~mmpol_init_from_mmp->module~mod_electrostatics module~mod_adjacency_mat mod_adjacency_mat proc~mmpol_init_from_mmp->module~mod_adjacency_mat module~mod_memory mod_memory proc~mmpol_init_from_mmp->module~mod_memory module~mod_io mod_io proc~mmpol_init_from_mmp->module~mod_io module~mod_constants mod_constants proc~mmpol_init_from_mmp->module~mod_constants module~mod_mmpol mod_mmpol proc~mmpol_init_from_mmp->module~mod_mmpol module~mod_utils mod_utils proc~mmpol_init_from_mmp->module~mod_utils module~mod_electrostatics->module~mod_adjacency_mat module~mod_electrostatics->module~mod_memory module~mod_electrostatics->module~mod_io module~mod_electrostatics->module~mod_constants module~mod_topology mod_topology module~mod_electrostatics->module~mod_topology module~mod_profiling mod_profiling module~mod_electrostatics->module~mod_profiling module~mod_adjacency_mat->module~mod_memory module~mod_memory->module~mod_io module~mod_memory->module~mod_constants iso_c_binding iso_c_binding module~mod_memory->iso_c_binding module~mod_io->module~mod_constants module~mod_constants->iso_c_binding module~mod_mmpol->module~mod_electrostatics module~mod_mmpol->module~mod_adjacency_mat module~mod_mmpol->module~mod_memory module~mod_mmpol->module~mod_io module~mod_mmpol->module~mod_constants module~mod_mmpol->module~mod_topology module~mod_nonbonded mod_nonbonded module~mod_mmpol->module~mod_nonbonded module~mod_bonded mod_bonded module~mod_mmpol->module~mod_bonded module~mod_link_atom mod_link_atom module~mod_mmpol->module~mod_link_atom module~mod_utils->module~mod_memory module~mod_utils->module~mod_constants module~mod_topology->module~mod_adjacency_mat module~mod_topology->module~mod_memory module~mod_profiling->module~mod_memory module~mod_profiling->module~mod_io module~mod_profiling->module~mod_constants module~mod_nonbonded->module~mod_adjacency_mat module~mod_nonbonded->module~mod_memory module~mod_nonbonded->module~mod_constants module~mod_nonbonded->module~mod_topology module~mod_neighbor_list mod_neighbor_list module~mod_nonbonded->module~mod_neighbor_list module~mod_bonded->module~mod_memory module~mod_bonded->module~mod_io module~mod_bonded->module~mod_topology module~mod_link_atom->module~mod_memory module~mod_link_atom->module~mod_io module~mod_link_atom->module~mod_constants module~mod_link_atom->module~mod_utils module~mod_link_atom->module~mod_topology module~mod_link_atom->module~mod_nonbonded module~mod_link_atom->module~mod_bonded module~mod_neighbor_list->module~mod_adjacency_mat module~mod_neighbor_list->module~mod_memory module~mod_neighbor_list->module~mod_io

This function read a .mmp file (revision 2 and 3) are supported and initialize all the quantities need to describe the environment within this library.

Arguments

Type IntentOptional Attributes Name
character(len=*), intent(in) :: input_file

name of the input file
type(ommp_system), intent(inout) :: sys_obj

System data structure to be initialized

Calls

proc~~mmpol_init_from_mmp~~CallsGraph proc~mmpol_init_from_mmp mmpol_init_from_mmp proc~adj_mat_from_conn adj_mat_from_conn proc~mmpol_init_from_mmp->proc~adj_mat_from_conn proc~mmpol_prepare mmpol_prepare proc~mmpol_init_from_mmp->proc~mmpol_prepare proc~fatal_error fatal_error proc~mmpol_init_from_mmp->proc~fatal_error proc~time_push time_push proc~mmpol_init_from_mmp->proc~time_push proc~skip_lines skip_lines proc~mmpol_init_from_mmp->proc~skip_lines proc~ommp_message ommp_message proc~mmpol_init_from_mmp->proc~ommp_message proc~memory_init memory_init proc~mmpol_init_from_mmp->proc~memory_init interface~mallocate mallocate proc~mmpol_init_from_mmp->interface~mallocate proc~mmpol_init mmpol_init proc~mmpol_init_from_mmp->proc~mmpol_init proc~set_screening_parameters set_screening_parameters proc~mmpol_init_from_mmp->proc~set_screening_parameters interface~mfree mfree proc~mmpol_init_from_mmp->interface~mfree proc~time_pull time_pull proc~mmpol_init_from_mmp->proc~time_pull proc~polgroup11_to_mm2pg polgroup11_to_mm2pg proc~mmpol_init_from_mmp->proc~polgroup11_to_mm2pg proc~adj_mat_from_conn->interface~mallocate proc~adj_mat_from_conn->interface~mfree proc~compress_list compress_list proc~adj_mat_from_conn->proc~compress_list proc~sort_ivec_inplace sort_ivec_inplace proc~adj_mat_from_conn->proc~sort_ivec_inplace proc~mmpol_prepare->proc~time_push proc~mmpol_prepare->proc~ommp_message proc~mmpol_prepare->proc~time_pull proc~build_conn_upto_n build_conn_upto_n proc~mmpol_prepare->proc~build_conn_upto_n proc~remove_null_pol remove_null_pol proc~mmpol_prepare->proc~remove_null_pol proc~reverse_grp_tab reverse_grp_tab proc~mmpol_prepare->proc~reverse_grp_tab proc~make_screening_lists make_screening_lists proc~mmpol_prepare->proc~make_screening_lists proc~matcpy matcpy proc~mmpol_prepare->proc~matcpy proc~thole_init thole_init proc~mmpol_prepare->proc~thole_init proc~build_pg_adjacency_matrix build_pg_adjacency_matrix proc~mmpol_prepare->proc~build_pg_adjacency_matrix proc~rotate_multipoles rotate_multipoles proc~mmpol_prepare->proc~rotate_multipoles proc~fatal_error->proc~ommp_message proc~close_output close_output proc~fatal_error->proc~close_output proc~time_push->proc~fatal_error proc~mem_stat mem_stat proc~time_push->proc~mem_stat proc~r_alloc1 r_alloc1 interface~mallocate->proc~r_alloc1 proc~i_alloc3 i_alloc3 interface~mallocate->proc~i_alloc3 proc~r_alloc3 r_alloc3 interface~mallocate->proc~r_alloc3 proc~i_alloc2 i_alloc2 interface~mallocate->proc~i_alloc2 proc~l_alloc2 l_alloc2 interface~mallocate->proc~l_alloc2 proc~i_alloc1 i_alloc1 interface~mallocate->proc~i_alloc1 proc~l_alloc1 l_alloc1 interface~mallocate->proc~l_alloc1 proc~r_alloc2 r_alloc2 interface~mallocate->proc~r_alloc2 proc~mmpol_init->proc~fatal_error proc~mmpol_init->proc~time_push proc~mmpol_init->proc~time_pull proc~topology_init topology_init proc~mmpol_init->proc~topology_init proc~electrostatics_init electrostatics_init proc~mmpol_init->proc~electrostatics_init proc~set_screening_parameters->proc~fatal_error proc~r_free1 r_free1 interface~mfree->proc~r_free1 proc~i_free3 i_free3 interface~mfree->proc~i_free3 proc~l_free2 l_free2 interface~mfree->proc~l_free2 proc~i_free1 i_free1 interface~mfree->proc~i_free1 proc~l_free1 l_free1 interface~mfree->proc~l_free1 proc~r_free3 r_free3 interface~mfree->proc~r_free3 proc~r_free2 r_free2 interface~mfree->proc~r_free2 proc~i_free2 i_free2 interface~mfree->proc~i_free2 proc~time_pull->proc~fatal_error proc~time_pull->proc~ommp_message proc~time_pull->proc~mem_stat proc~polgroup11_to_mm2pg->proc~fatal_error proc~chk_free chk_free proc~r_free1->proc~chk_free proc~r_alloc1->proc~memory_init proc~chk_alloc chk_alloc proc~r_alloc1->proc~chk_alloc proc~i_alloc3->proc~memory_init proc~i_alloc3->proc~chk_alloc proc~r_alloc3->proc~memory_init proc~r_alloc3->proc~chk_alloc proc~build_conn_upto_n->proc~matcpy proc~mat_mult mat_mult proc~build_conn_upto_n->proc~mat_mult proc~matfree matfree proc~build_conn_upto_n->proc~matfree proc~sparse_identity sparse_identity proc~build_conn_upto_n->proc~sparse_identity proc~mat_andnot mat_andnot proc~build_conn_upto_n->proc~mat_andnot proc~remove_null_pol->proc~ommp_message proc~remove_null_pol->interface~mallocate proc~remove_null_pol->interface~mfree proc~compress_list->interface~mallocate proc~compress_list->interface~mfree proc~reverse_grp_tab->interface~mallocate proc~reverse_grp_tab->interface~mfree proc~reverse_grp_tab->proc~compress_list proc~make_screening_lists->interface~mallocate proc~make_screening_lists->interface~mfree proc~make_screening_lists->proc~compress_list proc~compress_data compress_data proc~make_screening_lists->proc~compress_data proc~screening_rules screening_rules proc~make_screening_lists->proc~screening_rules proc~i_alloc2->proc~memory_init proc~i_alloc2->proc~chk_alloc proc~l_alloc2->proc~memory_init proc~l_alloc2->proc~chk_alloc proc~i_free3->proc~chk_free proc~l_free2->proc~chk_free proc~i_free1->proc~chk_free proc~topology_init->interface~mallocate proc~mem_stat->proc~memory_init proc~i_alloc1->proc~memory_init proc~i_alloc1->proc~chk_alloc proc~l_alloc1->proc~memory_init proc~l_alloc1->proc~chk_alloc proc~r_alloc2->proc~memory_init proc~r_alloc2->proc~chk_alloc proc~sort_ivec_inplace->interface~mfree proc~sort_ivec sort_ivec proc~sort_ivec_inplace->proc~sort_ivec proc~thole_init->proc~fatal_error proc~thole_init->proc~ommp_message proc~reallocate_mat reallocate_mat proc~build_pg_adjacency_matrix->proc~reallocate_mat proc~rotation_matrix rotation_matrix proc~rotate_multipoles->proc~rotation_matrix proc~close_output->proc~ommp_message proc~l_free1->proc~chk_free proc~r_free3->proc~chk_free proc~r_free2->proc~chk_free proc~electrostatics_init->interface~mallocate proc~i_free2->proc~chk_free proc~chk_free->proc~fatal_error proc~mat_mult->proc~reallocate_mat proc~chk_alloc->proc~fatal_error proc~mat_andnot->proc~reallocate_mat proc~compress_data->interface~mallocate proc~screening_rules->proc~fatal_error proc~sort_ivec->interface~mallocate proc~sort_ivec->interface~mfree proc~rotation_matrix->proc~fatal_error

Called by

proc~~mmpol_init_from_mmp~~CalledByGraph proc~mmpol_init_from_mmp mmpol_init_from_mmp proc~ommp_init_mmp ommp_init_mmp proc~ommp_init_mmp->proc~mmpol_init_from_mmp proc~c_ommp_init_mmp C_ommp_init_mmp proc~c_ommp_init_mmp->proc~ommp_init_mmp

Contents

Source Code

mmpol_init_from_mmp

Source Code

    subroutine mmpol_init_from_mmp(input_file, sys_obj)
        !! This function read a .mmp file (revision 2 and 3) are supported
        !! and initialize all the quantities need to describe the environment
        !! within this library.

        use mod_mmpol, only: mmpol_prepare, mmpol_init
        use mod_electrostatics, only: set_screening_parameters
        
        use mod_memory, only: ip, rp, mfree, mallocate, memory_init
        use mod_io, only: iof_mmpinp
        use mod_adjacency_mat, only: adj_mat_from_conn
        use mod_utils, only: skip_lines
        use mod_constants, only: angstrom2au, OMMP_VERBOSE_DEBUG
        use mod_constants, only: mscale_wang_al, &
                                 pscale_wang_al, &
                                 dscale_wang_al, &
                                 uscale_wang_al, & 
                                 mscale_wang_dl, &
                                 pscale_wang_dl, &
                                 dscale_wang_dl, &
                                 uscale_wang_dl, & 
                                 mscale_amoeba, &
                                 pscale_amoeba, &
                                 dscale_amoeba, &
                                 uscale_amoeba, &
                                 pscale_intra_amoeba, &
                                 thole_scale_wang_dl, &
                                 thole_scale_wang_al

        implicit none
    
        type(ommp_system), intent(inout) :: sys_obj
        !! System data structure to be initialized

        character(len=*), intent(in) :: input_file
        !! name of the input file

        ! read the input for the mmpol calculation and process it.
        integer(ip) :: input_revision
        integer(ip) :: my_mm_atoms, my_pol_atoms, my_ff_type, my_ff_rules
        integer(ip) :: my_ld_cart
        
        integer(ip) :: i, ist
        
        real(rp), allocatable :: my_pol(:), my_cmm(:,:), my_q(:,:)
        integer(ip), allocatable :: pol_atoms_list(:), my_ip11(:,:)
        
        integer(ip), parameter :: maxn12 = 8
        !! maximum number of neighbour atoms
        integer(ip), parameter :: maxpgp = 120
        !! maximum number of members for the same polarization group
        real(rp), parameter :: thres = 1e-8
        !! Threshold used to decide if a polarizability is zero
        ! TODO why we do not use eps_rp directly?

        integer(ip), allocatable :: i12(:,:)
        character(len=OMMP_STR_CHAR_MAX) :: msg
        logical :: fex
       
        call time_push()
        write(msg, "(A)") "Reading MMP file: "//input_file(1:len(trim(input_file)))
        call ommp_message(msg, OMMP_VERBOSE_DEBUG)

        inquire(file=input_file(1:len(trim(input_file))), exist=fex)
        if(.not. fex) then
            call fatal_error("Input file '"//&
                &input_file(1:len(trim(input_file)))//"' does not exist.")
        end if

        open(unit=iof_mmpinp, &
             file=input_file(1:len(trim(input_file))), &
             form='formatted', &
             access='sequential', &
             iostat=ist, &
             action='read')

        if(ist /= 0) then
            call fatal_error('Error while opening MMP input file')
        end if

        call ommp_message("Reading input parameters", OMMP_VERBOSE_DEBUG)

        ! Read input revision, supported revisions are 2 and 3.
        read(iof_mmpinp,*) input_revision
        
        if (input_revision /= 3 .and. input_revision /= 2) &
            call fatal_error('MMP file is only supported with version 2 or 3.')

        !TODO
        call memory_init(.false., 0.0_rp)
        
        call skip_lines(iof_mmpinp, 2)
        read(iof_mmpinp,*) my_ff_type
        read(iof_mmpinp,*) my_ff_rules
        
        if(input_revision == 3) then
            call skip_lines(iof_mmpinp, 17)
        else if(input_revision == 2) then
            call skip_lines(iof_mmpinp, 15)
        end if
        
        read(iof_mmpinp,*) my_mm_atoms
        
        if(my_ff_type == 1) then
            my_ld_cart = 10
        else
            my_ld_cart = 1
        end if
        
        call ommp_message("Allocating memory", OMMP_VERBOSE_DEBUG)

        call mallocate('mmpol_init_from_mmp [my_cmm]', 3_ip, my_mm_atoms, my_cmm)
        call mallocate('mmpol_init_from_mmp [my_q]', my_ld_cart, my_mm_atoms, my_q)
        call mallocate('mmpol_init_from_mmp [my_pol]', my_mm_atoms, my_pol)
        call mallocate('mmpol_init_from_mmp [my_pol_atoms]', my_mm_atoms, pol_atoms_list)
        call mallocate('mmpol_init_from_mmp [my_ip11]', maxpgp, my_mm_atoms, my_ip11)
       
        call skip_lines(iof_mmpinp, my_mm_atoms+1) ! Skip a zero and the section of atomic numbers
        
        call ommp_message("Reading coordinates", OMMP_VERBOSE_DEBUG)
        ! coordinates:
        do i = 1, my_mm_atoms
            read(iof_mmpinp,*) my_cmm(1:3,i)
        end do

        call skip_lines(iof_mmpinp, my_mm_atoms) ! Skip section of residues number

        call ommp_message("Reading fixed multipoles", OMMP_VERBOSE_DEBUG)
        ! charges/multipoles:
        do i = 1, my_mm_atoms
            read(iof_mmpinp,*) my_q(1:my_ld_cart,i)
        end do

        call ommp_message("Reading polarizabilities", OMMP_VERBOSE_DEBUG)
        ! polarizabilities:
        my_pol = 0.0_rp
        do i = 1, my_mm_atoms
            read(iof_mmpinp,*) my_pol(i)
        end do

        call ommp_message("Processing polarizabilities", OMMP_VERBOSE_DEBUG)
        ! count how many atoms are polarizable:
        ! TODO this is more efficiently done with pack and count
        my_pol_atoms = 0
        pol_atoms_list(:) = 0
        do i = 1, my_mm_atoms
            if (my_pol(i) > thres) then
                my_pol_atoms = my_pol_atoms + 1
                pol_atoms_list(my_pol_atoms) = i
            end if
        end do
        
        ! remove null polarizabilities from the list
        do i = 1, my_pol_atoms
            my_pol(i) = my_pol(pol_atoms_list(i))
        end do
        my_pol(my_pol_atoms+1:my_mm_atoms) = 0.0_rp
        
        call ommp_message("Initializing open mmpol module", OMMP_VERBOSE_DEBUG)
        ! mmpol module initialization
        call mmpol_init(sys_obj, my_ff_type, my_mm_atoms, my_pol_atoms)
        if(sys_obj%amoeba) then
                call set_screening_parameters(sys_obj%eel, &
                                          mscale_amoeba, &
                                          pscale_amoeba, &
                                          dscale_amoeba, &
                                          uscale_amoeba, &
                                          pscale_intra_amoeba)
            sys_obj%eel%thole_scale = 0.0_rp
        else
            if(my_ff_rules == 1) then
                call set_screening_parameters(sys_obj%eel, &
                                              mscale_wang_dl, &
                                              pscale_wang_dl, &
                                              dscale_wang_dl, &
                                              uscale_wang_dl)
                sys_obj%eel%thole_scale = thole_scale_wang_dl
            else if(my_ff_rules == 0) then
                call set_screening_parameters(sys_obj%eel, &
                                              mscale_wang_al, &
                                              pscale_wang_al, &
                                              dscale_wang_al, &
                                              uscale_wang_al)
                sys_obj%eel%thole_scale = thole_scale_wang_al
            end if
        end if 
        
        call ommp_message("Converting input units to A.U.", OMMP_VERBOSE_DEBUG)
        ! Copy data in the correct units (this means AU)
        sys_obj%top%cmm = my_cmm * angstrom2au
        sys_obj%eel%q = my_q
        if(sys_obj%amoeba) then
            my_pol = my_pol*angstrom2au**3
        end if
        sys_obj%eel%pol = my_pol(1:my_pol_atoms)
        sys_obj%eel%polar_mm = pol_atoms_list(1:my_pol_atoms)
        
        call mfree('mmpol_init_from_mmp [my_cmm]', my_cmm)
        call mfree('mmpol_init_from_mmp [my_q]', my_q)
        call mfree('mmpol_init_from_mmp [pol]', my_pol)

        call ommp_message("Processing connectivity informations", OMMP_VERBOSE_DEBUG)

        ! 1-2 connectivity:
        call mallocate('mmpol_init_from_mmp [i12]', maxn12, my_mm_atoms, i12)
        do i = 1, my_mm_atoms
            read(iof_mmpinp,*) i12(1:maxn12,i)
        end do
        
        ! Writes the adjacency matrix in Yale sparse format in conn(1)
        call adj_mat_from_conn(i12, sys_obj%top%conn(1))
        call mfree('mmpol_init_from_mmp [i12]', i12)
         
        if(sys_obj%amoeba) then
            ! group 11 connectivity:
            ! (to be replaced with polarization group)
            do i = 1, my_mm_atoms
                read(iof_mmpinp,*) my_ip11(1:maxpgp,i)
            end do

            call polgroup11_to_mm2pg(my_ip11, sys_obj%eel%mmat_polgrp)
            call mfree('mmpol_init_from_mmp [my_ip11]', my_ip11)

            ! information to rotate the multipoles to the lab frame.
            ! mol_frame, iz, ix, iy:
            do i = 1, my_mm_atoms
                read(iof_mmpinp,*) sys_obj%eel%mol_frame(i), &
                                   sys_obj%eel%iz(i), sys_obj%eel%ix(i), &
                                   sys_obj%eel%iy(i)
            end do
        end if
        ! now, process the input, create all the required arrays 
        ! and the correspondence lists:
     
        call ommp_message("Populating utility arrays", OMMP_VERBOSE_DEBUG)
        call mmpol_prepare(sys_obj)
        close(iof_mmpinp) 
        call ommp_message("Initialization from MMP file done.", OMMP_VERBOSE_DEBUG)
        call time_pull('MMPol initialization from .mmp file')
    end subroutine mmpol_init_from_mmp