fdurastante/psblas3
1
0
Fork 0
mirror of https://github.com/sfilippone/psblas3.git synced 2026-10-19 20:32:12 +00:00
Code Issues Packages Projects Releases Wiki Activity

Merge branch 'unify_aggr_bld' into development

Browse Source
This commit is contained in:
Salvatore Filippone
2020-06-04 10:39:03 +02:00
parent 78f42c8e22 76ff6c824c
commit 1493be1de6
16 changed files with 887 additions and 311 deletions
Download Patch File Download Diff File Expand all files Collapse all files
base/internals/Makefile
+4 -3
View File
@@ -1,11 +1,12 @@
include ../../Make.inc
FOBJS = psi_compute_size.o psi_crea_bnd_elem.o psi_crea_index.o \
psi_crea_ovr_elem.o psi_bld_tmpovrl.o psi_dl_check.o \
psi_crea_ovr_elem.o psi_bld_tmpovrl.o \
psi_bld_tmphalo.o psi_sort_dl.o \
psi_indx_map_fnd_owner.o \
psi_desc_impl.o psi_hash_impl.o psi_list_search.o psi_srtlist.o
psi_desc_impl.o psi_hash_impl.o psi_srtlist.o \
psi_bld_glb_dep_list.o psi_xtr_loc_dl.o
#psi_list_search.o psi_dl_check.o
MPFOBJS = psi_desc_index.o psi_extrct_dl.o psi_fnd_owner.o psi_a2a_fnd_owner.o \
psi_graph_fnd_owner.o psi_adjcncy_fnd_owner.o psi_symm_dep_list.o
base/internals/psi_bld_glb_dep_list.F90
+215
View File
@@ -0,0 +1,215 @@
!
! 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.
!
!
subroutine psi_i_bld_glb_dep_list(ictxt,loc_dl,length_dl,dep_list,dl_lda,info)
use psi_mod, psb_protect_name => psi_i_bld_glb_dep_list
#ifdef MPI_MOD
use mpi
#endif
use psb_penv_mod
use psb_const_mod
use psb_error_mod
use psb_desc_mod
use psb_sort_mod
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
! ....scalar parameters...
integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(out) :: dl_lda
integer(psb_ipk_), intent(in) :: loc_dl(:), length_dl(0:)
integer(psb_ipk_), allocatable, intent(out) :: dep_list(:,:)
integer(psb_ipk_), intent(out) :: info
! .....local arrays....
integer(psb_ipk_) :: int_err(5)
! .....local scalars...
integer(psb_ipk_) :: i, proc,j,err_act
integer(psb_ipk_) :: err
integer(psb_ipk_) :: debug_level, debug_unit
integer(psb_mpk_) :: iictxt, icomm, me, np, minfo
logical, parameter :: dist_symm_list=.false., print_dl=.false.
character name*20
name='psi_bld_glb_dep_list'
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
iictxt = ictxt
info = psb_success_
call psb_info(iictxt,me,np)
dl_lda = length_dl(me)
call psb_max(iictxt, dl_lda)
if (debug_level >= psb_debug_inner_) &
& write(debug_unit,*) me,' ',trim(name),': Dep_list length ',length_dl(me),dl_lda
dl_lda = max(dl_lda,1)
allocate(dep_list(dl_lda,0:np),stat=info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='Allocate')
goto 9999
end if
icomm = psb_get_mpi_comm(iictxt)
call mpi_allgather(loc_dl,dl_lda,psb_mpi_ipk_,&
& dep_list,dl_lda,psb_mpi_ipk_,icomm,minfo)
info = minfo
if (info /= psb_success_) then
info=psb_err_internal_error_
goto 9999
endif
if (print_dl) then
if (me == 0) then
write(0,*) ' Dep_list '
do i=0,np-1
j = length_dl(i)
write(0,*) 'Proc ',i,':',dep_list(1:j,i)
end do
flush(0)
end if
call psb_barrier(ictxt)
end if
call psb_erractionrestore(err_act)
return
9999 continue
call psb_errpush(info,name,i_err=int_err)
call psb_error_handler(err_act)
return
end subroutine psi_i_bld_glb_dep_list
subroutine psi_i_bld_glb_csr_dep_list(ictxt,loc_dl,length_dl,c_dep_list,dl_ptr,info)
use psi_mod, psb_protect_name => psi_i_bld_glb_csr_dep_list
#ifdef MPI_MOD
use mpi
#endif
use psb_penv_mod
use psb_const_mod
use psb_error_mod
use psb_desc_mod
use psb_sort_mod
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
! ....scalar parameters...
integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(in) :: loc_dl(:), length_dl(0:)
integer(psb_ipk_), allocatable, intent(out) :: c_dep_list(:), dl_ptr(:)
integer(psb_ipk_), intent(out) :: info
! .....local arrays....
integer(psb_ipk_) :: int_err(5)
! .....local scalars...
integer(psb_ipk_) :: i, proc,j,err_act, length, myld
integer(psb_ipk_) :: err
integer(psb_ipk_) :: debug_level, debug_unit
integer(psb_mpk_) :: iictxt, icomm, me, np, minfo
logical, parameter :: dist_symm_list=.false., print_dl=.false.
character name*20
name='psi_bld_glb_csr_dep_list'
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
iictxt = ictxt
info = psb_success_
call psb_info(iictxt,me,np)
myld = length_dl(me)
length = sum(length_dl(0:np-1))
allocate(dl_ptr(0:np),stat=info)
dl_ptr(0) = 0
do i=1, np
dl_ptr(i) = dl_ptr(i-1) + length_dl(i-1)
end do
if (length /= dl_ptr(np)) then
write(0,*) me,trim(name),' Inconsistency: ',length,dl_ptr(np)
end if
if (debug_level >= psb_debug_inner_) &
& write(debug_unit,*) me,' ',trim(name),': Dep_list length ',length_dl(me)
allocate(c_dep_list(length),stat=info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='Allocate')
goto 9999
end if
icomm = psb_get_mpi_comm(iictxt)
call mpi_allgatherv(loc_dl,myld,psb_mpi_ipk_,&
& c_dep_list,length_dl,dl_ptr,psb_mpi_ipk_,icomm,minfo)
info = minfo
if (info /= psb_success_) then
info=psb_err_internal_error_
goto 9999
endif
dl_ptr = dl_ptr + 1
if (print_dl) then
if (me == 0) then
write(0,*) ' Dep_list '
do i=0,np-1
write(0,*) 'Proc ',i,':',c_dep_list(dl_ptr(i):dl_ptr(i+1)-1)
end do
flush(0)
end if
call psb_barrier(ictxt)
end if
call psb_erractionrestore(err_act)
return
9999 continue
call psb_errpush(info,name,i_err=int_err)
call psb_error_handler(err_act)
return
end subroutine psi_i_bld_glb_csr_dep_list
base/internals/psi_crea_index.f90
+56 -26
View File
@@ -64,9 +64,10 @@ subroutine psi_i_crea_index(desc_a,index_in,index_out,nxch,nsnd,nrcv,info)
integer(psb_ipk_), allocatable, intent(inout) :: index_out(:)
! ....local scalars...
integer(psb_ipk_) :: ictxt, me, np, mode, err_act, dl_lda
integer(psb_ipk_) :: ictxt, me, np, mode, err_act, dl_lda, ldl
! ...parameters...
integer(psb_ipk_), allocatable :: dep_list(:,:), length_dl(:)
integer(psb_ipk_), allocatable :: dep_list(:,:), length_dl(:), loc_dl(:), c_dep_list(:), dl_ptr(:)
integer(psb_ipk_) :: dlmax, dlavg
integer(psb_ipk_),parameter :: root=psb_root_,no_comm=-1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
@@ -109,50 +110,69 @@ subroutine psi_i_crea_index(desc_a,index_in,index_out,nxch,nsnd,nrcv,info)
mode = 1
if (do_timings) call psb_tic(idx_phase1)
call psi_extract_dep_list(ictxt,&
call psi_extract_loc_dl(ictxt,&
& desc_a%is_bld(), desc_a%is_upd(),&
& index_in, dep_list,length_dl,dl_lda,mode,info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='extrct_dl')
goto 9999
end if
& index_in, loc_dl,length_dl,info)
dlmax = maxval(length_dl(:))
dlavg = (sum(length_dl(:))+np-1)/np
!!$ if ((dlmax>0).and.(me==0)) write(0,*) 'Dependency list : max:',dlmax,&
!!$ & ' avg:',dlavg, choose_sorting(dlmax,dlavg,np)
if (debug_level >= psb_debug_inner_) &
& write(debug_unit,*) me,' ',trim(name),': from extract_dep_list',&
& me,length_dl(0),index_in(1), ':',dep_list(:length_dl(me),me)
! ...now process root contains dependence list of all processes...
if (debug_level >= psb_debug_inner_) &
& write(debug_unit,*) me,' ',trim(name),': root sorting dep list'
if (do_timings) call psb_toc(idx_phase1)
if (do_timings) call psb_tic(idx_phase2)
call psi_dl_check(dep_list,dl_lda,np,length_dl)
if (choose_sorting(dlmax,dlavg,np)) then
call psi_bld_glb_dep_list(ictxt,&
& loc_dl,length_dl,c_dep_list,dl_ptr,info)
if (info /= 0) then
write(0,*) me,trim(name),' From bld_glb_list ',info
end if
!!$ call psi_dl_check(dep_list,dl_lda,np,length_dl)
!!$
!!$ ! ....now i can sort dependency lists.
call psi_sort_dl(dl_ptr,c_dep_list,length_dl,ictxt,info)
if (info /= 0) then
write(0,*) me,trim(name),' From sort_dl ',info
end if
ldl = length_dl(me)
loc_dl = c_dep_list(dl_ptr(me):dl_ptr(me)+ldl-1)
! ....now i can sort dependency lists.
call psi_sort_dl(dep_list,length_dl,np,info)
if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psi_sort_dl')
goto 9999
!!$ if(info /= psb_success_) then
!!$ call psb_errpush(psb_err_from_subroutine_,name,a_err='psi_sort_dl')
!!$ goto 9999
!!$ end if
else
! Do nothing
ldl = length_dl(me)
loc_dl = loc_dl(1:ldl)
end if
if (do_timings) call psb_toc(idx_phase2)
if (do_timings) call psb_tic(idx_phase3)
if(debug_level >= psb_debug_inner_)&
& write(debug_unit,*) me,' ',trim(name),': calling psi_desc_index'
& write(debug_unit,*) me,' ',trim(name),': calling psi_desc_index',ldl,':',loc_dl(1:ldl)
! Do the actual format conversion.
call psi_desc_index(desc_a,index_in,dep_list(1:,me),&
& length_dl(me),nsnd,nrcv, index_out,info)
call psi_desc_index(desc_a,index_in,loc_dl,ldl,nsnd,nrcv,index_out,info)
if(debug_level >= psb_debug_inner_) &
& write(debug_unit,*) me,' ',trim(name),': out of psi_desc_index',&
& size(index_out)
nxch = length_dl(me)
nxch = ldl
if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psi_desc_index')
goto 9999
end if
if (do_timings) call psb_toc(idx_phase3)
deallocate(dep_list,length_dl)
if (allocated(dep_list)) deallocate(dep_list,stat=info)
if ((info==0).and.allocated(length_dl)) deallocate(length_dl,stat=info)
if (info /= 0) then
info = psb_err_alloc_dealloc_
goto 9999
end if
if(debug_level >= psb_debug_inner_) &
& write(debug_unit,*) me,' ',trim(name),': done'
@@ -162,4 +182,14 @@ subroutine psi_i_crea_index(desc_a,index_in,index_out,nxch,nsnd,nrcv,info)
9999 call psb_error_handler(ictxt,err_act)
return
contains
function choose_sorting(dlmax,dlavg,np) result(val)
implicit none
integer(psb_ipk_), intent(in) :: dlmax,dlavg,np
logical :: val
val = .not.(((dlmax>(26*4)).or.((dlavg>=(26*2)).and.(np>=128))))
!val = .true.
end function choose_sorting
end subroutine psi_i_crea_index
base/internals/psi_desc_index.F90
+20 -2
View File
@@ -105,6 +105,7 @@ subroutine psi_i_desc_index(desc,index_in,dep_list,&
use mpi
#endif
use psb_penv_mod
use psb_timers_mod
use psi_mod, psb_protect_name => psi_i_desc_index
implicit none
#ifdef MPI_H
@@ -136,6 +137,8 @@ subroutine psi_i_desc_index(desc,index_in,dep_list,&
& idxr, idxs, iszs, iszr, nesd, nerv, ixp, idx
integer(psb_mpk_) :: icomm, minfo
logical, parameter :: do_timings=.true., oldstyle=.false., debug=.false.
integer(psb_ipk_), save :: idx_phase1=-1, idx_phase2=-1, idx_phase3=-1, idx_phase4=-1
logical, parameter :: usempi=.false.
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
@@ -159,13 +162,22 @@ subroutine psi_i_desc_index(desc,index_in,dep_list,&
write(debug_unit,*) me,' ',trim(name),': start'
call psb_barrier(ictxt)
endif
if ((do_timings).and.(idx_phase1==-1)) &
& idx_phase1 = psb_get_timer_idx("I_DSC_IDX: phase1 ")
if ((do_timings).and.(idx_phase2==-1)) &
& idx_phase2 = psb_get_timer_idx("I_DSC_IDX: phase2 ")
if ((do_timings).and.(idx_phase3==-1)) &
& idx_phase3 = psb_get_timer_idx("I_DSC_IDX: phase3 ")
if ((do_timings).and.(idx_phase4==-1)) &
& idx_phase4 = psb_get_timer_idx("I_DSC_IDX: phase4 ")
if (do_timings) call psb_tic(idx_phase1)
!
! first, find out the sizes to be exchanged.
! note: things marked here as sndbuf/rcvbuf (for mpi) corresponds to things
! to be received/sent (in the final psblas descriptor).
! be careful of the inversion
!
!
allocate(sdsz(np),rvsz(np),bsdindx(np),brvindx(np),stat=info)
if(info /= psb_success_) then
info=psb_err_alloc_dealloc_
@@ -239,6 +251,8 @@ subroutine psi_i_desc_index(desc,index_in,dep_list,&
call psb_errpush(info,name)
goto 9999
end if
if (do_timings) call psb_toc(idx_phase1)
if (do_timings) call psb_tic(idx_phase2)
!
! Second build the lists of requests
@@ -285,6 +299,8 @@ subroutine psi_i_desc_index(desc,index_in,dep_list,&
brvindx(proc+1) = idxr
idxr = idxr + rvsz(proc+1)
end do
if (do_timings) call psb_toc(idx_phase2)
if (do_timings) call psb_tic(idx_phase3)
if (usempi) then
call mpi_alltoallv(sndbuf,sdsz,bsdindx,psb_mpi_lpk_,&
@@ -357,6 +373,8 @@ subroutine psi_i_desc_index(desc,index_in,dep_list,&
end if
end if
if (do_timings) call psb_toc(idx_phase3)
if (do_timings) call psb_tic(idx_phase4)
!
! at this point we can finally build the output desc_index. beware
@@ -386,7 +404,7 @@ subroutine psi_i_desc_index(desc,index_in,dep_list,&
call psb_errpush(info,name)
goto 9999
end if
if (do_timings) call psb_toc(idx_phase4)
if (debug_level >= psb_debug_inner_) then
write(debug_unit,*) me,' ',trim(name),': done'
call psb_barrier(ictxt)
base/internals/psi_dl_check.f90
-95
View File
@@ -1,95 +0,0 @@
!
! 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: psi_dl_check.f90
!
! Subroutine: psi_dl_check
! Make sure a dependency list is symmetric, i.e. if process i depends on j
! then process j should depend on i (even if the data to be sent in one of the
! directions happens to be empty)
!
! Arguments:
! dep_list(:,:) - integer Initial dependency lists
! dl_lda - integer Allocated size of dep_list
! np - integer Total number of processes.
! length_dl(:) - integer Items in dependency lists; updated on
! exit
!
subroutine psi_i_dl_check(dep_list,dl_lda,np,length_dl)
use psi_mod, psb_protect_name => psi_i_dl_check
use psb_const_mod
use psb_desc_mod
implicit none
integer(psb_ipk_) :: np,dl_lda,length_dl(0:np)
integer(psb_ipk_) :: dep_list(dl_lda,0:np)
! locals
integer(psb_ipk_) :: proc, proc2, i, j
! ...if j is in dep_list of process i
! and i is not in dep_list of process j
! fix it.
do proc=0,np-1
i=1
outer: do
if (i >length_dl(proc)) exit outer
proc2=dep_list(i,proc)
if ((proc2 /= -1).and.(proc2 /= proc)) then
! ...search proc in proc2's dep_list....
j=1
p2loop:do
if (j > length_dl(proc2)) exit p2loop
if (dep_list(j,proc2) == proc) exit p2loop
j=j+1
enddo p2loop
if (j > length_dl(proc2)) then
! ...add proc to proc2 s dep_list.....',proc,proc2
length_dl(proc2) = length_dl(proc2)+1
if (length_dl(proc2) > size(dep_list,1)) then
write(psb_err_unit,*)'error in dl_check', proc2,proc,&
& length_dl(proc2),'>',size(dep_list,1)
endif
dep_list(length_dl(proc2),proc2) = proc
else if (dep_list(j,proc2) /= proc) then
write(psb_err_unit,*) 'PSI_DL_CHECK This should not happen!!! ',&
& j,proc2,dep_list(j,proc2),proc
endif
endif
i=i+1
enddo outer
enddo
end subroutine psi_i_dl_check
base/internals/psi_extrct_dl.F90
+31 -1
View File
@@ -127,6 +127,7 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
use psb_error_mod
use psb_desc_mod
use psb_sort_mod
use psb_timers_mod
implicit none
#ifdef MPI_H
include 'mpif.h'
@@ -147,7 +148,9 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
integer(psb_ipk_) :: err
integer(psb_ipk_) :: debug_level, debug_unit
integer(psb_mpk_) :: iictxt, icomm, me, np, minfo
logical, parameter :: dist_symm_list=.false., print_dl=.false.
logical, parameter :: dist_symm_list=.false., print_dl=.false., profile=.true.
logical, parameter :: do_timings=.false.
integer(psb_ipk_), save :: idx_phase1=-1, idx_phase2=-1, idx_phase3=-1
character name*20
name='psi_extrct_dl'
@@ -156,8 +159,16 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
debug_level = psb_get_debug_level()
iictxt = ictxt
info = psb_success_
if ((do_timings).and.(idx_phase1==-1)) &
& idx_phase1 = psb_get_timer_idx("PSI_XTR_DL: phase1 ")
if ((do_timings).and.(idx_phase2==-1)) &
& idx_phase2 = psb_get_timer_idx("PSI_XTR_DL: phase2")
!!$ if ((do_timings).and.(idx_phase3==-1)) &
!!$ & idx_phase3 = psb_get_timer_idx("PSI_XTR_DL: phase3")
call psb_info(iictxt,me,np)
if (do_timings) call psb_tic(idx_phase1)
allocate(itmp(2*np+1),length_dl(0:np),stat=info)
if (info /= psb_success_) then
info=psb_err_alloc_dealloc_
@@ -252,11 +263,17 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
endif
length_dl(me)=pointer_dep_list-1
if (do_timings) call psb_toc(idx_phase1)
if (do_timings) call psb_tic(idx_phase2)
if (dist_symm_list) then
call psb_realloc(length_dl(me),itmp,info)
call psi_symm_dep_list(itmp,ictxt,info)
dl_lda = max(size(itmp),1)
call psb_max(iictxt, dl_lda)
if (debug_level >= psb_debug_inner_) &
& write(debug_unit,*) me,' ',trim(name),': Dep_list length ',length_dl(me),dl_lda
call psb_realloc(dl_lda,itmp,info)
! dl_lda = min(np,2*dl_lda)
allocate(dep_list(dl_lda,0:np),stat=info)
@@ -282,6 +299,8 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
dl_lda = max(length_dl(me),1)
call psb_max(iictxt, dl_lda)
if (debug_level >= psb_debug_inner_) &
& write(debug_unit,*) me,' ',trim(name),': Dep_list length ',length_dl(me),dl_lda
allocate(dep_list(dl_lda,0:np),stat=info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='Allocate')
@@ -362,6 +381,17 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
end if
call psb_barrier(ictxt)
end if
if (do_timings) call psb_toc(idx_phase2)
if ((profile).and.(me==0)) then
block
integer(psb_ipk_) :: dlmax, dlavg
dlmax = maxval(length_dl(:))
dlavg = (sum(length_dl(:))+np-1)/np
if (dlmax>0) write(0,*) 'Dependency list : max:',dlmax,&
& ' avg:',dlavg, ((dlmax>np/2).or.((dlavg>=np/4).and.(np>128)))
end block
end if
call psb_erractionrestore(err_act)
return
base/internals/psi_graph_fnd_owner.F90
+27 -27
View File
@@ -47,7 +47,7 @@
! This is the method to find out who owns a set of indices.
! In principle we could do the following:
! 1. Do an allgatherv of IDX
! 2. For each of the collected indices figure if current proces owns it
! 2. For each of the collected indices figure out if current proces owns it
! 3. Scatter the results
! 4. Loop through the answers
! This method is guaranteed to find the owner, unless an input index has
@@ -101,8 +101,8 @@ subroutine psi_graph_fnd_owner(idx,iprc,idxmap,info)
integer(psb_ipk_), allocatable :: tprc(:), tsmpl(:), ladj(:)
integer(psb_mpk_) :: icomm, minfo, iictxt
integer(psb_ipk_) :: i,n_row,n_col,err_act,ip,j,ipnt, nsampl_out,&
& nv, n_answers, nreqst, nsampl_in, locr_max, &
& nreqst_max, nadj, maxspace, mxnsin
& nv, n_answers, nqries, nsampl_in, locr_max, &
& nqries_max, nadj, maxspace, mxnsin
integer(psb_lpk_) :: mglob, ih
integer(psb_ipk_) :: ictxt,np,me, nresp
integer(psb_ipk_), parameter :: nt=4
@@ -165,22 +165,22 @@ subroutine psi_graph_fnd_owner(idx,iprc,idxmap,info)
! This makes ladj allocated with size 0 if needed, as opposed to unallocated
call psb_realloc(nadj,ladj,info)
!
! Throughout the subroutine, nreqst is the number of local inquiries
! Throughout the subroutine, nqries is the number of local inquiries
! that have not been answered yet
!
nreqst = nv - n_answers
nreqst_max = nreqst
nqries = nv - n_answers
nqries_max = nqries
!
! Choice of maxspace should be adjusted to account for a default
! "sensible" size and/or a user-specified value
!
tmpv(1) = nadj
tmpv(2) = nreqst_max
tmpv(2) = nqries_max
tmpv(3) = n_row
tmpv(4) = psb_cd_get_maxspace()
call psb_max(ictxt,tmpv)
nreqst_max = tmpv(2)
nqries_max = tmpv(2)
locr_max = tmpv(3)
maxspace = nt*locr_max
if (tmpv(4) > 0) maxspace = min(maxspace,tmpv(4))
@@ -192,21 +192,21 @@ subroutine psi_graph_fnd_owner(idx,iprc,idxmap,info)
! Do a preliminary run on the user-defined adjacency lists
!
if (trace.and.(me == 0)) write(0,*) ' Initial sweep on user-defined topology'
if (debugsz) write(0,*) me,' Initial sweep on user-defined topology',nreqst
nsampl_in = min(nreqst,max(1,(maxspace+max(1,nadj)-1))/(max(1,nadj)))
if (debugsz) write(0,*) me,' Initial sweep on user-defined topology',nqries
nsampl_in = min(nqries,max(1,(maxspace+max(1,nadj)-1))/(max(1,nadj)))
call psi_adj_fnd_sweep(idx,iprc,ladj,idxmap,nsampl_in,n_answers)
call idxmap%xtnd_p_adjcncy(ladj)
nreqst = nv - n_answers
nreqst_max = nreqst
call psb_max(ictxt,nreqst_max)
if (trace.and.(me == 0)) write(0,*) ' After initial sweep:',nreqst_max
if (debugsz) write(0,*) me,' After sweep on user-defined topology',nreqst_max
nqries = nv - n_answers
nqries_max = nqries
call psb_max(ictxt,nqries_max)
if (trace.and.(me == 0)) write(0,*) ' After initial sweep:',nqries_max
if (debugsz) write(0,*) me,' After sweep on user-defined topology',nqries_max
end if
if (do_timings) call psb_toc(idx_sweep0)
fnd_owner_loop: do while (nreqst_max>0)
fnd_owner_loop: do while (nqries_max>0)
if (do_timings) call psb_tic(idx_loop_a2a)
if (debugsz) write(0,*) me,' fnd_owner_loop',nreqst_max
if (debugsz) write(0,*) me,' fnd_owner_loop',nqries_max
!
! The basic idea of this loop is to alternate between
! searching through all processes and searching
@@ -215,8 +215,8 @@ subroutine psi_graph_fnd_owner(idx,iprc,idxmap,info)
! 1. Select a sample such that the total size is <= maxspace
! sample query is then sent to all processes
!
! if (trace.and.(me == 0)) write(0,*) 'Looping in graph_fnd_owner: ', nreqst_max
nsampl_in = psb_cd_get_samplesize()
! if (trace.and.(me == 0)) write(0,*) 'Looping in graph_fnd_owner: ', nqries_max
nsampl_in = nqries
nsampl_in = min(max(1,(maxspace+np-1)/np),nsampl_in)
!
! Choose a sample, should it be done in this simplistic way?
@@ -236,13 +236,13 @@ subroutine psi_graph_fnd_owner(idx,iprc,idxmap,info)
! We might have padded when looking for owners, so the actual samples
! could be less than they appear. Should be explained better.
!
nsampl_in = min(nreqst,nsampl_in)
nsampl_in = min(nqries,nsampl_in)
call psi_cpy_out(iprc,tprc,tsmpl,nsampl_in,nsampl_out)
if (nsampl_out /= nsampl_in) then
write(0,*) me,'Warning: indices not found by a2a_fnd_owner ',nsampl_out,nsampl_in
end if
n_answers = n_answers + nsampl_out
nreqst = nv - n_answers
nqries = nv - n_answers
!
! 3. Extract the resulting adjacency list and add it to the
! indxmap;
@@ -259,18 +259,18 @@ subroutine psi_graph_fnd_owner(idx,iprc,idxmap,info)
! Need to set up a proper loop here to have a complete
! sweep over the input vector. Done inside adj_fnd_sweep.
!
!!$ write(0,*) me,' After a2a ',nreqst
nsampl_in = min(nreqst,max(1,(maxspace+max(1,nadj)-1))/(max(1,nadj)))
!!$ write(0,*) me,' After a2a ',nqries
nsampl_in = min(nqries,max(1,(maxspace+max(1,nadj)-1))/(max(1,nadj)))
mxnsin = nsampl_in
call psb_max(ictxt,mxnsin)
!!$ write(0,*) me, ' mxnsin ',mxnsin
if (mxnsin>0) call psi_adj_fnd_sweep(idx,iprc,ladj,idxmap,nsampl_in,n_answers)
call idxmap%xtnd_p_adjcncy(ladj)
nreqst = nv - n_answers
nreqst_max = nreqst
call psb_max(ictxt,nreqst_max)
if (trace.and.(me == 0)) write(0,*) ' fnd_owner_loop remaining:',nreqst_max
nqries = nv - n_answers
nqries_max = nqries
call psb_max(ictxt,nqries_max)
if (trace.and.(me == 0)) write(0,*) ' fnd_owner_loop remaining:',nqries_max
if (do_timings) call psb_toc(idx_loop_neigh)
end do fnd_owner_loop
base/internals/psi_list_search.f90
-58
View File
@@ -1,58 +0,0 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006, 2010, 2015, 2017
! Salvatore Filippone
! Alfredo Buttari CNRS-IRIT, Toulouse
!
! 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.
!
!
integer function psi_list_search(list,lenght_list,elem)
use psb_const_mod
implicit none
!returns position of elem in a array list
!of lenght lenght_list, if this element does not exist
!returns -1
integer(psb_ipk_) :: list(*)
integer(psb_ipk_) :: lenght_list
integer(psb_ipk_) :: elem
integer(psb_ipk_) :: i
i=1
do while ((i.le.lenght_list).and.(list(i).ne.elem))
i=i+1
enddo
if (i.le.lenght_list) then
if (list(i).eq.elem) then
psi_list_search=i
else
psi_list_search=-1
endif
else
psi_list_search=-1
endif
end function psi_list_search
base/internals/psi_sort_dl.f90
+199 -55
View File
@@ -29,65 +29,209 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psi_i_sort_dl(dep_list,l_dep_list,np,info)
!
! interface between former sort_dep_list subroutine
! and new srtlist
!
use psi_mod, psb_protect_name => psi_i_sort_dl
!**********************************************************************
! *
! The communication step among processors at each *
! matrix-vector product is a variable all-to-all *
! collective communication that we reimplement *
! in terms of point-to-point communications. *
! The data in input is a list of dependencies: *
! for each node a list of all the nodes it has to *
! communicate with. The lists are guaranteed to be *
! symmetric, i.e. for each pair (I,J) there is a *
! pair (J,I). The idea is to organize the ordering *
! so that at each communication step as many *
! processors as possible are communicating at the *
! same time, i.e. a step is defined by the fact *
! that all edges (I,J) in it have no common node. *
! *
! Formulation of the problem is: *
! Given an undirected graph (forest): *
! Find the shortest series of steps to cancel all *
! graph edges, where at each step all edges belonging *
! to a matching in the graph are canceled. *
! *
! An obvious lower bound to the optimum number of steps *
! is the largest degree of any node in the graph. *
! *
! The algorithm proceeds as follows: *
! 1. Build a list of all edges, e.g. copy the *
! dependencies lists keeping only (I,J) with I<J *
! 2. Compute an auxiliary vector with the degree of *
! each node of the graph. *
! 3. While there are edges in the graph do: *
! 4. Weight the edges with the sum of the degrees *
! of their nodes and sort them into descending order *
! 5. Scan the list of edges; if neither node of the *
! edge has been marked yet, cancel the edge and mark *
! the two nodes *
! 6. If no edge was chosen but the graph is nonempty *
! raise an error condition *
! 7. Queue the edges in the matchin to the output *
! sequence; *
! 8. Decrease by 1 the degree of all marked nodes, *
! then clear all marks *
! 9. Cycle to 3. *
! 10. For each node: scan the edge sequence; if an *
! edge has the node as an endpoint, queue the other *
! node in the dependency list for the current one *
! *
!**********************************************************************
subroutine psi_i_csr_sort_dl(dl_ptr,c_dep_list,l_dep_list,ictxt,info)
use psi_mod, psb_protect_name => psi_i_csr_sort_dl
use psb_const_mod
use psb_error_mod
use psb_sort_mod
implicit none
integer(psb_ipk_) :: np,dep_list(:,:), l_dep_list(:)
integer(psb_ipk_) :: idg, iupd, idgp, iedges, iidx, iich,ndgmx, isz, err_act
integer(psb_ipk_) :: i, info
integer(psb_ipk_), allocatable :: work(:)
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name='psi_sort_dl'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
integer(psb_ipk_), intent(inout) :: c_dep_list(:), dl_ptr(0:), l_dep_list(0:)
integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(out) :: info
! Local variables
integer(psb_ipk_), allocatable :: dg(:), dgp(:),&
& idx(:), upd(:), edges(:,:), ich(:)
integer(psb_ipk_) :: i, j, nedges, ip1, ip2, nch, ip, iedge,&
& i1, ix, ist, iswap(2)
logical :: internal_error
integer(psb_ipk_) :: me, np
info = 0
call psb_info(ictxt,me,np)
nedges = size(c_dep_list)
info = psb_success_
ndgmx = 0
do i=1,np
ndgmx = ndgmx + l_dep_list(i)
if (debug_level >= psb_debug_inner_)&
& write(debug_unit,*) name,': ',i,l_dep_list(i)
enddo
idg = 1
iupd = idg+np
idgp = iupd+np
iedges = idgp + ndgmx
iidx = iedges + 2*ndgmx
iich = iidx + ndgmx
isz = iich + ndgmx
if (debug_level >= psb_debug_inner_)&
& write(debug_unit,*) name,': ndgmx ',ndgmx,isz
allocate(work(isz))
! call srtlist(dep_list, dl_lda, l_dep_list, np, info)
call srtlist(dep_list,size(dep_list,1,kind=psb_ipk_),l_dep_list,np,work(idg),&
& work(idgp),work(iupd),work(iedges),work(iidx),work(iich),info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='srtlist')
goto 9999
endif
allocate(dg(0:np-1),dgp(nedges),edges(2,nedges),upd(0:np-1),&
& idx(nedges),ich(nedges),stat = info)
deallocate(work)
call psb_erractionrestore(err_act)
return
if (info /= 0) then
info = -9
return
end if
!
! 1. Compute an auxiliary vector with the degree of
! each node of the graph.
dg(0:np-1) = l_dep_list(0:np-1)
!
! 2. Build a list of all edges, e.g. copy the
! dependencies lists keeping only (I,J) with I<J
!
nedges = 0
do i = 0, np-1
do j = dl_ptr(i),dl_ptr(i+1) - 1
ip = c_dep_list(j)
if (i<=ip) then
nedges = nedges + 1
edges(1,nedges) = i
edges(2,nedges) = ip
end if
end do
end do
!
! 3. Loop over all edges
!
ist = 1
do while (ist <= nedges)
!
! 4. Weight the edges with the sum of the degrees
! of their nodes and sort them into descending order
upd(:) = 0
do i = ist, nedges
dgp(i) = (dg(edges(1,i)) + dg(edges(2,i)))
end do
call psb_msort(dgp(ist:nedges),ix=idx(ist:nedges),dir=psb_sort_down_)
! 5. Scan the list of edges; if neither node of the
! edge has been marked yet, take out the edge and mark
! the two nodes
i1 = ist
nch = 0
do i = ist, nedges
ix = idx(i)+ist-1
ip1 = edges(1,ix)
ip2 = edges(2,ix)
if ((upd(ip1)==0).and.(upd(ip2)==0)) then
upd(ip1) = -1
upd(ip2) = -1
nch = nch + 1
ich(nch) = ix
end if
end do
!
! 6. If no edge was chosen but the graph is nonempty
! raise an error condition
if (nch == 0) then
write(psb_err_unit,*)&
& 'srtlist ?????? impossible error !!!!!?????',&
& nedges,ist
do i=ist, nedges
ix = idx(i)+ist-1
write(psb_err_unit,*)&
& 'SRTLIST: Edge:',ix,edges(1,ix),&
& edges(2,ix),dgp(ix)
end do
info = psb_err_input_value_invalid_i_
return
end if
!
! 7. Queue the edges in the matching to the output
! sequence; decrease by 1 the degree of all marked
! nodes, then clear all marks
!
call psb_msort(ich(1:nch))
do i=1, nch
iswap(1:2) = edges(1:2,ist)
edges(1:2,ist) = edges(1:2,ich(i))
edges(1:2,ich(i)) = iswap(1:2)
ist = ist + 1
end do
do i=0, np-1
dg(i) = dg(i) + upd(i)
end do
end do
internal_error = .false.
do i=0, np-1
if (dg(i) /= 0) then
internal_error = .true.
if (me == 0) write(psb_err_unit,*)&
& 'csr_SRTLIST Error on exit:',i,dg(i)
end if
dg(i) = 0
end do
if (internal_error .and. (me==0)) then
write(0,*) 'Error on srt_list. Input:'
do i = 0, np-1
write(0,*) 'Proc: ',i,' list: '
write(0,*) c_dep_list(dl_ptr(i):dl_ptr(i+1) - 1)
end do
end if
!
! 10. Scan the edge sequence;
! for each edge, take each one of its
! endpoints and queue the other
! node in the endpoint dependency list
!
do j=1,nedges
i = edges(1,j)
ix = dl_ptr(i)
c_dep_list(ix+dg(i)) = edges(2,j)
dg(i) = dg(i)+1
i = edges(2,j)
ix = dl_ptr(i)
c_dep_list(ix+dg(i)) = edges(1,j)
dg(i) = dg(i)+1
!
! If there are any self loops, adjust for error condition
! check
!
if (edges(1,j) == edges(2,j)) dg(i) = dg(i) -1
end do
do i=0, np-1
if (dg(i) /= l_dep_list(i)) then
if (me == 0) write(psb_err_unit,*) &
& 'SRTLIST Mismatch on output',i,dg(i),l_dep_list(i)
end if
end do
9999 call psb_error_handler(err_act)
return
end subroutine psi_i_sort_dl
end subroutine psi_i_csr_sort_dl
base/internals/psi_srtlist.f90
+2 -2
View File
@@ -114,7 +114,7 @@ subroutine srtlist(dep_list,dl_lda,ldl,np,dg,dgp,upd, edges,idx,ich,info)
do i=1, np
do j=1, dg(i)
ip = dep_list(j,i) + 1
if (ip.gt.i) then
if (ip >= i) then
iedge = iedge + 1
edges(1,iedge) = i
edges(2,iedge) = ip
@@ -124,7 +124,6 @@ subroutine srtlist(dep_list,dl_lda,ldl,np,dg,dgp,upd, edges,idx,ich,info)
ist = 1
do while (ist.le.nedges)
do i=1, np
upd(i) = 0
enddo
@@ -188,6 +187,7 @@ subroutine srtlist(dep_list,dl_lda,ldl,np,dg,dgp,upd, edges,idx,ich,info)
i = edges(2,j)
dg(i) = dg(i)+1
dep_list(dg(i),i) = edges(1,j)-1
if (edges(1,j) == edges(2,j)) dg(i) = dg(i) -1
enddo
do i=1, np
if (dg(i).ne.ldl(i)) then
base/internals/psi_xtr_loc_dl.F90
+230
View File
@@ -0,0 +1,230 @@
!
! 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.
!
!
subroutine psi_i_xtr_loc_dl(ictxt,is_bld,is_upd,desc_str,loc_dl,length_dl,info)
! internal routine
! == = =============
!
! _____called by psi_crea_halo and psi_crea_ovrlap ______
!
! purpose
! == = ====
! process root (pid=0) extracts for each process "k" the ordered list of process
! to which "k" must communicate. this list with its order is extracted from
! desc_str list
!
!
! input
! == = ====
! desc_data :integer array
! explanation:
! name explanation
! ------------------ -------------------------------------------------------
! desc_data array of integer that contains some local and global
! information of matrix.
!
!
! now we explain each of the above vectors.
!
! let a be a generic sparse matrix. we denote with matdata_a the matrix_data
! array for matrix a.
! data stored in matrix_data array are:
!
! notation stored in explanation
! --------------- ---------------------- -------------------------------------
! dec_type matdata_a[psb_dec_type_] decomposition type
! m matdata_a[m_] total number of equations
! n matdata_a[n_] total number of variables
! n_row matdata_a[psb_n_row_] number of local equations
! n_col matdata_a[psb_n_col_] number of local variables
! psb_ctxt_a matdata_a[ctxt_] the blacs context handle, indicating
! the global context of the operation
! on the matrix.
! the context itself is global.
! desc_str integer array
! explanation:
! let desc_str_p be the array desc_str for local process.
!! this is composed of variable dimension blocks for each process to
! communicate to.
! each block contain indexes of local halo elements to exchange with other
! process.
! let p be the pointer to the first element of a block in desc_str_p.
! this block is stored in desc_str_p as :
!
! notation stored in explanation
! --------------- --------------------------- -----------------------------------
! process_id desc_str_p[p+psb_proc_id_] identifier of process which exchange
! data with.
! n_elements_recv desc_str_p[p+n_elem_recv_] number of elements to receive.
! elements_recv desc_str_p[p+elem_recv_+i] indexes of local elements to
! receive. these are stored in the
! array from location p+elem_recv_ to
! location p+elem_recv_+
! desc_str_p[p+n_elem_recv_]-1.
! if desc_data(psb_dec_type_) == 0
! then also will be:
! n_elements_send desc_str_p[p+n_elem_send_] number of elements to send.
! elements_send desc_str_p[p+elem_send_+i] indexes of local elements to
! send. these are stored in the
! array from location p+elem_send_ to
! location p+elem_send_+
! desc_str_p[p+n_elem_send_]-1.
! list is ended by -1 value
!
! np integer (global input)
! number of grid process.
!
! output
! == = ==
! loc_dl integer array(:)
! dependence list of current process
!
use psi_mod, psb_protect_name => psi_i_xtr_loc_dl
#ifdef MPI_MOD
use mpi
#endif
use psb_penv_mod
use psb_const_mod
use psb_error_mod
use psb_desc_mod
use psb_sort_mod
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
! ....scalar parameters...
logical, intent(in) :: is_bld, is_upd
integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(in) :: desc_str(:)
integer(psb_ipk_), allocatable, intent(out) :: loc_dl(:), length_dl(:)
integer(psb_ipk_), intent(out) :: info
! .....local arrays....
integer(psb_ipk_) :: int_err(5)
! .....local scalars...
integer(psb_ipk_) :: i,pdl,proc,j,err_act, ldu
integer(psb_ipk_) :: err
integer(psb_ipk_) :: debug_level, debug_unit
integer(psb_mpk_) :: iictxt, icomm, me, np, minfo
logical, parameter :: dist_symm_list=.true., print_dl=.false., profile=.true.
character name*20
name='psi_extrct_dl'
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
iictxt = ictxt
info = psb_success_
call psb_info(iictxt,me,np)
pdl = size(desc_str)
allocate(loc_dl(pdl+1),length_dl(0:np),stat=info)
if (info /= psb_success_) then
info=psb_err_alloc_dealloc_
goto 9999
end if
if (debug_level >= psb_debug_inner_)&
& write(debug_unit,*) me,' ',trim(name),': start ',info
loc_dl = -1
i = 1
pdl = 0
if (is_bld) then
do while (desc_str(i) /= -1)
if (debug_level >= psb_debug_inner_)&
& write(debug_unit,*) me,' ',trim(name),' : looping ',i,&
& desc_str(i),desc_str(i+1),desc_str(i+2)
! ...with different decomposition type we have different
! structure of indices lists............................
if ((desc_str(i+1) /= 0).or.(desc_str(i+2) /= 0)) then
! ..if number of element to be exchanged !=0
proc=desc_str(i)
if ((proc < 0).or.(proc >= np)) then
if (debug_level >= psb_debug_inner_)&
& write(debug_unit,*) me,' ',trim(name),': error ',i,desc_str(i)
info = 9999
int_err(1) = i
int_err(2) = desc_str(i)
goto 9999
endif
! if((me == 1).and.(proc == 3))write(psb_err_unit,*)'found 3'
pdl=pdl+1
loc_dl(pdl)=proc
endif
i=i+desc_str(i+1)+2
enddo
else if (is_upd) then
do while (desc_str(i) /= -1)
if (debug_level >= psb_debug_inner_) &
& write(debug_unit,*) me,' ',trim(name),': looping ',i,desc_str(i)
! ...with different decomposition type we have different
! structure of indices lists............................
if (desc_str(i+1) /= 0) then
proc=desc_str(i)
! ..if number of element to be exchanged !=0
pdl=pdl+1
loc_dl(pdl)=proc
endif
i=i+desc_str(i+1)+2
enddo
else
info = 2020
goto 9999
endif
call psb_msort_unique(loc_dl(1:pdl),ldu)
pdl = ldu
call psb_realloc(pdl,loc_dl,info)
call psi_symm_dep_list(loc_dl,ictxt,info)
pdl = size(loc_dl)
length_dl = 0
length_dl(me) = pdl
call psb_sum(ictxt,length_dl)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_errpush(info,name,i_err=int_err)
call psb_error_handler(err_act)
return
end subroutine psi_i_xtr_loc_dl
base/modules/desc/psb_gen_block_map_mod.F90
+35 -27
View File
@@ -2192,16 +2192,7 @@ contains
integer(psb_ipk_) :: lb, ub, m
if (n < 5) then
! don't bother with binary search for very
! small vectors
ipos = 0
do
if (ipos == n) return
if (key < v(ipos+1)) return
ipos = ipos + 1
end do
else
choice: if (n >5) then
lb = 1
ub = n
ipos = -1
@@ -2210,7 +2201,7 @@ contains
m = (lb+ub)/2
if (key==v(m)) then
ipos = m
return
exit choice
else if (key < v(m)) then
ub = m-1
else
@@ -2220,8 +2211,21 @@ contains
if (v(ub) > key) then
ub = ub - 1
end if
ipos = ub
endif
ipos = ub
else
! No binary search, do everything in the final cleanup
ipos = 0
end if choice
! Final cleanup
! This is needed because V may contain repeated entries
! i.e. there may be processes that own 0 indices
do
if (ipos == n) exit
if (key < v(ipos+1) ) exit
ipos = ipos + 1
end do
return
end function i_gen_block_search
@@ -2234,17 +2238,8 @@ contains
integer(psb_lpk_) :: v(:)
integer(psb_ipk_) :: lb, ub, m
if (n < 5) then
! don't bother with binary search for very
! small vectors
ipos = 0
do
if (ipos == n) return
if (key < v(ipos+1)) return
ipos = ipos + 1
end do
else
choice: if (n >5) then
lb = 1
ub = n
ipos = -1
@@ -2253,7 +2248,7 @@ contains
m = (lb+ub)/2
if (key==v(m)) then
ipos = m
return
exit choice
else if (key < v(m)) then
ub = m-1
else
@@ -2263,8 +2258,21 @@ contains
if (v(ub) > key) then
ub = ub - 1
end if
ipos = ub
endif
ipos = ub
else
! No binary search, do everything in the final cleanup
ipos = 0
end if choice
! Final cleanup
! This is needed because V may contain repeated entries
! i.e. there may be processes that own 0 indices
do
if (ipos == n) exit
if (key < v(ipos+1) ) exit
ipos = ipos + 1
end do
return
end function l_gen_block_search
#endif
base/modules/psi_i_mod.F90
+34 -13
View File
@@ -94,23 +94,15 @@ module psi_i_mod
end subroutine psi_i_desc_index
end interface
interface psi_dl_check
subroutine psi_i_dl_check(dep_list,dl_lda,np,length_dl)
import
implicit none
integer(psb_ipk_) :: np,dl_lda,length_dl(0:np)
integer(psb_ipk_) :: dep_list(dl_lda,0:np)
end subroutine psi_i_dl_check
end interface
interface psi_sort_dl
subroutine psi_i_sort_dl(dep_list,l_dep_list,np,info)
subroutine psi_i_csr_sort_dl(dl_ptr,c_dep_list,l_dep_list,ictxt,info)
import
implicit none
integer(psb_ipk_) :: dep_list(:,:), l_dep_list(:)
integer(psb_ipk_) :: np
integer(psb_ipk_), intent(in) :: c_dep_list(:), dl_ptr(0:)
integer(psb_ipk_), intent(inout) :: l_dep_list(0:)
integer(psb_ipk_) :: ictxt
integer(psb_ipk_) :: info
end subroutine psi_i_sort_dl
end subroutine psi_i_csr_sort_dl
end interface
interface psi_extract_dep_list
@@ -127,6 +119,35 @@ module psi_i_mod
end subroutine psi_i_extract_dep_list
end interface
interface psi_bld_glb_dep_list
subroutine psi_i_bld_glb_dep_list(ictxt,loc_dl,length_dl,dep_list,dl_lda,info)
import
integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(out) :: dl_lda
integer(psb_ipk_), intent(in) :: loc_dl(:), length_dl(0:)
integer(psb_ipk_), allocatable, intent(out) :: dep_list(:,:)
integer(psb_ipk_), intent(out) :: info
end subroutine psi_i_bld_glb_dep_list
subroutine psi_i_bld_glb_csr_dep_list(ictxt,loc_dl,length_dl,c_dep_list,dl_ptr,info)
import
integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(in) :: loc_dl(:), length_dl(0:)
integer(psb_ipk_), allocatable, intent(out) :: c_dep_list(:), dl_ptr(:)
integer(psb_ipk_), intent(out) :: info
end subroutine psi_i_bld_glb_csr_dep_list
end interface
interface psi_extract_loc_dl
subroutine psi_i_xtr_loc_dl(ictxt,is_bld,is_upd,desc_str,loc_dl,length_dl,info)
import
logical, intent(in) :: is_bld, is_upd
integer(psb_ipk_), intent(in) :: ictxt
integer(psb_ipk_), intent(in) :: desc_str(:)
integer(psb_ipk_), allocatable, intent(out) :: loc_dl(:), length_dl(:)
integer(psb_ipk_), intent(out) :: info
end subroutine psi_i_xtr_loc_dl
end interface
interface psi_fnd_owner
subroutine psi_i_fnd_owner(nv,idx,iprc,desc,info)
import
base/modules/serial/psb_base_mat_mod.F90
+32
View File
@@ -588,8 +588,16 @@ module psb_base_mat_mod
procedure, pass(a) :: reinit => psb_lbase_reinit
procedure, pass(a) :: allocate_mnnz => psb_lbase_allocate_mnnz
procedure, pass(a) :: reallocate_nz => psb_lbase_reallocate_nz
#if defined(IPK4) && defined(LPK8)
procedure, pass(a) :: allocate_imnnz => psb_lbase_allocate_imnnz
procedure, pass(a) :: reallocate_inz => psb_lbase_reallocate_inz
generic, public :: allocate => allocate_mnnz, allocate_imnnz
generic, public :: reallocate => reallocate_nz, reallocate_inz
#else
generic, public :: allocate => allocate_mnnz
generic, public :: reallocate => reallocate_nz
#endif
procedure, pass(a) :: csgetptn => psb_lbase_csgetptn
@@ -1413,6 +1421,30 @@ contains
end subroutine psb_lbase_set_lncols
#if defined(IPK4) && defined(LPK8)
subroutine psb_lbase_allocate_imnnz(m,n,a,nz)
implicit none
integer(psb_ipk_), intent(in) :: m,n
class(psb_lbase_sparse_mat), intent(inout) :: a
integer(psb_ipk_), intent(in), optional :: nz
integer(psb_lpk_) :: lm, ln, lnz
lm = m
ln = n
if (present(nz)) then
lnz = nz
call a%allocate(lm,ln,lnz)
else
call a%allocate(lm,ln)
end if
end subroutine psb_lbase_allocate_imnnz
subroutine psb_lbase_reallocate_inz(nz,a)
integer(psb_ipk_), intent(in) :: nz
class(psb_lbase_sparse_mat), intent(inout) :: a
integer(psb_lpk_) :: lnz
lnz = nz
call a%reallocate(lnz)
end subroutine psb_lbase_reallocate_inz
subroutine psb_lbase_set_inrows(m,a)
implicit none
class(psb_lbase_sparse_mat), intent(inout) :: a
base/tools/psb_icdasb.F90
+1 -1
View File
@@ -65,7 +65,7 @@ subroutine psb_icdasb(desc,info,ext_hv,mold)
integer(psb_ipk_) :: i, n_col, dectype, err_act, n_row
integer(psb_mpk_) :: np,me, icomm, ictxt
logical :: ext_hv_
logical, parameter :: do_timings=.false.
logical, parameter :: do_timings=.true.
integer(psb_ipk_), save :: idx_phase1=-1, idx_phase2=-1, idx_phase3=-1
integer(psb_ipk_), save :: idx_phase11=-1, idx_phase12=-1, idx_phase13=-1
integer(psb_ipk_), save :: idx_total=-1
test/pargen/runs/ppde.inp
+1 -1
View File
@@ -6,7 +6,7 @@ CSR Storage format for matrix A: CSR COO
3 Partition: 1 BLOCK 3 3D
2 Stopping criterion 1 2
0100 MAXIT
-1 ITRACE
05 ITRACE
002 IRST restart for RGMRES and BiCGSTABL