Extract_dep_list new symm

New dependency list symmetrizatoin, two alternatives: distributed, or
replicated after allgather.
fnd_owner
Salvatore Filippone 5 years ago
parent c845a7881e
commit 43c1ccfc87

@ -79,7 +79,7 @@ subroutine psi_i_dl_check(dep_list,dl_lda,np,length_dl)
! ...add proc to proc2 s dep_list.....',proc,proc2 ! ...add proc to proc2 s dep_list.....',proc,proc2
length_dl(proc2) = length_dl(proc2)+1 length_dl(proc2) = length_dl(proc2)+1
if (length_dl(proc2) > size(dep_list,1)) then if (length_dl(proc2) > size(dep_list,1)) then
write(psb_err_unit,*)'error in crea_halo', proc2,proc,& write(psb_err_unit,*)'error in dl_check', proc2,proc,&
& length_dl(proc2),'>',size(dep_list,1) & length_dl(proc2),'>',size(dep_list,1)
endif endif
dep_list(length_dl(proc2),proc2) = proc dep_list(length_dl(proc2),proc2) = proc

@ -126,6 +126,7 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
use psb_const_mod use psb_const_mod
use psb_error_mod use psb_error_mod
use psb_desc_mod use psb_desc_mod
use psb_sort_mod
implicit none implicit none
#ifdef MPI_H #ifdef MPI_H
include 'mpif.h' include 'mpif.h'
@ -145,7 +146,8 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
integer(psb_ipk_) :: i,pointer_dep_list,proc,j,err_act integer(psb_ipk_) :: i,pointer_dep_list,proc,j,err_act
integer(psb_ipk_) :: err integer(psb_ipk_) :: err
integer(psb_ipk_) :: debug_level, debug_unit integer(psb_ipk_) :: debug_level, debug_unit
integer(psb_mpk_) :: iictxt, icomm, me, np, dl_mpi, minfo integer(psb_mpk_) :: iictxt, icomm, me, np, minfo
logical, parameter :: dist_symm_list=.false., print_dl=.false.
character name*20 character name*20
name='psi_extrct_dl' name='psi_extrct_dl'
@ -156,13 +158,14 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
info = psb_success_ info = psb_success_
call psb_info(iictxt,me,np) call psb_info(iictxt,me,np)
allocate(itmp(np+1),length_dl(0:np),stat=info) allocate(itmp(2*np+1),length_dl(0:np),stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
info=psb_err_alloc_dealloc_ info=psb_err_alloc_dealloc_
goto 9999 goto 9999
end if end if
do i=0,np do i=0,np
length_dl(i) = 0 length_dl(i) = 0
itmp(i+1) = -1
enddo enddo
i=1 i=1
if (debug_level >= psb_debug_inner_)& if (debug_level >= psb_debug_inner_)&
@ -249,32 +252,116 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
endif endif
length_dl(me)=pointer_dep_list-1 length_dl(me)=pointer_dep_list-1
dl_lda = max(length_dl(me),1) if (dist_symm_list) then
call psb_max(iictxt, dl_lda) call psb_realloc(length_dl(me),itmp,info)
! call psi_symm_dep_list(itmp,ictxt,info)
! This doubling of DL_LDA is not 100% safe, dl_lda = max(size(itmp),1)
! but should work most of the time. call psb_max(iictxt, dl_lda)
! Will need to be improved later, perhaps move call psb_realloc(dl_lda,itmp,info)
! from a 2D allocation (ellpack style) to ! dl_lda = min(np,2*dl_lda)
! a 1D allocation (csr like). allocate(dep_list(dl_lda,0:np),stat=info)
! if (info /= psb_success_) then
dl_lda = min(2*dl_lda,np+1) call psb_errpush(psb_err_from_subroutine_,name,a_err='Allocate')
allocate(dep_list(dl_lda,0:np),stat=info) goto 9999
if (info /= psb_success_) then end if
call psb_errpush(psb_err_from_subroutine_,name,a_err='Allocate')
goto 9999
end if
call psb_sum(iictxt,length_dl(0:np)) call psb_sum(iictxt,length_dl(0:np))
icomm = psb_get_mpicomm(iictxt) icomm = psb_get_mpicomm(iictxt)
call mpi_allgather(itmp,dl_lda,psb_mpi_ipk_,& call mpi_allgather(itmp,dl_lda,psb_mpi_ipk_,&
& dep_list,dl_lda,psb_mpi_ipk_,icomm,minfo) & dep_list,dl_lda,psb_mpi_ipk_,icomm,minfo)
info = minfo info = minfo
if (info == 0) deallocate(itmp,stat=info) if (info == 0) deallocate(itmp,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
info=psb_err_alloc_dealloc_ info=psb_err_alloc_dealloc_
goto 9999 goto 9999
endif
else
block
integer(psb_ipk_), allocatable :: list1(:,:), ldl2(:), list2(:,:)
integer(psb_ipk_) :: i,j,ip,dlsym, ldu, mdl, l1, l2
dl_lda = max(length_dl(me),1)
call psb_max(iictxt, dl_lda)
allocate(dep_list(dl_lda,0:np),stat=info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='Allocate')
goto 9999
end if
call psb_sum(iictxt,length_dl(0:np))
icomm = psb_get_mpicomm(iictxt)
call mpi_allgather(itmp,dl_lda,psb_mpi_ipk_,&
& dep_list,dl_lda,psb_mpi_ipk_,icomm,minfo)
info = minfo
if (info /= psb_success_) then
info=psb_err_alloc_dealloc_
goto 9999
endif
allocate(ldl2(0:np),stat=info)
ldl2 = 0
do j=0, np-1
do i=1,length_dl(j)
ip = dep_list(i,j)
ldl2(ip) = ldl2(ip) + 1
end do
end do
dlsym = maxval(ldl2)
allocate(list2(dlsym,0:np),stat=info)
call move_alloc(dep_list,list1)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='Allocate')
goto 9999
end if
ldl2 = 0
do j=0, np-1
do i=1,length_dl(j)
ip = list1(i,j)
ldl2(ip) = ldl2(ip) + 1
list2(ldl2(ip),ip) = j
end do
end do
mdl = 0
do j = 0, np-1
l1 = length_dl(j)
l2 = ldl2(j)
itmp(1:l1) = list1(1:l1,j)
itmp(l1+1:l1+l2) = list2(1:l2,j)
ldu = l1 + l2
!if (me == 0) write(0,*) 'Iter ',j,':',l1,l2,':',itmp(1:l1),':',itmp(l1+1:l1+l2)
call psb_msort_unique(itmp(1:l1+l2),ldu)
mdl = max(mdl, ldu)
!if (me == 0) write(0,*) 'Iter ',j,':',ldu,':',itmp(1:ldu)
end do
dl_lda = mdl
allocate(dep_list(dl_lda,0:np),stat=info)
dep_list = -1
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='Allocate')
goto 9999
end if
do j = 0, np-1
l1 = length_dl(j)
l2 = ldl2(j)
itmp(1:l1) = list1(1:l1,j)
itmp(l1+1:l1+l2) = list2(1:l2,j)
ldu = l1 + l2
call psb_msort_unique(itmp(1:l1+l2),ldu)
length_dl(j) = ldu
dep_list(1:ldu,j) = itmp(1:ldu)
end do
end block
endif endif
if (print_dl) then
if (me == 0) then
write(0,*) ' Dep_list '
do i=0,np-1
j = length_dl(i)
write(0,*) 'Proc ',i,':',dep_list(1:j,i)
end do
flush(0)
end if
call psb_barrier(ictxt)
end if
call psb_erractionrestore(err_act) call psb_erractionrestore(err_act)
return return

