Optimized version of X_glob_transpose: do not go through LX_glob_transpose.

merge-paraggr-newops^2
Salvatore Filippone 5 years ago
parent dbc20d482e
commit 91f737475e

@ -113,7 +113,6 @@ subroutine psb_lc_coo_glob_transpose(ain,desc_r,info,atrans,desc_c,desc_rx)
integer(psb_ipk_) :: ictxt, np,me integer(psb_ipk_) :: ictxt, np,me
integer(psb_ipk_) :: counter,proc, err_act, j integer(psb_ipk_) :: counter,proc, err_act, j
integer(psb_lpk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,& integer(psb_lpk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,&
& irmin,icmin,irmax,icmax,&
& l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzt, nzd & l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzt, nzd
integer(psb_mpk_) :: icomm, minfo integer(psb_mpk_) :: icomm, minfo
integer(psb_mpk_), allocatable :: brvindx(:), & integer(psb_mpk_), allocatable :: brvindx(:), &
@ -407,19 +406,293 @@ subroutine psb_c_coo_glob_transpose(ain,desc_r,info,atrans,desc_c,desc_rx)
type(psb_desc_type), intent(out), optional :: desc_rx type(psb_desc_type), intent(out), optional :: desc_rx
integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), intent(out) :: info
type(psb_lc_coo_sparse_mat) :: atcoo integer(psb_ipk_) :: ictxt, np,me
integer(psb_ipk_) :: counter,proc, err_act, j
integer(psb_ipk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,&
& l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzd
integer(psb_lpk_) :: nzt, lszr
integer(psb_mpk_) :: icomm, minfo
integer(psb_mpk_), allocatable :: brvindx(:), &
& rvsz(:), bsdindx(:), sdsz(:), tsdx(:), trvx(:)
integer(psb_ipk_), allocatable :: halo_owner(:)
integer(psb_lpk_), allocatable :: iasnd(:), jasnd(:),iarcv(:),jarcv(:)
complex(psb_spk_), allocatable :: valsnd(:)
type(psb_c_coo_sparse_mat), allocatable :: acoo
logical :: rowcnv_,colcnv_,rowscale_,colscale_,outcol_glob_
type(psb_desc_type), pointer :: p_desc_c
character(len=5) :: outfmt_
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
if(psb_get_errstatus() /= 0) return
info=psb_success_
name='mld_glob_transpose'
call psb_erractionsave(err_act)
if (psb_errstatus_fatal()) then
info = psb_err_internal_error_ ; goto 9999
end if
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
ictxt = desc_r%get_context()
icomm = desc_r%get_mpic()
Call psb_info(ictxt, me, np)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),': Start'
if (present(desc_c)) then
p_desc_c => desc_c
else
p_desc_c => desc_r
end if
Allocate(brvindx(np+1),&
& rvsz(np),sdsz(np),bsdindx(np+1), acoo,stat=info)
if (info /= psb_success_) then
info=psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
sdsz(:)=0
rvsz(:)=0
l1 = 0
brvindx(1) = 0
bsdindx(1) = 0
counter=1
idx = 0
idxs = 0
idxr = 0
if (present(atrans)) then if (present(atrans)) then
call ain%cp_to_lcoo(atcoo,info) call ain%cp_to_coo(acoo,info)
else
call ain%mv_to_coo(acoo,info)
end if
!
! Compute number of entries in the
! halo part, sorted by destination process
!
nr = desc_r%get_local_rows()
nc = p_desc_c%get_local_cols()
nzl = acoo%get_nzeros()
hlstart = p_desc_c%get_local_rows()
do k=1, nzl
j = acoo%ja(k)
if (j > hlstart) then
call p_desc_c%indxmap%fnd_halo_owner(j,proc,info)
sdsz(proc+1) = sdsz(proc+1) +1
end if
end do
!
! Exchange sizes, so as to know sends/receives.
! This is different from the halo exchange because the
! number of entries was not precomputed in the descriptor,
! which was vector-oriented and not matrix-entry-oriented
!
call mpi_alltoall(sdsz,1,psb_mpi_mpk_,&
& rvsz,1,psb_mpi_mpk_,icomm,minfo)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mpi_alltoall')
goto 9999
end if
nsnds = count(sdsz /= 0)
nrcvs = count(rvsz /= 0)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done initial alltoall',nsnds,nrcvs
idxs = 0
idxr = 0
counter = 1
Do proc = 0, np-1
bsdindx(proc+1) = idxs
idxs = idxs + sdsz(proc+1)
brvindx(proc+1) = idxr
idxr = idxr + rvsz(proc+1)
Enddo
tsdx = bsdindx
trvx = brvindx
iszr = sum(rvsz)
iszs = sum(sdsz)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Sizes:',&
& ' Send:',sdsz(:),' Receive:',rvsz(:)
call psb_ensure_size(max(iszs,1),iasnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszs,1),jasnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszs,1),valsnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszr,1),iarcv,info)
if (info == psb_success_) call psb_ensure_size(max(iszr,1),jarcv,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='ensure_size')
goto 9999
end if
!
! Now, transpose the matrix, then split between itself
! and the send buffers
!
call acoo%transp()
if (acoo%get_nzeros()/= nzl) then
write(0,*) me,'Something strange upon transpose: ',nzl,acoo%get_nzeros()
end if
nzl = acoo%get_nzeros()
hlstart = p_desc_c%get_local_rows()
nzd = 0
do k = 1, nzl
j = acoo%ia(k)
if (j<=hlstart) then
nzd = nzd + 1
acoo%ia(nzd) = acoo%ia(k)
acoo%ja(nzd) = acoo%ja(k)
acoo%val(nzd) = acoo%val(k)
else
call p_desc_c%indxmap%fnd_halo_owner(j,proc,info)
tsdx(proc+1) = tsdx(proc+1) +1
iasnd(tsdx(proc+1)) = acoo%ia(k)
jasnd(tsdx(proc+1)) = acoo%ja(k)
valsnd(tsdx(proc+1)) = acoo%val(k)
end if
end do
!
! Put halo entries in global numbering
!
call desc_r%indxmap%l2gip(jasnd(1:iszs),info)
call p_desc_c%indxmap%l2gip(iasnd(1:iszs),info)
call acoo%set_nzeros(nzd)
! And exchange data.
! Normally we'll use our SIMPLE A2AV and not MPI, because
! the communication pattern is sparse, so ours is more
! efficient. Using ACOO for the receive buffers.
nzl = acoo%get_nzeros()
call acoo%ensure_size(nzl+iszr)
select case(psb_get_sp_a2av_alg())
case(psb_sp_a2av_smpl_triad_)
call psb_simple_triad_a2av(valsnd,iasnd,jasnd,sdsz,bsdindx,&
& acoo%val(nzl+1:nzl+iszr),iarcv(1:iszr),&
& jarcv(1:iszr),rvsz,brvindx,ictxt,info)
case(psb_sp_a2av_smpl_v_)
call psb_simple_a2av(valsnd,sdsz,bsdindx,&
& acoo%val(nzl+1:nzl+iszr),rvsz,brvindx,ictxt,info)
if (info == psb_success_) call psb_simple_a2av(iasnd,sdsz,bsdindx,&
& iarcv(1:iszr),rvsz,brvindx,ictxt,info)
if (info == psb_success_) call psb_simple_a2av(jasnd,sdsz,bsdindx,&
& jarcv(1:iszr),rvsz,brvindx,ictxt,info)
case(psb_sp_a2av_mpi_)
call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_c_spk_,&
& acoo%val(nzl+1:nzl+iszr),rvsz,brvindx,psb_mpi_c_spk_,icomm,minfo)
if (minfo == mpi_success) &
& call mpi_alltoallv(iasnd,sdsz,bsdindx,psb_mpi_lpk_,&
& iarcv(1:iszr),rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo == mpi_success) &
& call mpi_alltoallv(jasnd,sdsz,bsdindx,psb_mpi_lpk_,&
& jarcv(1:iszr),rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo /= mpi_success) info = minfo
case default
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='wrong A2AV alg selector')
goto 9999
end select
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='alltoallv')
goto 9999
end if
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done alltoallv'
if (present(desc_rx)) then
!
! Extend the appropriate descriptor; started as R but on
! transpose it now describes C
!
call desc_r%clone(desc_rx,info)
call psb_cd_reinit(desc_rx,info)
!
! Take out non-local rows
!
call psb_glob_to_loc(iarcv(1:iszr),p_desc_c,info,iact='I',owned=.true.)
lszr = iszr
call psb_coo_clean_negidx_inner(lszr,iarcv(1:iszr),jarcv(1:iszr),&
& acoo%val(nzl+1:nzl+iszr),nzt,info)
call desc_rx%g2lip_ins(jarcv(1:nzt),info)
call psb_cdasb(desc_rx,info)
acoo%ia(nzl+1:nzl+nzt) = iarcv(1:nzt)
acoo%ja(nzl+1:nzl+nzt) = jarcv(1:nzt)
nzl = nzl + nzt
call acoo%set_nzeros(nzl)
nzl = acoo%get_nzeros()
call acoo%set_sorted(.false.)
!
! Insert to extend descriptor
!
call acoo%set_nrows(p_desc_c%get_local_rows())
call acoo%set_ncols(desc_rx%get_local_cols())
!write(0,*) me,' Trans RX ',acoo%get_nrows(),acoo%get_ncols(),acoo%get_nzeros()
else else
call ain%mv_to_lcoo(atcoo,info) !
!
! Take out non-local rows
!
call psb_glob_to_loc(iarcv(1:iszr),p_desc_c,info,iact='I',owned=.true.)
call psb_glob_to_loc(jarcv(1:iszr),desc_r,info,iact='I')
lszr = iszr
call psb_coo_clean_negidx_inner(lszr,iarcv(1:iszr),jarcv(1:iszr),&
& acoo%val(nzl+1:nzl+iszr),nzt,info)
acoo%ia(nzl+1:nzl+nzt) = iarcv(1:nzt)
acoo%ja(nzl+1:nzl+nzt) = jarcv(1:nzt)
nzl = nzl + nzt
call acoo%set_nzeros(nzl)
nzl = acoo%get_nzeros()
call acoo%set_sorted(.false.)
call acoo%set_nrows(p_desc_c%get_local_rows())
call acoo%set_ncols(desc_r%get_local_cols())
!write(0,*) me,' Trans R- ',acoo%get_nrows(),acoo%get_ncols(),acoo%get_nzeros()
end if end if
if (info == 0) call psb_glob_transpose(atcoo,desc_r,info,desc_c=desc_c,desc_rx=desc_rx)
!!$ write(0,*) me,' Sanity check after rx%g2l :',count(acoo%ja(1:nzl)<0)
call acoo%fix(info)
nzl = acoo%get_nzeros()
if (present(atrans)) then if (present(atrans)) then
call atrans%mv_from_lcoo(atcoo,info) call atrans%mv_from_coo(acoo,info)
else else
call ain%mv_from_lcoo(atcoo,info) call ain%mv_from_coo(acoo,info)
end if end if
Deallocate(brvindx,bsdindx,rvsz,sdsz,&
& iasnd,jasnd,valsnd,&
& iarcv,jarcv,stat=info)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done'
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(ictxt,err_act)
return
end subroutine psb_c_coo_glob_transpose end subroutine psb_c_coo_glob_transpose
subroutine psb_c_simple_glob_transpose_ip(ain,desc_a,info) subroutine psb_c_simple_glob_transpose_ip(ain,desc_a,info)

