rotate_multipoles Subroutine

subroutine rotate_multipoles(eel)

Uses

  • proc~~rotate_multipoles~~UsesGraph proc~rotate_multipoles rotate_multipoles module~mod_electrostatics mod_electrostatics proc~rotate_multipoles->module~mod_electrostatics module~mod_memory mod_memory proc~rotate_multipoles->module~mod_memory module~mod_electrostatics->module~mod_memory module~mod_io mod_io module~mod_electrostatics->module~mod_io module~mod_profiling mod_profiling module~mod_electrostatics->module~mod_profiling module~mod_topology mod_topology module~mod_electrostatics->module~mod_topology module~mod_constants mod_constants module~mod_electrostatics->module~mod_constants module~mod_adjacency_mat mod_adjacency_mat module~mod_electrostatics->module~mod_adjacency_mat 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_profiling->module~mod_memory module~mod_profiling->module~mod_io module~mod_profiling->module~mod_constants module~mod_topology->module~mod_memory module~mod_topology->module~mod_adjacency_mat module~mod_constants->iso_c_binding module~mod_adjacency_mat->module~mod_memory

this routine rotates the atomic multipoles from the molecular frame where they are defined as force field parameters to the lab frame. if required, it also computes the contribution to the forces that stems from the derivatives of the rotation matrices, sometimes referred to as "torques". for the latter task, it uses the field and field gradient from at the multipoles, which is passed in def. for consistency reasons, def is dimensioned (ld_cder,mm_atoms).

Arguments

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

Calls

proc~~rotate_multipoles~~CallsGraph proc~rotate_multipoles rotate_multipoles proc~rotation_matrix rotation_matrix proc~rotate_multipoles->proc~rotation_matrix proc~fatal_error fatal_error proc~rotation_matrix->proc~fatal_error proc~ommp_message ommp_message proc~fatal_error->proc~ommp_message proc~close_output close_output proc~fatal_error->proc~close_output proc~close_output->proc~ommp_message

Called by

proc~~rotate_multipoles~~CalledByGraph proc~rotate_multipoles rotate_multipoles proc~mmpol_prepare mmpol_prepare proc~mmpol_prepare->proc~rotate_multipoles proc~init_eel_for_link_atom init_eel_for_link_atom proc~init_eel_for_link_atom->proc~rotate_multipoles proc~update_coordinates update_coordinates proc~update_coordinates->proc~rotate_multipoles proc~mmpol_init_from_mmp mmpol_init_from_mmp proc~mmpol_init_from_mmp->proc~mmpol_prepare proc~ommp_create_link_atom ommp_create_link_atom proc~ommp_create_link_atom->proc~init_eel_for_link_atom proc~mmpol_init_from_xyz mmpol_init_from_xyz proc~mmpol_init_from_xyz->proc~mmpol_prepare proc~c_ommp_update_coordinates C_ommp_update_coordinates proc~c_ommp_update_coordinates->proc~update_coordinates proc~ommp_system_from_qm_helper ommp_system_from_qm_helper proc~ommp_system_from_qm_helper->proc~mmpol_prepare proc~numerical_geomgrad numerical_geomgrad proc~numerical_geomgrad->proc~update_coordinates proc~update_coordinates_qmmm update_coordinates_qmmm proc~update_coordinates_qmmm->proc~update_coordinates proc~ommp_init_mmp ommp_init_mmp proc~ommp_init_mmp->proc~mmpol_init_from_mmp proc~c_ommp_system_from_qm_helper C_ommp_system_from_qm_helper proc~c_ommp_system_from_qm_helper->proc~ommp_system_from_qm_helper program~test_si_geomgrad test_SI_geomgrad program~test_si_geomgrad->proc~ommp_system_from_qm_helper proc~ommp_init_xyz ommp_init_xyz proc~ommp_init_xyz->proc~mmpol_init_from_xyz proc~c_ommp_create_link_atom C_ommp_create_link_atom proc~c_ommp_create_link_atom->proc~ommp_create_link_atom program~test_si_geomgrad_num test_SI_geomgrad_num program~test_si_geomgrad_num->proc~ommp_system_from_qm_helper proc~num_grd_print num_grd_print program~test_si_geomgrad_num->proc~num_grd_print proc~num_grd_print->proc~numerical_geomgrad proc~numerical_geomgrad_qmmm numerical_geomgrad_qmmm proc~numerical_geomgrad_qmmm->proc~update_coordinates_qmmm proc~c_ommp_init_mmp C_ommp_init_mmp proc~c_ommp_init_mmp->proc~ommp_init_mmp proc~c_ommp_init_xyz C_ommp_init_xyz proc~c_ommp_init_xyz->proc~ommp_init_xyz