@ -37,13 +37,13 @@
! !
! Arguments: ! Arguments:
! !
subroutine psi_symm_dep_list(rvsz,adj,idxmap,info,flag) subroutine psi_symm_dep_list_inrv(rvsz,adj,ictxt,info)
use psb_serial_mod use psb_serial_mod
use psb_const_mod use psb_const_mod
use psb_error_mod use psb_error_mod
use psb_penv_mod use psb_penv_mod
use psb_realloc_mod use psb_realloc_mod
use psb_indx_map_mod, psb_protect_name => psi_symm_dep_list use psb_indx_map_mod, psb_protect_name => psi_symm_dep_list_inrv
#ifdef MPI_MOD #ifdef MPI_MOD
use mpi use mpi
#endif #endif
@ -54,30 +54,23 @@ subroutine psi_symm_dep_list(rvsz,adj,idxmap,info,flag)
#endif #endif
integer(psb_mpk_), intent(inout) :: rvsz(0:) integer(psb_mpk_), intent(inout) :: rvsz(0:)
integer(psb_ipk_), allocatable, intent(inout) :: adj(:) integer(psb_ipk_), allocatable, intent(inout) :: adj(:)
class(psb_indx_map), intent(in) :: idxmap integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: flag
! !
integer(psb_ipk_), allocatable :: ladj(:) integer(psb_ipk_), allocatable :: ladj(:)
integer(psb_mpk_), allocatable :: sdsz(:) integer(psb_mpk_) :: icomm, minfo
integer(psb_mpk_) :: icomm, minfo, iictxt
integer(psb_ipk_) :: i,n_row,n_col,err_act,hsize,ip,& integer(psb_ipk_) :: i,n_row,n_col,err_act,hsize,ip,&
& last_ih, last_j, nidx, nrecv, nadj, flag_ & last_ih, last_j, nidx, nrecv, nadj, flag_
integer(psb_lpk_) :: mglob, ih integer(psb_lpk_) :: mglob, ih
integer(psb_ipk_) :: ictxt,np,me integer(psb_ipk_) :: np,me
character(len=20) :: name character(len=20) :: name
info = psb_success_ info = psb_success_
name = 'psi_symm_dep_list' name = 'psi_symm_dep_list'
call psb_erractionsave(err_act) call psb_erractionsave(err_act)
ictxt = idxmap%get_ctxt() icomm = psb_get_mpicomm(ictxt)
icomm = idxmap%get_mpic()
mglob = idxmap%get_gr()
n_row = idxmap%get_lr()
n_col = idxmap%get_lc()
iictxt = ictxt
call psb_info(ictxt, me, np) call psb_info(ictxt, me, np)
@ -87,51 +80,15 @@ subroutine psi_symm_dep_list(rvsz,adj,idxmap,info,flag)
goto 9999 goto 9999
endif endif
if (.not.(idxmap%is_valid())) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='invalid idxmap')
goto 9999
end if
nadj = size(adj) nadj = size(adj)
!
! Am getting this from the outside, possibly from the wrapper call
!
if (size(rvsz)<np) then if (size(rvsz)<np) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='invalid rvsz') call psb_errpush(psb_err_from_subroutine_,name,a_err='invalid rvsz')
goto 9999 goto 9999
end if end if
if (present(flag)) then
flag_ = flag
else
flag_ = psi_symm_flag_norv_
end if
select case(flag_)
case(psi_symm_flag_norv_, psi_symm_flag_inrv_)
! Ok, keep going
case default
call psb_errpush(psb_err_from_subroutine_,name,a_err='invalid flag')
goto 9999
end select
if (flag_ == psi_symm_flag_norv_) then
! write(0,*) me,name,' Going through ',nidx,nadj
Allocate(sdsz(0:np-1), stat=info)
!
! First, send sizes according to adjcncy list
!
sdsz = 0
do i=1, nadj
sdsz(adj(i)) = 1
end do
!write(0,*)me,' Check on sizes into a2a:',adj(:),nadj,':',sdsz(:)
call mpi_alltoall(sdsz,1,psb_mpi_mpk_,&
& rvsz,1,psb_mpi_mpk_,icomm,minfo)
end if
nrecv = 0 nrecv = 0
do ip=0, np-1 do ip=0, np-1
if (rvsz(ip) > 0) nrecv = nrecv + 1 if (rvsz(ip) > 0) nrecv = nrecv + 1
@ -145,12 +102,86 @@ subroutine psi_symm_dep_list(rvsz,adj,idxmap,info,flag)
do ip=0, np-1 do ip=0, np-1
if (rvsz(ip)>0) then if (rvsz(ip)>0) then
nadj = nadj + 1 nadj = nadj + 1
ladj(nadj+1:nadj+nrecv) = ip ladj(nadj) = ip
end if end if
end do end do
call psb_msort_unique(ladj,nadj) call psb_msort_unique(ladj,nadj)
call psb_realloc(nadj,adj,info) call psb_realloc(nadj,adj,info)
adj(1:nadj) = ladj(1:nadj) if (info == 0) adj(1:nadj) = ladj(1:nadj)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(ictxt,err_act)
return
end subroutine psi_symm_dep_list_inrv
subroutine psi_symm_dep_list_norv(adj,ictxt,info)
use psb_serial_mod
use psb_const_mod
use psb_error_mod
use psb_penv_mod
use psb_realloc_mod
use psb_indx_map_mod, psb_protect_name => psi_symm_dep_list_norv
#ifdef MPI_MOD
use mpi
#endif
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
integer(psb_ipk_), allocatable, intent(inout) :: adj(:)
integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(out) :: info
!
integer(psb_mpk_), allocatable :: rvsz(:), sdsz(:)
integer(psb_mpk_) :: icomm, minfo
integer(psb_ipk_) :: i,n_row,n_col,err_act,hsize,ip,&
& last_ih, last_j, nidx, nrecv, nadj
integer(psb_ipk_) :: mglob, ih
integer(psb_ipk_) :: np,me
character(len=20) :: name
info = psb_success_
name = 'psi_symm_dep_list'
call psb_erractionsave(err_act)
icomm = psb_get_mpicomm(ictxt)
call psb_info(ictxt, me, np)
if (np == -1) then
info = psb_err_context_error_
call psb_errpush(info,name)
goto 9999
endif
nadj = size(adj)
! write(0,*) me,name,' Going through ',nidx,nadj
Allocate(sdsz(0:np-1), rvsz(0:np-1), stat=info)
!
! First, send sizes according to adjcncy list
!
sdsz = 0
do i=1, nadj
sdsz(adj(i)) = 1
end do
!write(0,*)me,' Check on sizes into a2a:',adj(:),nadj,':',sdsz(:)
call mpi_alltoall(sdsz,1,psb_mpi_mpk_,&
& rvsz,1,psb_mpi_mpk_,icomm,minfo)
if (minfo == 0) call psi_symm_dep_list(rvsz,adj,ictxt,info)
if ((minfo /=0).or.(info /= 0)) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='inner call symm_dep')
goto 9999
end if
call psb_erractionrestore(err_act) call psb_erractionrestore(err_act)
return return
@ -159,4 +190,4 @@ subroutine psi_symm_dep_list(rvsz,adj,idxmap,info,flag)
return return
end subroutine psi_symm_dep_list end subroutine psi_symm_dep_list_norv

