diff --git a/base/modules/Makefile b/base/modules/Makefile index 9951bd89..cfc3c9a4 100644 --- a/base/modules/Makefile +++ b/base/modules/Makefile @@ -358,8 +358,8 @@ tools/psb_c_tools_mod.o tools/psb_z_tools_mod.o tools/psb_m_tools_a_mod.o tools/ tools/psb_s_tools_a_mod.o tools/psb_d_tools_a_mod.o\ tools/psb_c_tools_a_mod.o tools/psb_z_tools_a_mod.o: desc/psb_desc_mod.o psi_mod.o serial/psb_mat_mod.o -tools/psb_i_tools_mod.o: serial/psb_i_vect_mod.o -tools/psb_l_tools_mod.o: serial/psb_l_vect_mod.o +tools/psb_i_tools_mod.o: serial/psb_i_vect_mod.o tools/psb_m_tools_a_mod.o tools/psb_e_tools_a_mod.o +tools/psb_l_tools_mod.o: serial/psb_l_vect_mod.o tools/psb_m_tools_a_mod.o tools/psb_e_tools_a_mod.o tools/psb_s_tools_mod.o: serial/psb_s_vect_mod.o tools/psb_d_tools_mod.o: serial/psb_d_vect_mod.o tools/psb_c_tools_mod.o: serial/psb_c_vect_mod.o diff --git a/base/modules/tools/psb_c_tools_a_mod.f90 b/base/modules/tools/psb_c_tools_a_mod.f90 index 6c864ead..8dd592b3 100644 --- a/base/modules/tools/psb_c_tools_a_mod.f90 +++ b/base/modules/tools/psb_c_tools_a_mod.f90 @@ -116,4 +116,19 @@ Module psb_c_tools_a_mod end subroutine psb_cinsvi end interface + + interface psb_remote_vect + subroutine psb_c_remote_vect(n,v,iv,desc_a,x,ix, info) + import + implicit none + integer(psb_ipk_), intent(in) :: n + complex(psb_spk_), intent(in) :: v(:) + integer(psb_lpk_), intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + complex(psb_spk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + end subroutine psb_c_remote_vect + end interface psb_remote_vect + end module psb_c_tools_a_mod diff --git a/base/modules/tools/psb_c_tools_mod.F90 b/base/modules/tools/psb_c_tools_mod.F90 index 78ec8518..2de8f906 100644 --- a/base/modules/tools/psb_c_tools_mod.F90 +++ b/base/modules/tools/psb_c_tools_mod.F90 @@ -171,7 +171,7 @@ Module psb_c_tools_mod logical, intent(in), optional :: local end subroutine psb_cins_multivect end interface - + interface psb_cdbldext Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info,extype) import @@ -262,20 +262,6 @@ Module psb_c_tools_mod end subroutine psb_cspasb end interface - interface psb_remote_vect - subroutine psb_c_remote_vect(n,v,iv,desc_a,x,ix, info) - import - implicit none - integer(psb_ipk_), intent(in) :: n - complex(psb_spk_), intent(in) :: v(:) - integer(psb_lpk_), intent(in) :: iv(:) - type(psb_desc_type),intent(in) :: desc_a - complex(psb_spk_), allocatable, intent(out) :: x(:) - integer(psb_lpk_), allocatable, intent(out) :: ix(:) - integer(psb_ipk_), intent(out) :: info - end subroutine psb_c_remote_vect - end interface psb_remote_vect - interface psb_remote_mat subroutine psb_lc_remote_mat(a,desc_a,b, info) import diff --git a/base/modules/tools/psb_d_tools_a_mod.f90 b/base/modules/tools/psb_d_tools_a_mod.f90 index 1ce3d774..638953b4 100644 --- a/base/modules/tools/psb_d_tools_a_mod.f90 +++ b/base/modules/tools/psb_d_tools_a_mod.f90 @@ -116,4 +116,19 @@ Module psb_d_tools_a_mod end subroutine psb_dinsvi end interface + + interface psb_remote_vect + subroutine psb_d_remote_vect(n,v,iv,desc_a,x,ix, info) + import + implicit none + integer(psb_ipk_), intent(in) :: n + real(psb_dpk_), intent(in) :: v(:) + integer(psb_lpk_), intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + real(psb_dpk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + end subroutine psb_d_remote_vect + end interface psb_remote_vect + end module psb_d_tools_a_mod diff --git a/base/modules/tools/psb_d_tools_mod.F90 b/base/modules/tools/psb_d_tools_mod.F90 index 30510123..30e45d53 100644 --- a/base/modules/tools/psb_d_tools_mod.F90 +++ b/base/modules/tools/psb_d_tools_mod.F90 @@ -171,7 +171,7 @@ Module psb_d_tools_mod logical, intent(in), optional :: local end subroutine psb_dins_multivect end interface - + interface psb_cdbldext Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info,extype) import @@ -262,20 +262,6 @@ Module psb_d_tools_mod end subroutine psb_dspasb end interface - interface psb_remote_vect - subroutine psb_d_remote_vect(n,v,iv,desc_a,x,ix, info) - import - implicit none - integer(psb_ipk_), intent(in) :: n - real(psb_dpk_), intent(in) :: v(:) - integer(psb_lpk_), intent(in) :: iv(:) - type(psb_desc_type),intent(in) :: desc_a - real(psb_dpk_), allocatable, intent(out) :: x(:) - integer(psb_lpk_), allocatable, intent(out) :: ix(:) - integer(psb_ipk_), intent(out) :: info - end subroutine psb_d_remote_vect - end interface psb_remote_vect - interface psb_remote_mat subroutine psb_ld_remote_mat(a,desc_a,b, info) import diff --git a/base/modules/tools/psb_e_tools_a_mod.f90 b/base/modules/tools/psb_e_tools_a_mod.f90 index bce8cb40..f8a27cb5 100644 --- a/base/modules/tools/psb_e_tools_a_mod.f90 +++ b/base/modules/tools/psb_e_tools_a_mod.f90 @@ -116,4 +116,19 @@ Module psb_e_tools_a_mod end subroutine psb_einsvi end interface + + interface psb_remote_vect + subroutine psb_e_remote_vect(n,v,iv,desc_a,x,ix, info) + import + implicit none + integer(psb_ipk_), intent(in) :: n + integer(psb_epk_), intent(in) :: v(:) + integer(psb_lpk_), intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + integer(psb_epk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + end subroutine psb_e_remote_vect + end interface psb_remote_vect + end module psb_e_tools_a_mod diff --git a/base/modules/tools/psb_i2_tools_a_mod.f90 b/base/modules/tools/psb_i2_tools_a_mod.f90 index 860a55b1..b8d52bb4 100644 --- a/base/modules/tools/psb_i2_tools_a_mod.f90 +++ b/base/modules/tools/psb_i2_tools_a_mod.f90 @@ -116,4 +116,19 @@ Module psb_i2_tools_a_mod end subroutine psb_i2insvi end interface + + interface psb_remote_vect + subroutine psb_i2_remote_vect(n,v,iv,desc_a,x,ix, info) + import + implicit none + integer(psb_ipk_), intent(in) :: n + integer(psb_i2pk_), intent(in) :: v(:) + integer(psb_lpk_), intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + integer(psb_i2pk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + end subroutine psb_i2_remote_vect + end interface psb_remote_vect + end module psb_i2_tools_a_mod diff --git a/base/modules/tools/psb_i_tools_mod.F90 b/base/modules/tools/psb_i_tools_mod.F90 index 94fb04bc..1c207fac 100644 --- a/base/modules/tools/psb_i_tools_mod.F90 +++ b/base/modules/tools/psb_i_tools_mod.F90 @@ -32,6 +32,8 @@ Module psb_i_tools_mod use psb_desc_mod, only : psb_desc_type, psb_ipk_, psb_lpk_, psb_success_ use psb_i_vect_mod, only : psb_i_base_vect_type, psb_i_vect_type + use psb_m_tools_a_mod + use psb_e_tools_a_mod use psb_l_vect_mod, only : psb_l_vect_type use psb_i_multivect_mod, only : psb_i_base_multivect_type, psb_i_multivect_type use psi_mod, only : psb_snd, psb_rcv ! Needed only for psb_getelem @@ -168,5 +170,5 @@ Module psb_i_tools_mod logical, intent(in), optional :: local end subroutine psb_iins_multivect end interface - + end module psb_i_tools_mod diff --git a/base/modules/tools/psb_l_tools_mod.F90 b/base/modules/tools/psb_l_tools_mod.F90 index 61840af6..058403d6 100644 --- a/base/modules/tools/psb_l_tools_mod.F90 +++ b/base/modules/tools/psb_l_tools_mod.F90 @@ -32,6 +32,8 @@ Module psb_l_tools_mod use psb_desc_mod, only : psb_desc_type, psb_ipk_, psb_lpk_, psb_success_ use psb_l_vect_mod, only : psb_l_base_vect_type, psb_l_vect_type + use psb_m_tools_a_mod + use psb_e_tools_a_mod ! use psb_i_vect_mod, only : psb_i_vect_type use psb_l_multivect_mod, only : psb_l_base_multivect_type, psb_l_multivect_type use psi_mod, only : psb_snd, psb_rcv ! Needed only for psb_getelem @@ -168,5 +170,5 @@ Module psb_l_tools_mod logical, intent(in), optional :: local end subroutine psb_lins_multivect end interface - + end module psb_l_tools_mod diff --git a/base/modules/tools/psb_m_tools_a_mod.f90 b/base/modules/tools/psb_m_tools_a_mod.f90 index a5dfdd72..a1de6caa 100644 --- a/base/modules/tools/psb_m_tools_a_mod.f90 +++ b/base/modules/tools/psb_m_tools_a_mod.f90 @@ -116,4 +116,19 @@ Module psb_m_tools_a_mod end subroutine psb_minsvi end interface + + interface psb_remote_vect + subroutine psb_m_remote_vect(n,v,iv,desc_a,x,ix, info) + import + implicit none + integer(psb_ipk_), intent(in) :: n + integer(psb_mpk_), intent(in) :: v(:) + integer(psb_lpk_), intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + integer(psb_mpk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + end subroutine psb_m_remote_vect + end interface psb_remote_vect + end module psb_m_tools_a_mod diff --git a/base/modules/tools/psb_s_tools_a_mod.f90 b/base/modules/tools/psb_s_tools_a_mod.f90 index 32f445cb..85a2b029 100644 --- a/base/modules/tools/psb_s_tools_a_mod.f90 +++ b/base/modules/tools/psb_s_tools_a_mod.f90 @@ -116,4 +116,19 @@ Module psb_s_tools_a_mod end subroutine psb_sinsvi end interface + + interface psb_remote_vect + subroutine psb_s_remote_vect(n,v,iv,desc_a,x,ix, info) + import + implicit none + integer(psb_ipk_), intent(in) :: n + real(psb_spk_), intent(in) :: v(:) + integer(psb_lpk_), intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + real(psb_spk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + end subroutine psb_s_remote_vect + end interface psb_remote_vect + end module psb_s_tools_a_mod diff --git a/base/modules/tools/psb_s_tools_mod.F90 b/base/modules/tools/psb_s_tools_mod.F90 index 72033781..5d2f8d00 100644 --- a/base/modules/tools/psb_s_tools_mod.F90 +++ b/base/modules/tools/psb_s_tools_mod.F90 @@ -171,7 +171,7 @@ Module psb_s_tools_mod logical, intent(in), optional :: local end subroutine psb_sins_multivect end interface - + interface psb_cdbldext Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info,extype) import @@ -262,20 +262,6 @@ Module psb_s_tools_mod end subroutine psb_sspasb end interface - interface psb_remote_vect - subroutine psb_s_remote_vect(n,v,iv,desc_a,x,ix, info) - import - implicit none - integer(psb_ipk_), intent(in) :: n - real(psb_spk_), intent(in) :: v(:) - integer(psb_lpk_), intent(in) :: iv(:) - type(psb_desc_type),intent(in) :: desc_a - real(psb_spk_), allocatable, intent(out) :: x(:) - integer(psb_lpk_), allocatable, intent(out) :: ix(:) - integer(psb_ipk_), intent(out) :: info - end subroutine psb_s_remote_vect - end interface psb_remote_vect - interface psb_remote_mat subroutine psb_ls_remote_mat(a,desc_a,b, info) import diff --git a/base/modules/tools/psb_z_tools_a_mod.f90 b/base/modules/tools/psb_z_tools_a_mod.f90 index 21f7ff0f..4c421f26 100644 --- a/base/modules/tools/psb_z_tools_a_mod.f90 +++ b/base/modules/tools/psb_z_tools_a_mod.f90 @@ -116,4 +116,19 @@ Module psb_z_tools_a_mod end subroutine psb_zinsvi end interface + + interface psb_remote_vect + subroutine psb_z_remote_vect(n,v,iv,desc_a,x,ix, info) + import + implicit none + integer(psb_ipk_), intent(in) :: n + complex(psb_dpk_), intent(in) :: v(:) + integer(psb_lpk_), intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + complex(psb_dpk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + end subroutine psb_z_remote_vect + end interface psb_remote_vect + end module psb_z_tools_a_mod diff --git a/base/modules/tools/psb_z_tools_mod.F90 b/base/modules/tools/psb_z_tools_mod.F90 index c96737b1..9d6bd77b 100644 --- a/base/modules/tools/psb_z_tools_mod.F90 +++ b/base/modules/tools/psb_z_tools_mod.F90 @@ -171,7 +171,7 @@ Module psb_z_tools_mod logical, intent(in), optional :: local end subroutine psb_zins_multivect end interface - + interface psb_cdbldext Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info,extype) import @@ -262,20 +262,6 @@ Module psb_z_tools_mod end subroutine psb_zspasb end interface - interface psb_remote_vect - subroutine psb_z_remote_vect(n,v,iv,desc_a,x,ix, info) - import - implicit none - integer(psb_ipk_), intent(in) :: n - complex(psb_dpk_), intent(in) :: v(:) - integer(psb_lpk_), intent(in) :: iv(:) - type(psb_desc_type),intent(in) :: desc_a - complex(psb_dpk_), allocatable, intent(out) :: x(:) - integer(psb_lpk_), allocatable, intent(out) :: ix(:) - integer(psb_ipk_), intent(out) :: info - end subroutine psb_z_remote_vect - end interface psb_remote_vect - interface psb_remote_mat subroutine psb_lz_remote_mat(a,desc_a,b, info) import diff --git a/base/tools/Makefile b/base/tools/Makefile index ce13caed..968aec18 100644 --- a/base/tools/Makefile +++ b/base/tools/Makefile @@ -31,7 +31,9 @@ FOBJS = psb_cdall.o psb_cdals.o psb_cdalv.o psb_cd_inloc.o psb_cdins.o psb_cdprt MPFOBJS = psb_icdasb.o psb_ssphalo.o psb_dsphalo.o psb_csphalo.o psb_zsphalo.o \ psb_dcdbldext.o psb_zcdbldext.o psb_scdbldext.o psb_ccdbldext.o \ - psb_s_remote_mat.o psb_d_remote_mat.o psb_c_remote_mat.o psb_z_remote_mat.o + psb_s_remote_mat.o psb_d_remote_mat.o psb_c_remote_mat.o psb_z_remote_mat.o \ + psb_s_remote_vect.o psb_d_remote_vect.o psb_c_remote_vect.o psb_z_remote_vect.o \ + psb_e_remote_vect.o psb_m_remote_vect.o LIBDIR=.. INCDIR=.. diff --git a/base/tools/psb_c_remote_mat.F90 b/base/tools/psb_c_remote_mat.F90 index 18c7a91d..ae2eaaf2 100644 --- a/base/tools/psb_c_remote_mat.F90 +++ b/base/tools/psb_c_remote_mat.F90 @@ -274,163 +274,3 @@ Subroutine psb_lc_remote_mat(a,desc_a,b,info) return End Subroutine psb_lc_remote_mat - - -subroutine psb_c_remote_vect(n,v,iv,desc_a,x,ix, info) - use psb_base_mod, psb_protect_name => psb_c_remote_vect - -#ifdef MPI_MOD - use mpi -#endif - Implicit None -#ifdef MPI_H - include 'mpif.h' -#endif - integer(psb_ipk_), intent(in) :: n - complex(psb_spk_), Intent(in) :: v(:) - integer(psb_lpk_), Intent(in) :: iv(:) - type(psb_desc_type),intent(in) :: desc_a - complex(psb_spk_), allocatable, intent(out) :: x(:) - integer(psb_lpk_), allocatable, intent(out) :: ix(:) - integer(psb_ipk_), intent(out) :: info - ! ...local scalars.... - type(psb_ctxt_type) :: ctxt - integer(psb_ipk_) :: np,me - integer(psb_ipk_) :: counter, proc, i, & - & j, idxs,idxr, k, iszs, iszr - integer(psb_ipk_) :: nrcvs, nsnds - integer(psb_mpk_) :: icomm, minfo - integer(psb_mpk_), allocatable :: brvindx(:), & - & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) - integer(psb_lpk_), allocatable :: lsnd(:) - complex(psb_spk_), allocatable :: valsnd(:) - integer(psb_ipk_), allocatable :: iprc(:) - integer(psb_ipk_) :: debug_level, debug_unit, err_act - character(len=20) :: name, ch_err - - info=psb_success_ - name='psb_c_remote_vect' - 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() - - ctxt = desc_a%get_context() - icomm = desc_a%get_mpic() - - Call psb_info(ctxt, me, np) - - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),': Start' - - Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& - & bsdindx(np+1), stat=info) - - if (info /= psb_success_) then - info=psb_err_alloc_dealloc_ - call psb_errpush(info,name) - goto 9999 - end if - - call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) - - icomm = desc_a%get_mpic() - sdsz(:) = 0 - rvsz(:) = 0 - sdsi(:) = 0 - rvsi(:) = 0 - brvindx(:) = 0 - bsdindx(:) = 0 - counter = 1 - idxs = 0 - idxr = 0 - do i=1,n - if (iprc(i) >=0) then - sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 - else - write(0,*)me,name,' Error from fnd_owner: ',iprc(i) - end if - end do - call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& - & rvsz,1,psb_mpi_mpk_,icomm,minfo) - if (minfo /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='mpi_alltoall') - goto 9999 - end if - !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) - nsnds = count(sdsz /= 0) - nrcvs = count(rvsz /= 0) - 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 - - iszs = sum(sdsz) - iszr = sum(rvsz) - call psb_realloc(iszs,lsnd,info) - if (info == 0) call psb_realloc(iszs,valsnd,info) - if (info == 0) call psb_realloc(iszr,x,info) - if (info == 0) call psb_realloc(iszr,ix,info) - if (info /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='realloc') - goto 9999 - end if - do k=1, n - proc = iprc(k) - sdsi(proc+1) = sdsi(proc+1) + 1 - lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) - valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) - end do - do proc=0,np-1 - if (sdsi(proc+1) /= sdsz(proc+1)) & - & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) - end do - - select case(psb_get_sp_a2av_alg()) - case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) - call psb_simple_a2av(valsnd,sdsz,bsdindx,& - & x,rvsz,brvindx,ctxt,info) - if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& - & ix,rvsz,brvindx,ctxt,info) - case(psb_sp_a2av_mpi_) - - call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_c_spk_,& - & x,rvsz,brvindx,psb_mpi_c_spk_,icomm,minfo) - if (minfo == mpi_success) & - & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& - & ix,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 - - Deallocate(brvindx,bsdindx,rvsz,sdsz,& - & lsnd,valsnd,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(ctxt,err_act) - - return - -End Subroutine psb_c_remote_vect diff --git a/base/tools/psb_c_remote_vect.F90 b/base/tools/psb_c_remote_vect.F90 new file mode 100644 index 00000000..bd5286fa --- /dev/null +++ b/base/tools/psb_c_remote_vect.F90 @@ -0,0 +1,223 @@ +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006-2018 +! Salvatore Filippone +! Alfredo Buttari +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: psb_c_remote_vect.f90 +! +! Subroutine: +! This routine does the retrieval of remote vector entries. +! +! There are three possible exchange algorithms: +! 1. Use MPI_Alltoallv +! 2. Use psb_simple_a2av +! 3. Use psb_simple_triad_a2av +! Default choice is 3. The MPI variant has proved to be inefficient; +! that is because it is not persistent, therefore you pay the initialization price +! every time, and it is not optimized for a sparse communication pattern, +! most MPI implementations assume that all communications are non-empty. +! The PSB_SIMPLE variants reuse the same communicator, and go for a simplistic +! sequence of sends/receive that is quite efficient for a sparse communication +! pattern. To be refined/reviewed in the future to compare with neighbour +! persistent collectives. +! +! Arguments: +! desc_a - type(psb_desc_type). The communication descriptor. +! info - integer. Return code +! rowcnv - logical Should row/col indices be converted +! colcnv - logical to/from global numbering when sent/received? +! default is .TRUE. +! rowscale - logical Should row/col indices on output be remapped +! colscale - logical from MIN:MAX to 1:(MAX-MIN+1) ? +! default is .FALSE. +! (commmon use is ROWSCALE=.TRUE., COLSCALE=.FALSE.) +! data - integer Which index list in desc_a should be used to retrieve +! rows, default psb_comm_halo_ +! psb_comm_halo_ use halo_index +! psb_comm_ext_ use ext_index +! psb_comm_ovrl_ DISABLED for this routine. +! +subroutine psb_c_remote_vect(n,v,iv,desc_a,x,ix, info) + use psb_base_mod, psb_protect_name => psb_c_remote_vect + +#ifdef MPI_MOD + use mpi +#endif + Implicit None +#ifdef MPI_H + include 'mpif.h' +#endif + integer(psb_ipk_), intent(in) :: n + complex(psb_spk_), Intent(in) :: v(:) + integer(psb_lpk_), Intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + complex(psb_spk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + ! ...local scalars.... + type(psb_ctxt_type) :: ctxt + integer(psb_ipk_) :: np,me + integer(psb_ipk_) :: counter, proc, i, & + & j, idxs,idxr, k, iszs, iszr + integer(psb_ipk_) :: nrcvs, nsnds + integer(psb_mpk_) :: icomm, minfo + integer(psb_mpk_), allocatable :: brvindx(:), & + & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) + integer(psb_lpk_), allocatable :: lsnd(:) + complex(psb_spk_), allocatable :: valsnd(:) + integer(psb_ipk_), allocatable :: iprc(:) + integer(psb_ipk_) :: debug_level, debug_unit, err_act + character(len=20) :: name, ch_err + + info=psb_success_ + name='psb_c_remote_vect' + 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() + + ctxt = desc_a%get_context() + icomm = desc_a%get_mpic() + + Call psb_info(ctxt, me, np) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),': Start' + + Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& + & bsdindx(np+1), stat=info) + + if (info /= psb_success_) then + info=psb_err_alloc_dealloc_ + call psb_errpush(info,name) + goto 9999 + end if + + call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) + + icomm = desc_a%get_mpic() + sdsz(:) = 0 + rvsz(:) = 0 + sdsi(:) = 0 + rvsi(:) = 0 + brvindx(:) = 0 + bsdindx(:) = 0 + counter = 1 + idxs = 0 + idxr = 0 + do i=1,n + if (iprc(i) >=0) then + sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 + else + write(0,*)me,name,' Error from fnd_owner: ',iprc(i) + end if + end do + call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& + & rvsz,1,psb_mpi_mpk_,icomm,minfo) + if (minfo /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='mpi_alltoall') + goto 9999 + end if + !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) + nsnds = count(sdsz /= 0) + nrcvs = count(rvsz /= 0) + 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 + + iszs = sum(sdsz) + iszr = sum(rvsz) + call psb_realloc(iszs,lsnd,info) + if (info == 0) call psb_realloc(iszs,valsnd,info) + if (info == 0) call psb_realloc(iszr,x,info) + if (info == 0) call psb_realloc(iszr,ix,info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='realloc') + goto 9999 + end if + do k=1, n + proc = iprc(k) + sdsi(proc+1) = sdsi(proc+1) + 1 + lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) + valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) + end do + do proc=0,np-1 + if (sdsi(proc+1) /= sdsz(proc+1)) & + & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) + end do + + select case(psb_get_sp_a2av_alg()) + case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) + call psb_simple_a2av(valsnd,sdsz,bsdindx,& + & x,rvsz,brvindx,ctxt,info) + if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& + & ix,rvsz,brvindx,ctxt,info) + case(psb_sp_a2av_mpi_) + + call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_c_spk_,& + & x,rvsz,brvindx,psb_mpi_c_spk_,icomm,minfo) + if (minfo == mpi_success) & + & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& + & ix,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 + + Deallocate(brvindx,bsdindx,rvsz,sdsz,& + & lsnd,valsnd,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(ctxt,err_act) + + return + +End Subroutine psb_c_remote_vect diff --git a/base/tools/psb_casb.f90 b/base/tools/psb_casb.f90 index d43b062c..83e2715b 100644 --- a/base/tools/psb_casb.f90 +++ b/base/tools/psb_casb.f90 @@ -113,7 +113,7 @@ subroutine psb_casb_vect(x, desc_a, info, mold, scratch) integer(psb_ipk_) :: nrmv, nx, i nrmv = x%get_nrmv() - call psb_c_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) + call psb_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) nx = size(vx) call psb_realloc(nx,ivx,info) call desc_a%g2l(lvx,ivx,info,owned=.true.) diff --git a/base/tools/psb_d_remote_mat.F90 b/base/tools/psb_d_remote_mat.F90 index a5f12755..35116dc5 100644 --- a/base/tools/psb_d_remote_mat.F90 +++ b/base/tools/psb_d_remote_mat.F90 @@ -274,163 +274,3 @@ Subroutine psb_ld_remote_mat(a,desc_a,b,info) return End Subroutine psb_ld_remote_mat - - -subroutine psb_d_remote_vect(n,v,iv,desc_a,x,ix, info) - use psb_base_mod, psb_protect_name => psb_d_remote_vect - -#ifdef MPI_MOD - use mpi -#endif - Implicit None -#ifdef MPI_H - include 'mpif.h' -#endif - integer(psb_ipk_), intent(in) :: n - real(psb_dpk_), Intent(in) :: v(:) - integer(psb_lpk_), Intent(in) :: iv(:) - type(psb_desc_type),intent(in) :: desc_a - real(psb_dpk_), allocatable, intent(out) :: x(:) - integer(psb_lpk_), allocatable, intent(out) :: ix(:) - integer(psb_ipk_), intent(out) :: info - ! ...local scalars.... - type(psb_ctxt_type) :: ctxt - integer(psb_ipk_) :: np,me - integer(psb_ipk_) :: counter, proc, i, & - & j, idxs,idxr, k, iszs, iszr - integer(psb_ipk_) :: nrcvs, nsnds - integer(psb_mpk_) :: icomm, minfo - integer(psb_mpk_), allocatable :: brvindx(:), & - & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) - integer(psb_lpk_), allocatable :: lsnd(:) - real(psb_dpk_), allocatable :: valsnd(:) - integer(psb_ipk_), allocatable :: iprc(:) - integer(psb_ipk_) :: debug_level, debug_unit, err_act - character(len=20) :: name, ch_err - - info=psb_success_ - name='psb_d_remote_vect' - 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() - - ctxt = desc_a%get_context() - icomm = desc_a%get_mpic() - - Call psb_info(ctxt, me, np) - - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),': Start' - - Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& - & bsdindx(np+1), stat=info) - - if (info /= psb_success_) then - info=psb_err_alloc_dealloc_ - call psb_errpush(info,name) - goto 9999 - end if - - call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) - - icomm = desc_a%get_mpic() - sdsz(:) = 0 - rvsz(:) = 0 - sdsi(:) = 0 - rvsi(:) = 0 - brvindx(:) = 0 - bsdindx(:) = 0 - counter = 1 - idxs = 0 - idxr = 0 - do i=1,n - if (iprc(i) >=0) then - sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 - else - write(0,*)me,name,' Error from fnd_owner: ',iprc(i) - end if - end do - call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& - & rvsz,1,psb_mpi_mpk_,icomm,minfo) - if (minfo /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='mpi_alltoall') - goto 9999 - end if - !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) - nsnds = count(sdsz /= 0) - nrcvs = count(rvsz /= 0) - 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 - - iszs = sum(sdsz) - iszr = sum(rvsz) - call psb_realloc(iszs,lsnd,info) - if (info == 0) call psb_realloc(iszs,valsnd,info) - if (info == 0) call psb_realloc(iszr,x,info) - if (info == 0) call psb_realloc(iszr,ix,info) - if (info /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='realloc') - goto 9999 - end if - do k=1, n - proc = iprc(k) - sdsi(proc+1) = sdsi(proc+1) + 1 - lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) - valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) - end do - do proc=0,np-1 - if (sdsi(proc+1) /= sdsz(proc+1)) & - & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) - end do - - select case(psb_get_sp_a2av_alg()) - case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) - call psb_simple_a2av(valsnd,sdsz,bsdindx,& - & x,rvsz,brvindx,ctxt,info) - if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& - & ix,rvsz,brvindx,ctxt,info) - case(psb_sp_a2av_mpi_) - - call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_r_dpk_,& - & x,rvsz,brvindx,psb_mpi_r_dpk_,icomm,minfo) - if (minfo == mpi_success) & - & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& - & ix,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 - - Deallocate(brvindx,bsdindx,rvsz,sdsz,& - & lsnd,valsnd,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(ctxt,err_act) - - return - -End Subroutine psb_d_remote_vect diff --git a/base/tools/psb_d_remote_vect.F90 b/base/tools/psb_d_remote_vect.F90 new file mode 100644 index 00000000..4a409fa5 --- /dev/null +++ b/base/tools/psb_d_remote_vect.F90 @@ -0,0 +1,223 @@ +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006-2018 +! Salvatore Filippone +! Alfredo Buttari +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: psb_d_remote_vect.f90 +! +! Subroutine: +! This routine does the retrieval of remote vector entries. +! +! There are three possible exchange algorithms: +! 1. Use MPI_Alltoallv +! 2. Use psb_simple_a2av +! 3. Use psb_simple_triad_a2av +! Default choice is 3. The MPI variant has proved to be inefficient; +! that is because it is not persistent, therefore you pay the initialization price +! every time, and it is not optimized for a sparse communication pattern, +! most MPI implementations assume that all communications are non-empty. +! The PSB_SIMPLE variants reuse the same communicator, and go for a simplistic +! sequence of sends/receive that is quite efficient for a sparse communication +! pattern. To be refined/reviewed in the future to compare with neighbour +! persistent collectives. +! +! Arguments: +! desc_a - type(psb_desc_type). The communication descriptor. +! info - integer. Return code +! rowcnv - logical Should row/col indices be converted +! colcnv - logical to/from global numbering when sent/received? +! default is .TRUE. +! rowscale - logical Should row/col indices on output be remapped +! colscale - logical from MIN:MAX to 1:(MAX-MIN+1) ? +! default is .FALSE. +! (commmon use is ROWSCALE=.TRUE., COLSCALE=.FALSE.) +! data - integer Which index list in desc_a should be used to retrieve +! rows, default psb_comm_halo_ +! psb_comm_halo_ use halo_index +! psb_comm_ext_ use ext_index +! psb_comm_ovrl_ DISABLED for this routine. +! +subroutine psb_d_remote_vect(n,v,iv,desc_a,x,ix, info) + use psb_base_mod, psb_protect_name => psb_d_remote_vect + +#ifdef MPI_MOD + use mpi +#endif + Implicit None +#ifdef MPI_H + include 'mpif.h' +#endif + integer(psb_ipk_), intent(in) :: n + real(psb_dpk_), Intent(in) :: v(:) + integer(psb_lpk_), Intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + real(psb_dpk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + ! ...local scalars.... + type(psb_ctxt_type) :: ctxt + integer(psb_ipk_) :: np,me + integer(psb_ipk_) :: counter, proc, i, & + & j, idxs,idxr, k, iszs, iszr + integer(psb_ipk_) :: nrcvs, nsnds + integer(psb_mpk_) :: icomm, minfo + integer(psb_mpk_), allocatable :: brvindx(:), & + & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) + integer(psb_lpk_), allocatable :: lsnd(:) + real(psb_dpk_), allocatable :: valsnd(:) + integer(psb_ipk_), allocatable :: iprc(:) + integer(psb_ipk_) :: debug_level, debug_unit, err_act + character(len=20) :: name, ch_err + + info=psb_success_ + name='psb_d_remote_vect' + 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() + + ctxt = desc_a%get_context() + icomm = desc_a%get_mpic() + + Call psb_info(ctxt, me, np) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),': Start' + + Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& + & bsdindx(np+1), stat=info) + + if (info /= psb_success_) then + info=psb_err_alloc_dealloc_ + call psb_errpush(info,name) + goto 9999 + end if + + call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) + + icomm = desc_a%get_mpic() + sdsz(:) = 0 + rvsz(:) = 0 + sdsi(:) = 0 + rvsi(:) = 0 + brvindx(:) = 0 + bsdindx(:) = 0 + counter = 1 + idxs = 0 + idxr = 0 + do i=1,n + if (iprc(i) >=0) then + sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 + else + write(0,*)me,name,' Error from fnd_owner: ',iprc(i) + end if + end do + call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& + & rvsz,1,psb_mpi_mpk_,icomm,minfo) + if (minfo /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='mpi_alltoall') + goto 9999 + end if + !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) + nsnds = count(sdsz /= 0) + nrcvs = count(rvsz /= 0) + 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 + + iszs = sum(sdsz) + iszr = sum(rvsz) + call psb_realloc(iszs,lsnd,info) + if (info == 0) call psb_realloc(iszs,valsnd,info) + if (info == 0) call psb_realloc(iszr,x,info) + if (info == 0) call psb_realloc(iszr,ix,info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='realloc') + goto 9999 + end if + do k=1, n + proc = iprc(k) + sdsi(proc+1) = sdsi(proc+1) + 1 + lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) + valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) + end do + do proc=0,np-1 + if (sdsi(proc+1) /= sdsz(proc+1)) & + & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) + end do + + select case(psb_get_sp_a2av_alg()) + case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) + call psb_simple_a2av(valsnd,sdsz,bsdindx,& + & x,rvsz,brvindx,ctxt,info) + if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& + & ix,rvsz,brvindx,ctxt,info) + case(psb_sp_a2av_mpi_) + + call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_r_dpk_,& + & x,rvsz,brvindx,psb_mpi_r_dpk_,icomm,minfo) + if (minfo == mpi_success) & + & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& + & ix,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 + + Deallocate(brvindx,bsdindx,rvsz,sdsz,& + & lsnd,valsnd,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(ctxt,err_act) + + return + +End Subroutine psb_d_remote_vect diff --git a/base/tools/psb_dasb.f90 b/base/tools/psb_dasb.f90 index f8225496..19a19ff1 100644 --- a/base/tools/psb_dasb.f90 +++ b/base/tools/psb_dasb.f90 @@ -113,7 +113,7 @@ subroutine psb_dasb_vect(x, desc_a, info, mold, scratch) integer(psb_ipk_) :: nrmv, nx, i nrmv = x%get_nrmv() - call psb_d_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) + call psb_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) nx = size(vx) call psb_realloc(nx,ivx,info) call desc_a%g2l(lvx,ivx,info,owned=.true.) diff --git a/base/tools/psb_e_remote_vect.F90 b/base/tools/psb_e_remote_vect.F90 new file mode 100644 index 00000000..9fb15ff9 --- /dev/null +++ b/base/tools/psb_e_remote_vect.F90 @@ -0,0 +1,223 @@ +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006-2018 +! Salvatore Filippone +! Alfredo Buttari +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: psb_e_remote_vect.f90 +! +! Subroutine: +! This routine does the retrieval of remote vector entries. +! +! There are three possible exchange algorithms: +! 1. Use MPI_Alltoallv +! 2. Use psb_simple_a2av +! 3. Use psb_simple_triad_a2av +! Default choice is 3. The MPI variant has proved to be inefficient; +! that is because it is not persistent, therefore you pay the initialization price +! every time, and it is not optimized for a sparse communication pattern, +! most MPI implementations assume that all communications are non-empty. +! The PSB_SIMPLE variants reuse the same communicator, and go for a simplistic +! sequence of sends/receive that is quite efficient for a sparse communication +! pattern. To be refined/reviewed in the future to compare with neighbour +! persistent collectives. +! +! Arguments: +! desc_a - type(psb_desc_type). The communication descriptor. +! info - integer. Return code +! rowcnv - logical Should row/col indices be converted +! colcnv - logical to/from global numbering when sent/received? +! default is .TRUE. +! rowscale - logical Should row/col indices on output be remapped +! colscale - logical from MIN:MAX to 1:(MAX-MIN+1) ? +! default is .FALSE. +! (commmon use is ROWSCALE=.TRUE., COLSCALE=.FALSE.) +! data - integer Which index list in desc_a should be used to retrieve +! rows, default psb_comm_halo_ +! psb_comm_halo_ use halo_index +! psb_comm_ext_ use ext_index +! psb_comm_ovrl_ DISABLED for this routine. +! +subroutine psb_e_remote_vect(n,v,iv,desc_a,x,ix, info) + use psb_base_mod, psb_protect_name => psb_e_remote_vect + +#ifdef MPI_MOD + use mpi +#endif + Implicit None +#ifdef MPI_H + include 'mpif.h' +#endif + integer(psb_ipk_), intent(in) :: n + integer(psb_epk_), Intent(in) :: v(:) + integer(psb_lpk_), Intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + integer(psb_epk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + ! ...local scalars.... + type(psb_ctxt_type) :: ctxt + integer(psb_ipk_) :: np,me + integer(psb_ipk_) :: counter, proc, i, & + & j, idxs,idxr, k, iszs, iszr + integer(psb_ipk_) :: nrcvs, nsnds + integer(psb_mpk_) :: icomm, minfo + integer(psb_mpk_), allocatable :: brvindx(:), & + & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) + integer(psb_lpk_), allocatable :: lsnd(:) + integer(psb_epk_), allocatable :: valsnd(:) + integer(psb_ipk_), allocatable :: iprc(:) + integer(psb_ipk_) :: debug_level, debug_unit, err_act + character(len=20) :: name, ch_err + + info=psb_success_ + name='psb_e_remote_vect' + 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() + + ctxt = desc_a%get_context() + icomm = desc_a%get_mpic() + + Call psb_info(ctxt, me, np) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),': Start' + + Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& + & bsdindx(np+1), stat=info) + + if (info /= psb_success_) then + info=psb_err_alloc_dealloc_ + call psb_errpush(info,name) + goto 9999 + end if + + call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) + + icomm = desc_a%get_mpic() + sdsz(:) = 0 + rvsz(:) = 0 + sdsi(:) = 0 + rvsi(:) = 0 + brvindx(:) = 0 + bsdindx(:) = 0 + counter = 1 + idxs = 0 + idxr = 0 + do i=1,n + if (iprc(i) >=0) then + sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 + else + write(0,*)me,name,' Error from fnd_owner: ',iprc(i) + end if + end do + call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& + & rvsz,1,psb_mpi_mpk_,icomm,minfo) + if (minfo /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='mpi_alltoall') + goto 9999 + end if + !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) + nsnds = count(sdsz /= 0) + nrcvs = count(rvsz /= 0) + 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 + + iszs = sum(sdsz) + iszr = sum(rvsz) + call psb_realloc(iszs,lsnd,info) + if (info == 0) call psb_realloc(iszs,valsnd,info) + if (info == 0) call psb_realloc(iszr,x,info) + if (info == 0) call psb_realloc(iszr,ix,info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='realloc') + goto 9999 + end if + do k=1, n + proc = iprc(k) + sdsi(proc+1) = sdsi(proc+1) + 1 + lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) + valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) + end do + do proc=0,np-1 + if (sdsi(proc+1) /= sdsz(proc+1)) & + & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) + end do + + select case(psb_get_sp_a2av_alg()) + case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) + call psb_simple_a2av(valsnd,sdsz,bsdindx,& + & x,rvsz,brvindx,ctxt,info) + if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& + & ix,rvsz,brvindx,ctxt,info) + case(psb_sp_a2av_mpi_) + + call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_epk_,& + & x,rvsz,brvindx,psb_mpi_epk_,icomm,minfo) + if (minfo == mpi_success) & + & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& + & ix,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 + + Deallocate(brvindx,bsdindx,rvsz,sdsz,& + & lsnd,valsnd,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(ctxt,err_act) + + return + +End Subroutine psb_e_remote_vect diff --git a/base/tools/psb_i2_remote_vect.F90 b/base/tools/psb_i2_remote_vect.F90 new file mode 100644 index 00000000..3f6bffbd --- /dev/null +++ b/base/tools/psb_i2_remote_vect.F90 @@ -0,0 +1,223 @@ +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006-2018 +! Salvatore Filippone +! Alfredo Buttari +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: psb_i2_remote_vect.f90 +! +! Subroutine: +! This routine does the retrieval of remote vector entries. +! +! There are three possible exchange algorithms: +! 1. Use MPI_Alltoallv +! 2. Use psb_simple_a2av +! 3. Use psb_simple_triad_a2av +! Default choice is 3. The MPI variant has proved to be inefficient; +! that is because it is not persistent, therefore you pay the initialization price +! every time, and it is not optimized for a sparse communication pattern, +! most MPI implementations assume that all communications are non-empty. +! The PSB_SIMPLE variants reuse the same communicator, and go for a simplistic +! sequence of sends/receive that is quite efficient for a sparse communication +! pattern. To be refined/reviewed in the future to compare with neighbour +! persistent collectives. +! +! Arguments: +! desc_a - type(psb_desc_type). The communication descriptor. +! info - integer. Return code +! rowcnv - logical Should row/col indices be converted +! colcnv - logical to/from global numbering when sent/received? +! default is .TRUE. +! rowscale - logical Should row/col indices on output be remapped +! colscale - logical from MIN:MAX to 1:(MAX-MIN+1) ? +! default is .FALSE. +! (commmon use is ROWSCALE=.TRUE., COLSCALE=.FALSE.) +! data - integer Which index list in desc_a should be used to retrieve +! rows, default psb_comm_halo_ +! psb_comm_halo_ use halo_index +! psb_comm_ext_ use ext_index +! psb_comm_ovrl_ DISABLED for this routine. +! +subroutine psb_i2_remote_vect(n,v,iv,desc_a,x,ix, info) + use psb_base_mod, psb_protect_name => psb_i2_remote_vect + +#ifdef MPI_MOD + use mpi +#endif + Implicit None +#ifdef MPI_H + include 'mpif.h' +#endif + integer(psb_ipk_), intent(in) :: n + integer(psb_i2pk_), Intent(in) :: v(:) + integer(psb_lpk_), Intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + integer(psb_i2pk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + ! ...local scalars.... + type(psb_ctxt_type) :: ctxt + integer(psb_ipk_) :: np,me + integer(psb_ipk_) :: counter, proc, i, & + & j, idxs,idxr, k, iszs, iszr + integer(psb_ipk_) :: nrcvs, nsnds + integer(psb_mpk_) :: icomm, minfo + integer(psb_mpk_), allocatable :: brvindx(:), & + & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) + integer(psb_lpk_), allocatable :: lsnd(:) + integer(psb_i2pk_), allocatable :: valsnd(:) + integer(psb_ipk_), allocatable :: iprc(:) + integer(psb_ipk_) :: debug_level, debug_unit, err_act + character(len=20) :: name, ch_err + + info=psb_success_ + name='psb_i2_remote_vect' + 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() + + ctxt = desc_a%get_context() + icomm = desc_a%get_mpic() + + Call psb_info(ctxt, me, np) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),': Start' + + Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& + & bsdindx(np+1), stat=info) + + if (info /= psb_success_) then + info=psb_err_alloc_dealloc_ + call psb_errpush(info,name) + goto 9999 + end if + + call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) + + icomm = desc_a%get_mpic() + sdsz(:) = 0 + rvsz(:) = 0 + sdsi(:) = 0 + rvsi(:) = 0 + brvindx(:) = 0 + bsdindx(:) = 0 + counter = 1 + idxs = 0 + idxr = 0 + do i=1,n + if (iprc(i) >=0) then + sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 + else + write(0,*)me,name,' Error from fnd_owner: ',iprc(i) + end if + end do + call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& + & rvsz,1,psb_mpi_mpk_,icomm,minfo) + if (minfo /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='mpi_alltoall') + goto 9999 + end if + !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) + nsnds = count(sdsz /= 0) + nrcvs = count(rvsz /= 0) + 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 + + iszs = sum(sdsz) + iszr = sum(rvsz) + call psb_realloc(iszs,lsnd,info) + if (info == 0) call psb_realloc(iszs,valsnd,info) + if (info == 0) call psb_realloc(iszr,x,info) + if (info == 0) call psb_realloc(iszr,ix,info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='realloc') + goto 9999 + end if + do k=1, n + proc = iprc(k) + sdsi(proc+1) = sdsi(proc+1) + 1 + lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) + valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) + end do + do proc=0,np-1 + if (sdsi(proc+1) /= sdsz(proc+1)) & + & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) + end do + + select case(psb_get_sp_a2av_alg()) + case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) + call psb_simple_a2av(valsnd,sdsz,bsdindx,& + & x,rvsz,brvindx,ctxt,info) + if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& + & ix,rvsz,brvindx,ctxt,info) + case(psb_sp_a2av_mpi_) + + call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_i2pk_,& + & x,rvsz,brvindx,psb_mpi_i2pk_,icomm,minfo) + if (minfo == mpi_success) & + & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& + & ix,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 + + Deallocate(brvindx,bsdindx,rvsz,sdsz,& + & lsnd,valsnd,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(ctxt,err_act) + + return + +End Subroutine psb_i2_remote_vect diff --git a/base/tools/psb_iasb.f90 b/base/tools/psb_iasb.f90 index 5e09d331..f5e5669f 100644 --- a/base/tools/psb_iasb.f90 +++ b/base/tools/psb_iasb.f90 @@ -113,7 +113,7 @@ subroutine psb_iasb_vect(x, desc_a, info, mold, scratch) integer(psb_ipk_) :: nrmv, nx, i nrmv = x%get_nrmv() - call psb_i_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) + call psb_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) nx = size(vx) call psb_realloc(nx,ivx,info) call desc_a%g2l(lvx,ivx,info,owned=.true.) diff --git a/base/tools/psb_lasb.f90 b/base/tools/psb_lasb.f90 index 1110ee11..baf55320 100644 --- a/base/tools/psb_lasb.f90 +++ b/base/tools/psb_lasb.f90 @@ -113,7 +113,7 @@ subroutine psb_lasb_vect(x, desc_a, info, mold, scratch) integer(psb_ipk_) :: nrmv, nx, i nrmv = x%get_nrmv() - call psb_l_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) + call psb_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) nx = size(vx) call psb_realloc(nx,ivx,info) call desc_a%g2l(lvx,ivx,info,owned=.true.) diff --git a/base/tools/psb_m_remote_vect.F90 b/base/tools/psb_m_remote_vect.F90 new file mode 100644 index 00000000..01b5aeb3 --- /dev/null +++ b/base/tools/psb_m_remote_vect.F90 @@ -0,0 +1,223 @@ +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006-2018 +! Salvatore Filippone +! Alfredo Buttari +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: psb_m_remote_vect.f90 +! +! Subroutine: +! This routine does the retrieval of remote vector entries. +! +! There are three possible exchange algorithms: +! 1. Use MPI_Alltoallv +! 2. Use psb_simple_a2av +! 3. Use psb_simple_triad_a2av +! Default choice is 3. The MPI variant has proved to be inefficient; +! that is because it is not persistent, therefore you pay the initialization price +! every time, and it is not optimized for a sparse communication pattern, +! most MPI implementations assume that all communications are non-empty. +! The PSB_SIMPLE variants reuse the same communicator, and go for a simplistic +! sequence of sends/receive that is quite efficient for a sparse communication +! pattern. To be refined/reviewed in the future to compare with neighbour +! persistent collectives. +! +! Arguments: +! desc_a - type(psb_desc_type). The communication descriptor. +! info - integer. Return code +! rowcnv - logical Should row/col indices be converted +! colcnv - logical to/from global numbering when sent/received? +! default is .TRUE. +! rowscale - logical Should row/col indices on output be remapped +! colscale - logical from MIN:MAX to 1:(MAX-MIN+1) ? +! default is .FALSE. +! (commmon use is ROWSCALE=.TRUE., COLSCALE=.FALSE.) +! data - integer Which index list in desc_a should be used to retrieve +! rows, default psb_comm_halo_ +! psb_comm_halo_ use halo_index +! psb_comm_ext_ use ext_index +! psb_comm_ovrl_ DISABLED for this routine. +! +subroutine psb_m_remote_vect(n,v,iv,desc_a,x,ix, info) + use psb_base_mod, psb_protect_name => psb_m_remote_vect + +#ifdef MPI_MOD + use mpi +#endif + Implicit None +#ifdef MPI_H + include 'mpif.h' +#endif + integer(psb_ipk_), intent(in) :: n + integer(psb_mpk_), Intent(in) :: v(:) + integer(psb_lpk_), Intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + integer(psb_mpk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + ! ...local scalars.... + type(psb_ctxt_type) :: ctxt + integer(psb_ipk_) :: np,me + integer(psb_ipk_) :: counter, proc, i, & + & j, idxs,idxr, k, iszs, iszr + integer(psb_ipk_) :: nrcvs, nsnds + integer(psb_mpk_) :: icomm, minfo + integer(psb_mpk_), allocatable :: brvindx(:), & + & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) + integer(psb_lpk_), allocatable :: lsnd(:) + integer(psb_mpk_), allocatable :: valsnd(:) + integer(psb_ipk_), allocatable :: iprc(:) + integer(psb_ipk_) :: debug_level, debug_unit, err_act + character(len=20) :: name, ch_err + + info=psb_success_ + name='psb_m_remote_vect' + 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() + + ctxt = desc_a%get_context() + icomm = desc_a%get_mpic() + + Call psb_info(ctxt, me, np) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),': Start' + + Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& + & bsdindx(np+1), stat=info) + + if (info /= psb_success_) then + info=psb_err_alloc_dealloc_ + call psb_errpush(info,name) + goto 9999 + end if + + call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) + + icomm = desc_a%get_mpic() + sdsz(:) = 0 + rvsz(:) = 0 + sdsi(:) = 0 + rvsi(:) = 0 + brvindx(:) = 0 + bsdindx(:) = 0 + counter = 1 + idxs = 0 + idxr = 0 + do i=1,n + if (iprc(i) >=0) then + sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 + else + write(0,*)me,name,' Error from fnd_owner: ',iprc(i) + end if + end do + call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& + & rvsz,1,psb_mpi_mpk_,icomm,minfo) + if (minfo /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='mpi_alltoall') + goto 9999 + end if + !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) + nsnds = count(sdsz /= 0) + nrcvs = count(rvsz /= 0) + 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 + + iszs = sum(sdsz) + iszr = sum(rvsz) + call psb_realloc(iszs,lsnd,info) + if (info == 0) call psb_realloc(iszs,valsnd,info) + if (info == 0) call psb_realloc(iszr,x,info) + if (info == 0) call psb_realloc(iszr,ix,info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='realloc') + goto 9999 + end if + do k=1, n + proc = iprc(k) + sdsi(proc+1) = sdsi(proc+1) + 1 + lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) + valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) + end do + do proc=0,np-1 + if (sdsi(proc+1) /= sdsz(proc+1)) & + & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) + end do + + select case(psb_get_sp_a2av_alg()) + case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) + call psb_simple_a2av(valsnd,sdsz,bsdindx,& + & x,rvsz,brvindx,ctxt,info) + if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& + & ix,rvsz,brvindx,ctxt,info) + case(psb_sp_a2av_mpi_) + + call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_mpk_,& + & x,rvsz,brvindx,psb_mpi_mpk_,icomm,minfo) + if (minfo == mpi_success) & + & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& + & ix,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 + + Deallocate(brvindx,bsdindx,rvsz,sdsz,& + & lsnd,valsnd,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(ctxt,err_act) + + return + +End Subroutine psb_m_remote_vect diff --git a/base/tools/psb_s_remote_mat.F90 b/base/tools/psb_s_remote_mat.F90 index d0886304..df64266b 100644 --- a/base/tools/psb_s_remote_mat.F90 +++ b/base/tools/psb_s_remote_mat.F90 @@ -274,163 +274,3 @@ Subroutine psb_ls_remote_mat(a,desc_a,b,info) return End Subroutine psb_ls_remote_mat - - -subroutine psb_s_remote_vect(n,v,iv,desc_a,x,ix, info) - use psb_base_mod, psb_protect_name => psb_s_remote_vect - -#ifdef MPI_MOD - use mpi -#endif - Implicit None -#ifdef MPI_H - include 'mpif.h' -#endif - integer(psb_ipk_), intent(in) :: n - real(psb_spk_), Intent(in) :: v(:) - integer(psb_lpk_), Intent(in) :: iv(:) - type(psb_desc_type),intent(in) :: desc_a - real(psb_spk_), allocatable, intent(out) :: x(:) - integer(psb_lpk_), allocatable, intent(out) :: ix(:) - integer(psb_ipk_), intent(out) :: info - ! ...local scalars.... - type(psb_ctxt_type) :: ctxt - integer(psb_ipk_) :: np,me - integer(psb_ipk_) :: counter, proc, i, & - & j, idxs,idxr, k, iszs, iszr - integer(psb_ipk_) :: nrcvs, nsnds - integer(psb_mpk_) :: icomm, minfo - integer(psb_mpk_), allocatable :: brvindx(:), & - & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) - integer(psb_lpk_), allocatable :: lsnd(:) - real(psb_spk_), allocatable :: valsnd(:) - integer(psb_ipk_), allocatable :: iprc(:) - integer(psb_ipk_) :: debug_level, debug_unit, err_act - character(len=20) :: name, ch_err - - info=psb_success_ - name='psb_s_remote_vect' - 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() - - ctxt = desc_a%get_context() - icomm = desc_a%get_mpic() - - Call psb_info(ctxt, me, np) - - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),': Start' - - Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& - & bsdindx(np+1), stat=info) - - if (info /= psb_success_) then - info=psb_err_alloc_dealloc_ - call psb_errpush(info,name) - goto 9999 - end if - - call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) - - icomm = desc_a%get_mpic() - sdsz(:) = 0 - rvsz(:) = 0 - sdsi(:) = 0 - rvsi(:) = 0 - brvindx(:) = 0 - bsdindx(:) = 0 - counter = 1 - idxs = 0 - idxr = 0 - do i=1,n - if (iprc(i) >=0) then - sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 - else - write(0,*)me,name,' Error from fnd_owner: ',iprc(i) - end if - end do - call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& - & rvsz,1,psb_mpi_mpk_,icomm,minfo) - if (minfo /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='mpi_alltoall') - goto 9999 - end if - !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) - nsnds = count(sdsz /= 0) - nrcvs = count(rvsz /= 0) - 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 - - iszs = sum(sdsz) - iszr = sum(rvsz) - call psb_realloc(iszs,lsnd,info) - if (info == 0) call psb_realloc(iszs,valsnd,info) - if (info == 0) call psb_realloc(iszr,x,info) - if (info == 0) call psb_realloc(iszr,ix,info) - if (info /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='realloc') - goto 9999 - end if - do k=1, n - proc = iprc(k) - sdsi(proc+1) = sdsi(proc+1) + 1 - lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) - valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) - end do - do proc=0,np-1 - if (sdsi(proc+1) /= sdsz(proc+1)) & - & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) - end do - - select case(psb_get_sp_a2av_alg()) - case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) - call psb_simple_a2av(valsnd,sdsz,bsdindx,& - & x,rvsz,brvindx,ctxt,info) - if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& - & ix,rvsz,brvindx,ctxt,info) - case(psb_sp_a2av_mpi_) - - call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_r_spk_,& - & x,rvsz,brvindx,psb_mpi_r_spk_,icomm,minfo) - if (minfo == mpi_success) & - & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& - & ix,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 - - Deallocate(brvindx,bsdindx,rvsz,sdsz,& - & lsnd,valsnd,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(ctxt,err_act) - - return - -End Subroutine psb_s_remote_vect diff --git a/base/tools/psb_s_remote_vect.F90 b/base/tools/psb_s_remote_vect.F90 new file mode 100644 index 00000000..a8464663 --- /dev/null +++ b/base/tools/psb_s_remote_vect.F90 @@ -0,0 +1,223 @@ +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006-2018 +! Salvatore Filippone +! Alfredo Buttari +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: psb_s_remote_vect.f90 +! +! Subroutine: +! This routine does the retrieval of remote vector entries. +! +! There are three possible exchange algorithms: +! 1. Use MPI_Alltoallv +! 2. Use psb_simple_a2av +! 3. Use psb_simple_triad_a2av +! Default choice is 3. The MPI variant has proved to be inefficient; +! that is because it is not persistent, therefore you pay the initialization price +! every time, and it is not optimized for a sparse communication pattern, +! most MPI implementations assume that all communications are non-empty. +! The PSB_SIMPLE variants reuse the same communicator, and go for a simplistic +! sequence of sends/receive that is quite efficient for a sparse communication +! pattern. To be refined/reviewed in the future to compare with neighbour +! persistent collectives. +! +! Arguments: +! desc_a - type(psb_desc_type). The communication descriptor. +! info - integer. Return code +! rowcnv - logical Should row/col indices be converted +! colcnv - logical to/from global numbering when sent/received? +! default is .TRUE. +! rowscale - logical Should row/col indices on output be remapped +! colscale - logical from MIN:MAX to 1:(MAX-MIN+1) ? +! default is .FALSE. +! (commmon use is ROWSCALE=.TRUE., COLSCALE=.FALSE.) +! data - integer Which index list in desc_a should be used to retrieve +! rows, default psb_comm_halo_ +! psb_comm_halo_ use halo_index +! psb_comm_ext_ use ext_index +! psb_comm_ovrl_ DISABLED for this routine. +! +subroutine psb_s_remote_vect(n,v,iv,desc_a,x,ix, info) + use psb_base_mod, psb_protect_name => psb_s_remote_vect + +#ifdef MPI_MOD + use mpi +#endif + Implicit None +#ifdef MPI_H + include 'mpif.h' +#endif + integer(psb_ipk_), intent(in) :: n + real(psb_spk_), Intent(in) :: v(:) + integer(psb_lpk_), Intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + real(psb_spk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + ! ...local scalars.... + type(psb_ctxt_type) :: ctxt + integer(psb_ipk_) :: np,me + integer(psb_ipk_) :: counter, proc, i, & + & j, idxs,idxr, k, iszs, iszr + integer(psb_ipk_) :: nrcvs, nsnds + integer(psb_mpk_) :: icomm, minfo + integer(psb_mpk_), allocatable :: brvindx(:), & + & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) + integer(psb_lpk_), allocatable :: lsnd(:) + real(psb_spk_), allocatable :: valsnd(:) + integer(psb_ipk_), allocatable :: iprc(:) + integer(psb_ipk_) :: debug_level, debug_unit, err_act + character(len=20) :: name, ch_err + + info=psb_success_ + name='psb_s_remote_vect' + 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() + + ctxt = desc_a%get_context() + icomm = desc_a%get_mpic() + + Call psb_info(ctxt, me, np) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),': Start' + + Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& + & bsdindx(np+1), stat=info) + + if (info /= psb_success_) then + info=psb_err_alloc_dealloc_ + call psb_errpush(info,name) + goto 9999 + end if + + call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) + + icomm = desc_a%get_mpic() + sdsz(:) = 0 + rvsz(:) = 0 + sdsi(:) = 0 + rvsi(:) = 0 + brvindx(:) = 0 + bsdindx(:) = 0 + counter = 1 + idxs = 0 + idxr = 0 + do i=1,n + if (iprc(i) >=0) then + sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 + else + write(0,*)me,name,' Error from fnd_owner: ',iprc(i) + end if + end do + call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& + & rvsz,1,psb_mpi_mpk_,icomm,minfo) + if (minfo /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='mpi_alltoall') + goto 9999 + end if + !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) + nsnds = count(sdsz /= 0) + nrcvs = count(rvsz /= 0) + 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 + + iszs = sum(sdsz) + iszr = sum(rvsz) + call psb_realloc(iszs,lsnd,info) + if (info == 0) call psb_realloc(iszs,valsnd,info) + if (info == 0) call psb_realloc(iszr,x,info) + if (info == 0) call psb_realloc(iszr,ix,info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='realloc') + goto 9999 + end if + do k=1, n + proc = iprc(k) + sdsi(proc+1) = sdsi(proc+1) + 1 + lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) + valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) + end do + do proc=0,np-1 + if (sdsi(proc+1) /= sdsz(proc+1)) & + & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) + end do + + select case(psb_get_sp_a2av_alg()) + case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) + call psb_simple_a2av(valsnd,sdsz,bsdindx,& + & x,rvsz,brvindx,ctxt,info) + if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& + & ix,rvsz,brvindx,ctxt,info) + case(psb_sp_a2av_mpi_) + + call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_r_spk_,& + & x,rvsz,brvindx,psb_mpi_r_spk_,icomm,minfo) + if (minfo == mpi_success) & + & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& + & ix,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 + + Deallocate(brvindx,bsdindx,rvsz,sdsz,& + & lsnd,valsnd,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(ctxt,err_act) + + return + +End Subroutine psb_s_remote_vect diff --git a/base/tools/psb_sasb.f90 b/base/tools/psb_sasb.f90 index 7238fecf..315e24ff 100644 --- a/base/tools/psb_sasb.f90 +++ b/base/tools/psb_sasb.f90 @@ -113,7 +113,7 @@ subroutine psb_sasb_vect(x, desc_a, info, mold, scratch) integer(psb_ipk_) :: nrmv, nx, i nrmv = x%get_nrmv() - call psb_s_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) + call psb_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) nx = size(vx) call psb_realloc(nx,ivx,info) call desc_a%g2l(lvx,ivx,info,owned=.true.) diff --git a/base/tools/psb_z_remote_mat.F90 b/base/tools/psb_z_remote_mat.F90 index fc23bd4a..5461c5d5 100644 --- a/base/tools/psb_z_remote_mat.F90 +++ b/base/tools/psb_z_remote_mat.F90 @@ -274,163 +274,3 @@ Subroutine psb_lz_remote_mat(a,desc_a,b,info) return End Subroutine psb_lz_remote_mat - - -subroutine psb_z_remote_vect(n,v,iv,desc_a,x,ix, info) - use psb_base_mod, psb_protect_name => psb_z_remote_vect - -#ifdef MPI_MOD - use mpi -#endif - Implicit None -#ifdef MPI_H - include 'mpif.h' -#endif - integer(psb_ipk_), intent(in) :: n - complex(psb_dpk_), Intent(in) :: v(:) - integer(psb_lpk_), Intent(in) :: iv(:) - type(psb_desc_type),intent(in) :: desc_a - complex(psb_dpk_), allocatable, intent(out) :: x(:) - integer(psb_lpk_), allocatable, intent(out) :: ix(:) - integer(psb_ipk_), intent(out) :: info - ! ...local scalars.... - type(psb_ctxt_type) :: ctxt - integer(psb_ipk_) :: np,me - integer(psb_ipk_) :: counter, proc, i, & - & j, idxs,idxr, k, iszs, iszr - integer(psb_ipk_) :: nrcvs, nsnds - integer(psb_mpk_) :: icomm, minfo - integer(psb_mpk_), allocatable :: brvindx(:), & - & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) - integer(psb_lpk_), allocatable :: lsnd(:) - complex(psb_dpk_), allocatable :: valsnd(:) - integer(psb_ipk_), allocatable :: iprc(:) - integer(psb_ipk_) :: debug_level, debug_unit, err_act - character(len=20) :: name, ch_err - - info=psb_success_ - name='psb_z_remote_vect' - 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() - - ctxt = desc_a%get_context() - icomm = desc_a%get_mpic() - - Call psb_info(ctxt, me, np) - - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),': Start' - - Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& - & bsdindx(np+1), stat=info) - - if (info /= psb_success_) then - info=psb_err_alloc_dealloc_ - call psb_errpush(info,name) - goto 9999 - end if - - call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) - - icomm = desc_a%get_mpic() - sdsz(:) = 0 - rvsz(:) = 0 - sdsi(:) = 0 - rvsi(:) = 0 - brvindx(:) = 0 - bsdindx(:) = 0 - counter = 1 - idxs = 0 - idxr = 0 - do i=1,n - if (iprc(i) >=0) then - sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 - else - write(0,*)me,name,' Error from fnd_owner: ',iprc(i) - end if - end do - call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& - & rvsz,1,psb_mpi_mpk_,icomm,minfo) - if (minfo /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='mpi_alltoall') - goto 9999 - end if - !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) - nsnds = count(sdsz /= 0) - nrcvs = count(rvsz /= 0) - 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 - - iszs = sum(sdsz) - iszr = sum(rvsz) - call psb_realloc(iszs,lsnd,info) - if (info == 0) call psb_realloc(iszs,valsnd,info) - if (info == 0) call psb_realloc(iszr,x,info) - if (info == 0) call psb_realloc(iszr,ix,info) - if (info /= psb_success_) then - info=psb_err_from_subroutine_ - call psb_errpush(info,name,a_err='realloc') - goto 9999 - end if - do k=1, n - proc = iprc(k) - sdsi(proc+1) = sdsi(proc+1) + 1 - lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) - valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) - end do - do proc=0,np-1 - if (sdsi(proc+1) /= sdsz(proc+1)) & - & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) - end do - - select case(psb_get_sp_a2av_alg()) - case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) - call psb_simple_a2av(valsnd,sdsz,bsdindx,& - & x,rvsz,brvindx,ctxt,info) - if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& - & ix,rvsz,brvindx,ctxt,info) - case(psb_sp_a2av_mpi_) - - call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_c_dpk_,& - & x,rvsz,brvindx,psb_mpi_c_dpk_,icomm,minfo) - if (minfo == mpi_success) & - & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& - & ix,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 - - Deallocate(brvindx,bsdindx,rvsz,sdsz,& - & lsnd,valsnd,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(ctxt,err_act) - - return - -End Subroutine psb_z_remote_vect diff --git a/base/tools/psb_z_remote_vect.F90 b/base/tools/psb_z_remote_vect.F90 new file mode 100644 index 00000000..ed705bb5 --- /dev/null +++ b/base/tools/psb_z_remote_vect.F90 @@ -0,0 +1,223 @@ +! +! Parallel Sparse BLAS version 3.5 +! (C) Copyright 2006-2018 +! Salvatore Filippone +! Alfredo Buttari +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the PSBLAS group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: psb_z_remote_vect.f90 +! +! Subroutine: +! This routine does the retrieval of remote vector entries. +! +! There are three possible exchange algorithms: +! 1. Use MPI_Alltoallv +! 2. Use psb_simple_a2av +! 3. Use psb_simple_triad_a2av +! Default choice is 3. The MPI variant has proved to be inefficient; +! that is because it is not persistent, therefore you pay the initialization price +! every time, and it is not optimized for a sparse communication pattern, +! most MPI implementations assume that all communications are non-empty. +! The PSB_SIMPLE variants reuse the same communicator, and go for a simplistic +! sequence of sends/receive that is quite efficient for a sparse communication +! pattern. To be refined/reviewed in the future to compare with neighbour +! persistent collectives. +! +! Arguments: +! desc_a - type(psb_desc_type). The communication descriptor. +! info - integer. Return code +! rowcnv - logical Should row/col indices be converted +! colcnv - logical to/from global numbering when sent/received? +! default is .TRUE. +! rowscale - logical Should row/col indices on output be remapped +! colscale - logical from MIN:MAX to 1:(MAX-MIN+1) ? +! default is .FALSE. +! (commmon use is ROWSCALE=.TRUE., COLSCALE=.FALSE.) +! data - integer Which index list in desc_a should be used to retrieve +! rows, default psb_comm_halo_ +! psb_comm_halo_ use halo_index +! psb_comm_ext_ use ext_index +! psb_comm_ovrl_ DISABLED for this routine. +! +subroutine psb_z_remote_vect(n,v,iv,desc_a,x,ix, info) + use psb_base_mod, psb_protect_name => psb_z_remote_vect + +#ifdef MPI_MOD + use mpi +#endif + Implicit None +#ifdef MPI_H + include 'mpif.h' +#endif + integer(psb_ipk_), intent(in) :: n + complex(psb_dpk_), Intent(in) :: v(:) + integer(psb_lpk_), Intent(in) :: iv(:) + type(psb_desc_type),intent(in) :: desc_a + complex(psb_dpk_), allocatable, intent(out) :: x(:) + integer(psb_lpk_), allocatable, intent(out) :: ix(:) + integer(psb_ipk_), intent(out) :: info + ! ...local scalars.... + type(psb_ctxt_type) :: ctxt + integer(psb_ipk_) :: np,me + integer(psb_ipk_) :: counter, proc, i, & + & j, idxs,idxr, k, iszs, iszr + integer(psb_ipk_) :: nrcvs, nsnds + integer(psb_mpk_) :: icomm, minfo + integer(psb_mpk_), allocatable :: brvindx(:), & + & rvsz(:), bsdindx(:), sdsz(:), sdsi(:), rvsi(:) + integer(psb_lpk_), allocatable :: lsnd(:) + complex(psb_dpk_), allocatable :: valsnd(:) + integer(psb_ipk_), allocatable :: iprc(:) + integer(psb_ipk_) :: debug_level, debug_unit, err_act + character(len=20) :: name, ch_err + + info=psb_success_ + name='psb_z_remote_vect' + 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() + + ctxt = desc_a%get_context() + icomm = desc_a%get_mpic() + + Call psb_info(ctxt, me, np) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),': Start' + + Allocate(rvsz(np),sdsz(np),sdsi(np),rvsi(np),brvindx(np+1),& + & bsdindx(np+1), stat=info) + + if (info /= psb_success_) then + info=psb_err_alloc_dealloc_ + call psb_errpush(info,name) + goto 9999 + end if + + call desc_a%indxmap%fnd_owner(iv(1:n),iprc,info) + + icomm = desc_a%get_mpic() + sdsz(:) = 0 + rvsz(:) = 0 + sdsi(:) = 0 + rvsi(:) = 0 + brvindx(:) = 0 + bsdindx(:) = 0 + counter = 1 + idxs = 0 + idxr = 0 + do i=1,n + if (iprc(i) >=0) then + sdsz(iprc(i)+1) = sdsz(iprc(i)+1) +1 + else + write(0,*)me,name,' Error from fnd_owner: ',iprc(i) + end if + end do + call mpi_alltoall(sdsz,1,psb_mpi_mpk_,& + & rvsz,1,psb_mpi_mpk_,icomm,minfo) + if (minfo /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='mpi_alltoall') + goto 9999 + end if + !write(0,*)me,name,' sdsz ',sdsz(:),' rvsz:',rvsz(:) + nsnds = count(sdsz /= 0) + nrcvs = count(rvsz /= 0) + 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 + + iszs = sum(sdsz) + iszr = sum(rvsz) + call psb_realloc(iszs,lsnd,info) + if (info == 0) call psb_realloc(iszs,valsnd,info) + if (info == 0) call psb_realloc(iszr,x,info) + if (info == 0) call psb_realloc(iszr,ix,info) + if (info /= psb_success_) then + info=psb_err_from_subroutine_ + call psb_errpush(info,name,a_err='realloc') + goto 9999 + end if + do k=1, n + proc = iprc(k) + sdsi(proc+1) = sdsi(proc+1) + 1 + lsnd(bsdindx(proc+1)+sdsi(proc+1)) = iv(k) + valsnd(bsdindx(proc+1)+sdsi(proc+1)) = v(k) + end do + do proc=0,np-1 + if (sdsi(proc+1) /= sdsz(proc+1)) & + & write(0,*) me,name,'Send mismacth ',sdsi(proc+1),sdsz(proc+1) + end do + + select case(psb_get_sp_a2av_alg()) + case(psb_sp_a2av_smpl_triad_,psb_sp_a2av_smpl_v_) + call psb_simple_a2av(valsnd,sdsz,bsdindx,& + & x,rvsz,brvindx,ctxt,info) + if (info == psb_success_) call psb_simple_a2av(lsnd,sdsz,bsdindx,& + & ix,rvsz,brvindx,ctxt,info) + case(psb_sp_a2av_mpi_) + + call mpi_alltoallv(valsnd,sdsz,bsdindx,psb_mpi_c_dpk_,& + & x,rvsz,brvindx,psb_mpi_c_dpk_,icomm,minfo) + if (minfo == mpi_success) & + & call mpi_alltoallv(lsnd,sdsz,bsdindx,psb_mpi_lpk_,& + & ix,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 + + Deallocate(brvindx,bsdindx,rvsz,sdsz,& + & lsnd,valsnd,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(ctxt,err_act) + + return + +End Subroutine psb_z_remote_vect diff --git a/base/tools/psb_zasb.f90 b/base/tools/psb_zasb.f90 index 11c08d6d..decbfdec 100644 --- a/base/tools/psb_zasb.f90 +++ b/base/tools/psb_zasb.f90 @@ -113,7 +113,7 @@ subroutine psb_zasb_vect(x, desc_a, info, mold, scratch) integer(psb_ipk_) :: nrmv, nx, i nrmv = x%get_nrmv() - call psb_z_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) + call psb_remote_vect(nrmv,x%rmtv,x%rmidx,desc_a,vx,lvx,info) nx = size(vx) call psb_realloc(nx,ivx,info) call desc_a%g2l(lvx,ivx,info,owned=.true.)