@ -113,7 +113,6 @@ subroutine psb_ld_coo_glob_transpose(ain,desc_r,info,atrans,desc_c,desc_rx)
integer(psb_ipk_) :: ictxt, np,me integer(psb_ipk_) :: ictxt, np,me
integer(psb_ipk_) :: counter,proc, err_act, j integer(psb_ipk_) :: counter,proc, err_act, j
integer(psb_lpk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,& integer(psb_lpk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,&
& irmin,icmin,irmax,icmax,&
& l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzt, nzd & l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzt, nzd
integer(psb_mpk_) :: icomm, minfo integer(psb_mpk_) :: icomm, minfo
integer(psb_mpk_), allocatable :: brvindx(:), & integer(psb_mpk_), allocatable :: brvindx(:), &
@ -407,19 +406,293 @@ subroutine psb_d_coo_glob_transpose(ain,desc_r,info,atrans,desc_c,desc_rx)
type(psb_desc_type), intent(out), optional :: desc_rx type(psb_desc_type), intent(out), optional :: desc_rx
integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), intent(out) :: info
type(psb_ld_coo_sparse_mat) :: atcoo integer(psb_ipk_) :: ictxt, np,me
integer(psb_ipk_) :: counter,proc, err_act, j
integer(psb_ipk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,&
& l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzd
integer(psb_lpk_) :: nzt, lszr
integer(psb_mpk_) :: icomm, minfo
integer(psb_mpk_), allocatable :: brvindx(:), &
& rvsz(:), bsdindx(:), sdsz(:), tsdx(:), trvx(:)
integer(psb_ipk_), allocatable :: halo_owner(:)
integer(psb_lpk_), allocatable :: iasnd(:), jasnd(:),iarcv(:),jarcv(:)
real(psb_dpk_), allocatable :: valsnd(:)
type(psb_d_coo_sparse_mat), allocatable :: acoo
logical :: rowcnv_,colcnv_,rowscale_,colscale_,outcol_glob_
type(psb_desc_type), pointer :: p_desc_c
character(len=5) :: outfmt_
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
if(psb_get_errstatus() /= 0) return
info=psb_success_
name='mld_glob_transpose'
call psb_erractionsave(err_act)
if (psb_errstatus_fatal()) then
info = psb_err_internal_error_ ; goto 9999
end if
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
ictxt = desc_r%get_context()
icomm = desc_r%get_mpic()
Call psb_info(ictxt, me, np)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),': Start'
if (present(desc_c)) then
p_desc_c => desc_c
else
p_desc_c => desc_r
end if
Allocate(brvindx(np+1),&
& rvsz(np),sdsz(np),bsdindx(np+1), acoo,stat=info)
if (info /= psb_success_) then
info=psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
sdsz(:)=0
rvsz(:)=0
l1 = 0
brvindx(1) = 0
bsdindx(1) = 0
counter=1
idx = 0
idxs = 0
idxr = 0
if (present(atrans)) then if (present(atrans)) then
call ain%cp_to_lcoo(atcoo,info) call ain%cp_to_coo(acoo,info)
else
call ain%mv_to_coo(acoo,info)
end if
!
! Compute number of entries in the
! halo part, sorted by destination process
!
nr = desc_r%get_local_rows()
nc = p_desc_c%get_local_cols()
nzl = acoo%get_nzeros()
hlstart = p_desc_c%get_local_rows()
do k=1, nzl
j = acoo%ja(k)
if (j > hlstart) then
call p_desc_c%indxmap%fnd_halo_owner(j,proc,info)
sdsz(proc+1) = sdsz(proc+1) +1
end if
end do
!
! Exchange sizes, so as to know sends/receives.
! This is different from the halo exchange because the
! number of entries was not precomputed in the descriptor,
! which was vector-oriented and not matrix-entry-oriented
!
call mpi_alltoall(sdsz,1,psb_mpi_mpk_,&
& rvsz,1,psb_mpi_mpk_,icomm,minfo)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mpi_alltoall')
goto 9999
end if
nsnds = count(sdsz /= 0)
nrcvs = count(rvsz /= 0)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done initial alltoall',nsnds,nrcvs
idxs = 0
idxr = 0
counter = 1
Do proc = 0, np-1
bsdindx(proc+1) = idxs
idxs = idxs + sdsz(proc+1)
brvindx(proc+1) = idxr
idxr = idxr + rvsz(proc+1)
Enddo
tsdx = bsdindx
trvx = brvindx
iszr = sum(rvsz)
iszs = sum(sdsz)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Sizes:',&
& ' Send:',sdsz(:),' Receive:',rvsz(:)
call psb_ensure_size(max(iszs,1),iasnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszs,1),jasnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszs,1),valsnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszr,1),iarcv,info)
if (info == psb_success_) call psb_ensure_size(max(iszr,1),jarcv,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='ensure_size')
goto 9999
end if
!
! Now, transpose the matrix, then split between itself
! and the send buffers
!
call acoo%transp()
if (acoo%get_nzeros()/= nzl) then
write(0,*) me,'Something strange upon transpose: ',nzl,acoo%get_nzeros()
end if
nzl = acoo%get_nzeros()
hlstart = p_desc_c%get_local_rows()
nzd = 0
do k = 1, nzl
j = acoo%ia(k)
if (j<=hlstart) then
nzd = nzd + 1
acoo%ia(nzd) = acoo%ia(k)
acoo%ja(nzd) = acoo%ja(k)
acoo%val(nzd) = acoo%val(k)
else
call p_desc_c%indxmap%fnd_halo_owner(j,proc,info)
tsdx(proc+1) = tsdx(proc+1) +1
iasnd(tsdx(proc+1)) = acoo%ia(k)
jasnd(tsdx(proc+1)) = acoo%ja(k)
valsnd(tsdx(proc+1)) = acoo%val(k)
end if
end do
!
! Put halo entries in global numbering
!
call desc_r%indxmap%l2gip(jasnd(1:iszs),info)
call p_desc_c%indxmap%l2gip(iasnd(1:iszs),info)
call acoo%set_nzeros(nzd)
! And exchange data.
! Normally we'll use our SIMPLE A2AV and not MPI, because
! the communication pattern is sparse, so ours is more
! efficient. Using ACOO for the receive buffers.
nzl = acoo%get_nzeros()
call acoo%ensure_size(nzl+iszr)
select case(psb_get_sp_a2av_alg())
case(psb_sp_a2av_smpl_triad_)
call psb_simple_triad_a2av(valsnd,iasnd,jasnd,sdsz,bsdindx,&
& acoo%val(nzl+1:nzl+iszr),iarcv(1:iszr),&
& jarcv(1:iszr),rvsz,brvindx,ictxt,info)
case(psb_sp_a2av_smpl_v_)
call psb_simple_a2av(valsnd,sdsz,bsdindx,&
& acoo%val(nzl+1:nzl+iszr),rvsz,brvindx,ictxt,info)
if (info == psb_success_) call psb_simple_a2av(iasnd,sdsz,bsdindx,&
& iarcv(1:iszr),rvsz,brvindx,ictxt,info)
if (info == psb_success_) call psb_simple_a2av(jasnd,sdsz,bsdindx,&
& jarcv(1:iszr),rvsz,brvindx,ictxt,info)
case(psb_sp_a2av_mpi_)
call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_r_dpk_,&
& acoo%val(nzl+1:nzl+iszr),rvsz,brvindx,psb_mpi_r_dpk_,icomm,minfo)
if (minfo == mpi_success) &
& call mpi_alltoallv(iasnd,sdsz,bsdindx,psb_mpi_lpk_,&
& iarcv(1:iszr),rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo == mpi_success) &
& call mpi_alltoallv(jasnd,sdsz,bsdindx,psb_mpi_lpk_,&
& jarcv(1:iszr),rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo /= mpi_success) info = minfo
case default
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='wrong A2AV alg selector')
goto 9999
end select
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='alltoallv')
goto 9999
end if
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done alltoallv'
if (present(desc_rx)) then
!
! Extend the appropriate descriptor; started as R but on
! transpose it now describes C
!
call desc_r%clone(desc_rx,info)
call psb_cd_reinit(desc_rx,info)
!
! Take out non-local rows
!
call psb_glob_to_loc(iarcv(1:iszr),p_desc_c,info,iact='I',owned=.true.)
lszr = iszr
call psb_coo_clean_negidx_inner(lszr,iarcv(1:iszr),jarcv(1:iszr),&
& acoo%val(nzl+1:nzl+iszr),nzt,info)
call desc_rx%g2lip_ins(jarcv(1:nzt),info)
call psb_cdasb(desc_rx,info)
acoo%ia(nzl+1:nzl+nzt) = iarcv(1:nzt)
acoo%ja(nzl+1:nzl+nzt) = jarcv(1:nzt)
nzl = nzl + nzt
call acoo%set_nzeros(nzl)
nzl = acoo%get_nzeros()
call acoo%set_sorted(.false.)
!
! Insert to extend descriptor
!
call acoo%set_nrows(p_desc_c%get_local_rows())
call acoo%set_ncols(desc_rx%get_local_cols())
!write(0,*) me,' Trans RX ',acoo%get_nrows(),acoo%get_ncols(),acoo%get_nzeros()
else else
call ain%mv_to_lcoo(atcoo,info) !
!
! Take out non-local rows
!
call psb_glob_to_loc(iarcv(1:iszr),p_desc_c,info,iact='I',owned=.true.)
call psb_glob_to_loc(jarcv(1:iszr),desc_r,info,iact='I')
lszr = iszr
call psb_coo_clean_negidx_inner(lszr,iarcv(1:iszr),jarcv(1:iszr),&
& acoo%val(nzl+1:nzl+iszr),nzt,info)
acoo%ia(nzl+1:nzl+nzt) = iarcv(1:nzt)
acoo%ja(nzl+1:nzl+nzt) = jarcv(1:nzt)
nzl = nzl + nzt
call acoo%set_nzeros(nzl)
nzl = acoo%get_nzeros()
call acoo%set_sorted(.false.)
call acoo%set_nrows(p_desc_c%get_local_rows())
call acoo%set_ncols(desc_r%get_local_cols())
!write(0,*) me,' Trans R- ',acoo%get_nrows(),acoo%get_ncols(),acoo%get_nzeros()
end if end if
if (info == 0) call psb_glob_transpose(atcoo,desc_r,info,desc_c=desc_c,desc_rx=desc_rx)
!!$ write(0,*) me,' Sanity check after rx%g2l :',count(acoo%ja(1:nzl)<0)
call acoo%fix(info)
nzl = acoo%get_nzeros()
if (present(atrans)) then if (present(atrans)) then
call atrans%mv_from_lcoo(atcoo,info) call atrans%mv_from_coo(acoo,info)
else else
call ain%mv_from_lcoo(atcoo,info) call ain%mv_from_coo(acoo,info)
end if end if
Deallocate(brvindx,bsdindx,rvsz,sdsz,&
& iasnd,jasnd,valsnd,&
& iarcv,jarcv,stat=info)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done'
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(ictxt,err_act)
return
end subroutine psb_d_coo_glob_transpose end subroutine psb_d_coo_glob_transpose
subroutine psb_d_simple_glob_transpose_ip(ain,desc_a,info) subroutine psb_d_simple_glob_transpose_ip(ain,desc_a,info)