@ -313,17 +313,23 @@ module psb_indx_map_mod
integer, parameter :: psi_symm_flag_norv_ = 0 integer, parameter :: psi_symm_flag_norv_ = 0
integer, parameter :: psi_symm_flag_inrv_ = 1 integer, parameter :: psi_symm_flag_inrv_ = 1
interface interface psi_symm_dep_list
subroutine psi_symm_dep_list(rvsz,adj,idxmap,info,flag) subroutine psi_symm_dep_list_inrv(rvsz,adj,ictxt,info)
import :: psb_indx_map, psb_ipk_, psb_lpk_, psb_mpk_ import :: psb_indx_map, psb_ipk_, psb_lpk_, psb_mpk_
implicit none implicit none
integer(psb_mpk_), intent(inout) :: rvsz(:) integer(psb_mpk_), intent(inout) :: rvsz(0:)
integer(psb_ipk_), intent(in) :: adj(:) integer(psb_ipk_), allocatable, intent(inout) :: adj(:)
class(psb_indx_map), intent(in) :: idxmap integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: flag end subroutine psi_symm_dep_list_inrv
end subroutine psi_symm_dep_list subroutine psi_symm_dep_list_norv(adj,ictxt,info)
end interface import :: psb_indx_map, psb_ipk_, psb_lpk_, psb_mpk_
implicit none
integer(psb_ipk_), allocatable, intent(inout) :: adj(:)
integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(out) :: info
end subroutine psi_symm_dep_list_norv
end interface psi_symm_dep_list
contains contains

Loading…
Cancel
Save