Builds the adjacency matrix of polarization groups starting from atomic adjacency matrix and list of polarization groups indices.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(ommp_electrostatics_type), | intent(in) | :: | eel | |||
type(yale_sparse), | intent(out) | :: | adj |
The group adjacency matrix to be saved. |
subroutine build_pg_adjacency_matrix(eel, adj)
!! Builds the adjacency matrix of polarization groups starting from
!! atomic adjacency matrix and list of polarization groups indices.
use mod_adjacency_mat, only: reallocate_mat
implicit none
type(ommp_electrostatics_type), intent(in) :: eel
type(yale_sparse), intent(out) :: adj
!! The group adjacency matrix to be saved.
integer(ip) :: npg, pg1, atm1, atm2, i, j
npg = eel%polgrp_mmat%n
adj%n = npg
allocate(adj%ri(adj%n+1))
allocate(adj%ci(adj%n*2))
adj%ri(1) = 1
do pg1=1, npg
! For each polarization group
adj%ri(pg1+1) = adj%ri(pg1)
do i=eel%polgrp_mmat%ri(pg1), eel%polgrp_mmat%ri(pg1+1)-1
! Loop on every atom of the group
atm1 = eel%polgrp_mmat%ci(i)
do j=eel%top%conn(1)%ri(atm1), eel%top%conn(1)%ri(atm1+1)-1
! Loop on each connected atom...
atm2 = eel%top%conn(1)%ci(j)
! If the two atoms are in different PG, then the two
! polarization groups are connected.
if(eel%mmat_polgrp(atm1) /= eel%mmat_polgrp(atm2) .and. &
! if the group is not already present in the matrix
all(adj%ci(adj%ri(pg1):adj%ri(pg1+1)-1) /= eel%mmat_polgrp(atm2))) then
adj%ci(adj%ri(pg1+1)) = eel%mmat_polgrp(atm2)
adj%ri(pg1+1) = adj%ri(pg1+1) + 1
if(adj%ri(pg1+1) > size(adj%ci)) then
! If matrix is too small, it could be enlarged...
call reallocate_mat(adj, size(adj%ci)+adj%n)
end if
end if
end do
end do
end do
! Finally trim the output matrix
call reallocate_mat(adj, adj%ri(adj%n+1)-1)
end subroutine build_pg_adjacency_matrix