@ -113,7 +113,6 @@ subroutine psb_ls_coo_glob_transpose(ain,desc_r,info,atrans,desc_c,desc_rx)
integer(psb_ipk_) :: ictxt, np,me integer(psb_ipk_) :: ictxt, np,me
integer(psb_ipk_) :: counter,proc, err_act, j integer(psb_ipk_) :: counter,proc, err_act, j
integer(psb_lpk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,& integer(psb_lpk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,&
& irmin,icmin,irmax,icmax,&
& l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzt, nzd & l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzt, nzd
integer(psb_mpk_) :: icomm, minfo integer(psb_mpk_) :: icomm, minfo
integer(psb_mpk_), allocatable :: brvindx(:), & integer(psb_mpk_), allocatable :: brvindx(:), &
@ -407,19 +406,293 @@ subroutine psb_s_coo_glob_transpose(ain,desc_r,info,atrans,desc_c,desc_rx)
type(psb_desc_type), intent(out), optional :: desc_rx type(psb_desc_type), intent(out), optional :: desc_rx
integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), intent(out) :: info
type(psb_ls_coo_sparse_mat) :: atcoo integer(psb_ipk_) :: ictxt, np,me
integer(psb_ipk_) :: counter,proc, err_act, j
integer(psb_ipk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,&
& l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzd
integer(psb_lpk_) :: nzt, lszr
integer(psb_mpk_) :: icomm, minfo
integer(psb_mpk_), allocatable :: brvindx(:), &
& rvsz(:), bsdindx(:), sdsz(:), tsdx(:), trvx(:)
integer(psb_ipk_), allocatable :: halo_owner(:)
integer(psb_lpk_), allocatable :: iasnd(:), jasnd(:),iarcv(:),jarcv(:)
real(psb_spk_), allocatable :: valsnd(:)
type(psb_s_coo_sparse_mat), allocatable :: acoo
logical :: rowcnv_,colcnv_,rowscale_,colscale_,outcol_glob_
type(psb_desc_type), pointer :: p_desc_c
character(len=5) :: outfmt_
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
if(psb_get_errstatus() /= 0) return
info=psb_success_
name='mld_glob_transpose'
call psb_erractionsave(err_act)
if (psb_errstatus_fatal()) then
info = psb_err_internal_error_ ; goto 9999
end if
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
ictxt = desc_r%get_context()
icomm = desc_r%get_mpic()
Call psb_info(ictxt, me, np)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),': Start'
if (present(desc_c)) then
p_desc_c => desc_c
else
p_desc_c => desc_r
end if
Allocate(brvindx(np+1),&
& rvsz(np),sdsz(np),bsdindx(np+1), acoo,stat=info)
if (info /= psb_success_) then
info=psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
sdsz(:)=0
rvsz(:)=0
l1 = 0
brvindx(1) = 0
bsdindx(1) = 0
counter=1
idx = 0
idxs = 0
idxr = 0
if (present(atrans)) then if (present(atrans)) then
call ain%cp_to_lcoo(atcoo,info) call ain%cp_to_coo(acoo,info)
else
call ain%mv_to_coo(acoo,info)
end if
!
! Compute number of entries in the
! halo part, sorted by destination process
!
nr = desc_r%get_local_rows()
nc = p_desc_c%get_local_cols()
nzl = acoo%get_nzeros()
hlstart = p_desc_c%get_local_rows()
do k=1, nzl
j = acoo%ja(k)
if (j > hlstart) then
call p_desc_c%indxmap%fnd_halo_owner(j,proc,info)
sdsz(proc+1) = sdsz(proc+1) +1
end if
end do
!
! Exchange sizes, so as to know sends/receives.
! This is different from the halo exchange because the
! number of entries was not precomputed in the descriptor,
! which was vector-oriented and not matrix-entry-oriented
!
call mpi_alltoall(sdsz,1,psb_mpi_mpk_,&
& rvsz,1,psb_mpi_mpk_,icomm,minfo)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mpi_alltoall')
goto 9999
end if
nsnds = count(sdsz /= 0)
nrcvs = count(rvsz /= 0)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done initial alltoall',nsnds,nrcvs
idxs = 0
idxr = 0
counter = 1
Do proc = 0, np-1
bsdindx(proc+1) = idxs
idxs = idxs + sdsz(proc+1)
brvindx(proc+1) = idxr
idxr = idxr + rvsz(proc+1)
Enddo
tsdx = bsdindx
trvx = brvindx
iszr = sum(rvsz)
iszs = sum(sdsz)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Sizes:',&
& ' Send:',sdsz(:),' Receive:',rvsz(:)
call psb_ensure_size(max(iszs,1),iasnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszs,1),jasnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszs,1),valsnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszr,1),iarcv,info)
if (info == psb_success_) call psb_ensure_size(max(iszr,1),jarcv,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='ensure_size')
goto 9999
end if
!
! Now, transpose the matrix, then split between itself
! and the send buffers
!
call acoo%transp()
if (acoo%get_nzeros()/= nzl) then
write(0,*) me,'Something strange upon transpose: ',nzl,acoo%get_nzeros()
end if
nzl = acoo%get_nzeros()
hlstart = p_desc_c%get_local_rows()
nzd = 0
do k = 1, nzl
j = acoo%ia(k)
if (j<=hlstart) then
nzd = nzd + 1
acoo%ia(nzd) = acoo%ia(k)
acoo%ja(nzd) = acoo%ja(k)
acoo%val(nzd) = acoo%val(k)
else
call p_desc_c%indxmap%fnd_halo_owner(j,proc,info)
tsdx(proc+1) = tsdx(proc+1) +1
iasnd(tsdx(proc+1)) = acoo%ia(k)
jasnd(tsdx(proc+1)) = acoo%ja(k)
valsnd(tsdx(proc+1)) = acoo%val(k)
end if
end do
!
! Put halo entries in global numbering
!
call desc_r%indxmap%l2gip(jasnd(1:iszs),info)
call p_desc_c%indxmap%l2gip(iasnd(1:iszs),info)
call acoo%set_nzeros(nzd)
! And exchange data.
! Normally we'll use our SIMPLE A2AV and not MPI, because
! the communication pattern is sparse, so ours is more
! efficient. Using ACOO for the receive buffers.
nzl = acoo%get_nzeros()
call acoo%ensure_size(nzl+iszr)
select case(psb_get_sp_a2av_alg())
case(psb_sp_a2av_smpl_triad_)
call psb_simple_triad_a2av(valsnd,iasnd,jasnd,sdsz,bsdindx,&
& acoo%val(nzl+1:nzl+iszr),iarcv(1:iszr),&
& jarcv(1:iszr),rvsz,brvindx,ictxt,info)
case(psb_sp_a2av_smpl_v_)
call psb_simple_a2av(valsnd,sdsz,bsdindx,&
& acoo%val(nzl+1:nzl+iszr),rvsz,brvindx,ictxt,info)
if (info == psb_success_) call psb_simple_a2av(iasnd,sdsz,bsdindx,&
& iarcv(1:iszr),rvsz,brvindx,ictxt,info)
if (info == psb_success_) call psb_simple_a2av(jasnd,sdsz,bsdindx,&
& jarcv(1:iszr),rvsz,brvindx,ictxt,info)
case(psb_sp_a2av_mpi_)
call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_r_spk_,&
& acoo%val(nzl+1:nzl+iszr),rvsz,brvindx,psb_mpi_r_spk_,icomm,minfo)
if (minfo == mpi_success) &
& call mpi_alltoallv(iasnd,sdsz,bsdindx,psb_mpi_lpk_,&
& iarcv(1:iszr),rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo == mpi_success) &
& call mpi_alltoallv(jasnd,sdsz,bsdindx,psb_mpi_lpk_,&
& jarcv(1:iszr),rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo /= mpi_success) info = minfo
case default
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='wrong A2AV alg selector')
goto 9999
end select
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='alltoallv')
goto 9999
end if
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done alltoallv'
if (present(desc_rx)) then
!
! Extend the appropriate descriptor; started as R but on
! transpose it now describes C
!
call desc_r%clone(desc_rx,info)
call psb_cd_reinit(desc_rx,info)
!
! Take out non-local rows
!
call psb_glob_to_loc(iarcv(1:iszr),p_desc_c,info,iact='I',owned=.true.)
lszr = iszr
call psb_coo_clean_negidx_inner(lszr,iarcv(1:iszr),jarcv(1:iszr),&
& acoo%val(nzl+1:nzl+iszr),nzt,info)
call desc_rx%g2lip_ins(jarcv(1:nzt),info)
call psb_cdasb(desc_rx,info)
acoo%ia(nzl+1:nzl+nzt) = iarcv(1:nzt)
acoo%ja(nzl+1:nzl+nzt) = jarcv(1:nzt)
nzl = nzl + nzt
call acoo%set_nzeros(nzl)
nzl = acoo%get_nzeros()
call acoo%set_sorted(.false.)
!
! Insert to extend descriptor
!
call acoo%set_nrows(p_desc_c%get_local_rows())
call acoo%set_ncols(desc_rx%get_local_cols())
!write(0,*) me,' Trans RX ',acoo%get_nrows(),acoo%get_ncols(),acoo%get_nzeros()
else else
call ain%mv_to_lcoo(atcoo,info) !
!
! Take out non-local rows
!
call psb_glob_to_loc(iarcv(1:iszr),p_desc_c,info,iact='I',owned=.true.)
call psb_glob_to_loc(jarcv(1:iszr),desc_r,info,iact='I')
lszr = iszr
call psb_coo_clean_negidx_inner(lszr,iarcv(1:iszr),jarcv(1:iszr),&
& acoo%val(nzl+1:nzl+iszr),nzt,info)
acoo%ia(nzl+1:nzl+nzt) = iarcv(1:nzt)
acoo%ja(nzl+1:nzl+nzt) = jarcv(1:nzt)
nzl = nzl + nzt
call acoo%set_nzeros(nzl)
nzl = acoo%get_nzeros()
call acoo%set_sorted(.false.)
call acoo%set_nrows(p_desc_c%get_local_rows())
call acoo%set_ncols(desc_r%get_local_cols())
!write(0,*) me,' Trans R- ',acoo%get_nrows(),acoo%get_ncols(),acoo%get_nzeros()
end if end if
if (info == 0) call psb_glob_transpose(atcoo,desc_r,info,desc_c=desc_c,desc_rx=desc_rx)
!!$ write(0,*) me,' Sanity check after rx%g2l :',count(acoo%ja(1:nzl)<0)
call acoo%fix(info)
nzl = acoo%get_nzeros()
if (present(atrans)) then if (present(atrans)) then
call atrans%mv_from_lcoo(atcoo,info) call atrans%mv_from_coo(acoo,info)
else else
call ain%mv_from_lcoo(atcoo,info) call ain%mv_from_coo(acoo,info)
end if end if
Deallocate(brvindx,bsdindx,rvsz,sdsz,&
& iasnd,jasnd,valsnd,&
& iarcv,jarcv,stat=info)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done'
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(ictxt,err_act)
return
end subroutine psb_s_coo_glob_transpose end subroutine psb_s_coo_glob_transpose
subroutine psb_s_simple_glob_transpose_ip(ain,desc_a,info) subroutine psb_s_simple_glob_transpose_ip(ain,desc_a,info)