Contents

Source Code


Source Code

subroutine rotate_multipoles(eel)
    !! this routine rotates the atomic multipoles from the molecular frame
    !! where they are defined as force field parameters to the lab frame.
    !! if required, it also computes the contribution to the forces that
    !! stems from the derivatives of the rotation matrices, sometimes 
    !! referred to as "torques".
    !! for the latter task, it uses the field and field gradient from 
    !! at the multipoles, which is passed in def. 
    !! for consistency reasons, def is dimensioned (ld_cder,mm_atoms). 
    
    use mod_memory, only: ip, rp
    use mod_electrostatics, only: ommp_electrostatics_type
    
    implicit none
  
    type(ommp_electrostatics_type), intent(inout) :: eel

    integer(ip) :: j, jx, jy, jz 
    real(rp), dimension(3,3) :: r, rt, qua, rqua, tmp
    real(rp), dimension(3,3,3) :: dri, driz, drix, driy

    ! loop over the mm sites and build the rotation matrices.
    do j = 1, eel%top%mm_atoms
        jz = eel%iz(j)
        if(jz == 0) jz = j
        jx = eel%ix(j)
        if(jx == 0) jx = j
        jy = eel%iy(j)
        if(jy == 0) jy = j
        
        call rotation_matrix(.false., & 
                             eel%top%cmm(:,j), eel%top%cmm(:,jx), &
                             eel%top%cmm(:,jy), eel%top%cmm(:,jz), &
                             eel%mol_frame(j), &
                             r, dri, driz, drix, driy)
        
        rt = transpose(r)
        ! copy the monopole:
        eel%q(1,j) = eel%q0(1,j)

        ! rotate the dipole 
        eel%q(2:4,j) = matmul(r,eel%q0(2:4,j))

        ! exctract, rotate and put back the quadrupole:
        qua(_x_,_x_)  = eel%q0(4+_xx_,j)
        qua(_x_,_y_)  = eel%q0(4+_xy_,j)
        qua(_x_,_z_)  = eel%q0(4+_xz_,j)
        qua(_y_,_x_)  = eel%q0(4+_yx_,j)
        qua(_y_,_y_)  = eel%q0(4+_yy_,j)
        qua(_y_,_z_)  = eel%q0(4+_yz_,j)
        qua(_z_,_x_)  = eel%q0(4+_zx_,j)
        qua(_z_,_y_)  = eel%q0(4+_zy_,j)
        qua(_z_,_z_)  = eel%q0(4+_zz_,j)
        
        tmp  = matmul(r,qua)
        rqua = matmul(tmp,rt)
        eel%q(4+_xx_,j)  = rqua(_x_,_x_)
        eel%q(4+_yy_,j)  = rqua(_y_,_y_)
        eel%q(4+_zz_,j)  = rqua(_z_,_z_)
        eel%q(4+_xy_,j)  = rqua(_x_,_y_)
        eel%q(4+_xz_,j)  = rqua(_x_,_z_)
        eel%q(4+_yz_,j)  = rqua(_y_,_z_)
    end do
end subroutine rotate_multipoles