@ -113,7 +113,6 @@ subroutine psb_lz_coo_glob_transpose(ain,desc_r,info,atrans,desc_c,desc_rx)
integer(psb_ipk_) :: ictxt, np,me integer(psb_ipk_) :: ictxt, np,me
integer(psb_ipk_) :: counter,proc, err_act, j integer(psb_ipk_) :: counter,proc, err_act, j
integer(psb_lpk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,& integer(psb_lpk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,&
& irmin,icmin,irmax,icmax,&
& l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzt, nzd & l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzt, nzd
integer(psb_mpk_) :: icomm, minfo integer(psb_mpk_) :: icomm, minfo
integer(psb_mpk_), allocatable :: brvindx(:), & integer(psb_mpk_), allocatable :: brvindx(:), &
@ -407,19 +406,293 @@ subroutine psb_z_coo_glob_transpose(ain,desc_r,info,atrans,desc_c,desc_rx)
type(psb_desc_type), intent(out), optional :: desc_rx type(psb_desc_type), intent(out), optional :: desc_rx
integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), intent(out) :: info
type(psb_lz_coo_sparse_mat) :: atcoo integer(psb_ipk_) :: ictxt, np,me
integer(psb_ipk_) :: counter,proc, err_act, j
integer(psb_ipk_) :: i, k, idx, r, ipx,mat_recv, iszs, iszr,idxs,idxr,nz,&
& l1, nsnds, nrcvs, nr,nc,nzl, hlstart, nzd
integer(psb_lpk_) :: nzt, lszr
integer(psb_mpk_) :: icomm, minfo
integer(psb_mpk_), allocatable :: brvindx(:), &
& rvsz(:), bsdindx(:), sdsz(:), tsdx(:), trvx(:)
integer(psb_ipk_), allocatable :: halo_owner(:)
integer(psb_lpk_), allocatable :: iasnd(:), jasnd(:),iarcv(:),jarcv(:)
complex(psb_dpk_), allocatable :: valsnd(:)
type(psb_z_coo_sparse_mat), allocatable :: acoo
logical :: rowcnv_,colcnv_,rowscale_,colscale_,outcol_glob_
type(psb_desc_type), pointer :: p_desc_c
character(len=5) :: outfmt_
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
if(psb_get_errstatus() /= 0) return
info=psb_success_
name='mld_glob_transpose'
call psb_erractionsave(err_act)
if (psb_errstatus_fatal()) then
info = psb_err_internal_error_ ; goto 9999
end if
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
ictxt = desc_r%get_context()
icomm = desc_r%get_mpic()
Call psb_info(ictxt, me, np)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),': Start'
if (present(desc_c)) then
p_desc_c => desc_c
else
p_desc_c => desc_r
end if
Allocate(brvindx(np+1),&
& rvsz(np),sdsz(np),bsdindx(np+1), acoo,stat=info)
if (info /= psb_success_) then
info=psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
sdsz(:)=0
rvsz(:)=0
l1 = 0
brvindx(1) = 0
bsdindx(1) = 0
counter=1
idx = 0
idxs = 0
idxr = 0
if (present(atrans)) then if (present(atrans)) then
call ain%cp_to_lcoo(atcoo,info) call ain%cp_to_coo(acoo,info)
else
call ain%mv_to_coo(acoo,info)
end if
!
! Compute number of entries in the
! halo part, sorted by destination process
!
nr = desc_r%get_local_rows()
nc = p_desc_c%get_local_cols()
nzl = acoo%get_nzeros()
hlstart = p_desc_c%get_local_rows()
do k=1, nzl
j = acoo%ja(k)
if (j > hlstart) then
call p_desc_c%indxmap%fnd_halo_owner(j,proc,info)
sdsz(proc+1) = sdsz(proc+1) +1
end if
end do
!
! Exchange sizes, so as to know sends/receives.
! This is different from the halo exchange because the
! number of entries was not precomputed in the descriptor,
! which was vector-oriented and not matrix-entry-oriented
!
call mpi_alltoall(sdsz,1,psb_mpi_mpk_,&
& rvsz,1,psb_mpi_mpk_,icomm,minfo)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mpi_alltoall')
goto 9999
end if
nsnds = count(sdsz /= 0)
nrcvs = count(rvsz /= 0)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done initial alltoall',nsnds,nrcvs
idxs = 0
idxr = 0
counter = 1
Do proc = 0, np-1
bsdindx(proc+1) = idxs
idxs = idxs + sdsz(proc+1)
brvindx(proc+1) = idxr
idxr = idxr + rvsz(proc+1)
Enddo
tsdx = bsdindx
trvx = brvindx
iszr = sum(rvsz)
iszs = sum(sdsz)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Sizes:',&
& ' Send:',sdsz(:),' Receive:',rvsz(:)
call psb_ensure_size(max(iszs,1),iasnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszs,1),jasnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszs,1),valsnd,info)
if (info == psb_success_) call psb_ensure_size(max(iszr,1),iarcv,info)
if (info == psb_success_) call psb_ensure_size(max(iszr,1),jarcv,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='ensure_size')
goto 9999
end if
!
! Now, transpose the matrix, then split between itself
! and the send buffers
!
call acoo%transp()
if (acoo%get_nzeros()/= nzl) then
write(0,*) me,'Something strange upon transpose: ',nzl,acoo%get_nzeros()
end if
nzl = acoo%get_nzeros()
hlstart = p_desc_c%get_local_rows()
nzd = 0
do k = 1, nzl
j = acoo%ia(k)
if (j<=hlstart) then
nzd = nzd + 1
acoo%ia(nzd) = acoo%ia(k)
acoo%ja(nzd) = acoo%ja(k)
acoo%val(nzd) = acoo%val(k)
else
call p_desc_c%indxmap%fnd_halo_owner(j,proc,info)
tsdx(proc+1) = tsdx(proc+1) +1
iasnd(tsdx(proc+1)) = acoo%ia(k)
jasnd(tsdx(proc+1)) = acoo%ja(k)
valsnd(tsdx(proc+1)) = acoo%val(k)
end if
end do
!
! Put halo entries in global numbering
!
call desc_r%indxmap%l2gip(jasnd(1:iszs),info)
call p_desc_c%indxmap%l2gip(iasnd(1:iszs),info)
call acoo%set_nzeros(nzd)
! And exchange data.
! Normally we'll use our SIMPLE A2AV and not MPI, because
! the communication pattern is sparse, so ours is more
! efficient. Using ACOO for the receive buffers.
nzl = acoo%get_nzeros()
call acoo%ensure_size(nzl+iszr)
select case(psb_get_sp_a2av_alg())
case(psb_sp_a2av_smpl_triad_)
call psb_simple_triad_a2av(valsnd,iasnd,jasnd,sdsz,bsdindx,&
& acoo%val(nzl+1:nzl+iszr),iarcv(1:iszr),&
& jarcv(1:iszr),rvsz,brvindx,ictxt,info)
case(psb_sp_a2av_smpl_v_)
call psb_simple_a2av(valsnd,sdsz,bsdindx,&
& acoo%val(nzl+1:nzl+iszr),rvsz,brvindx,ictxt,info)
if (info == psb_success_) call psb_simple_a2av(iasnd,sdsz,bsdindx,&
& iarcv(1:iszr),rvsz,brvindx,ictxt,info)
if (info == psb_success_) call psb_simple_a2av(jasnd,sdsz,bsdindx,&
& jarcv(1:iszr),rvsz,brvindx,ictxt,info)
case(psb_sp_a2av_mpi_)
call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_c_dpk_,&
& acoo%val(nzl+1:nzl+iszr),rvsz,brvindx,psb_mpi_c_dpk_,icomm,minfo)
if (minfo == mpi_success) &
& call mpi_alltoallv(iasnd,sdsz,bsdindx,psb_mpi_lpk_,&
& iarcv(1:iszr),rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo == mpi_success) &
& call mpi_alltoallv(jasnd,sdsz,bsdindx,psb_mpi_lpk_,&
& jarcv(1:iszr),rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo /= mpi_success) info = minfo
case default
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='wrong A2AV alg selector')
goto 9999
end select
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='alltoallv')
goto 9999
end if
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done alltoallv'
if (present(desc_rx)) then
!
! Extend the appropriate descriptor; started as R but on
! transpose it now describes C
!
call desc_r%clone(desc_rx,info)
call psb_cd_reinit(desc_rx,info)
!
! Take out non-local rows
!
call psb_glob_to_loc(iarcv(1:iszr),p_desc_c,info,iact='I',owned=.true.)
lszr = iszr
call psb_coo_clean_negidx_inner(lszr,iarcv(1:iszr),jarcv(1:iszr),&
& acoo%val(nzl+1:nzl+iszr),nzt,info)
call desc_rx%g2lip_ins(jarcv(1:nzt),info)
call psb_cdasb(desc_rx,info)
acoo%ia(nzl+1:nzl+nzt) = iarcv(1:nzt)
acoo%ja(nzl+1:nzl+nzt) = jarcv(1:nzt)
nzl = nzl + nzt
call acoo%set_nzeros(nzl)
nzl = acoo%get_nzeros()
call acoo%set_sorted(.false.)
!
! Insert to extend descriptor
!
call acoo%set_nrows(p_desc_c%get_local_rows())
call acoo%set_ncols(desc_rx%get_local_cols())
!write(0,*) me,' Trans RX ',acoo%get_nrows(),acoo%get_ncols(),acoo%get_nzeros()
else else
call ain%mv_to_lcoo(atcoo,info) !
!
! Take out non-local rows
!
call psb_glob_to_loc(iarcv(1:iszr),p_desc_c,info,iact='I',owned=.true.)
call psb_glob_to_loc(jarcv(1:iszr),desc_r,info,iact='I')
lszr = iszr
call psb_coo_clean_negidx_inner(lszr,iarcv(1:iszr),jarcv(1:iszr),&
& acoo%val(nzl+1:nzl+iszr),nzt,info)
acoo%ia(nzl+1:nzl+nzt) = iarcv(1:nzt)
acoo%ja(nzl+1:nzl+nzt) = jarcv(1:nzt)
nzl = nzl + nzt
call acoo%set_nzeros(nzl)
nzl = acoo%get_nzeros()
call acoo%set_sorted(.false.)
call acoo%set_nrows(p_desc_c%get_local_rows())
call acoo%set_ncols(desc_r%get_local_cols())
!write(0,*) me,' Trans R- ',acoo%get_nrows(),acoo%get_ncols(),acoo%get_nzeros()
end if end if
if (info == 0) call psb_glob_transpose(atcoo,desc_r,info,desc_c=desc_c,desc_rx=desc_rx)
!!$ write(0,*) me,' Sanity check after rx%g2l :',count(acoo%ja(1:nzl)<0)
call acoo%fix(info)
nzl = acoo%get_nzeros()
if (present(atrans)) then if (present(atrans)) then
call atrans%mv_from_lcoo(atcoo,info) call atrans%mv_from_coo(acoo,info)
else else
call ain%mv_from_lcoo(atcoo,info) call ain%mv_from_coo(acoo,info)
end if end if
Deallocate(brvindx,bsdindx,rvsz,sdsz,&
& iasnd,jasnd,valsnd,&
& iarcv,jarcv,stat=info)
if (debug_level >= psb_debug_outer_)&
& write(debug_unit,*) me,' ',trim(name),': Done'
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(ictxt,err_act)
return
end subroutine psb_z_coo_glob_transpose end subroutine psb_z_coo_glob_transpose
subroutine psb_z_simple_glob_transpose_ip(ain,desc_a,info) subroutine psb_z_simple_glob_transpose_ip(ain,desc_a,info)

Loading…
Cancel
Save