Second step of major reorg: tested compilation.

psblas3-type-indexed
Salvatore Filippone 20 years ago
parent 0383e80618
commit af5b65606c

@ -1,6 +1,15 @@
Changelog. A lot less detailed than usual, at least for past
history.
2007/01/11: Migrated repository to SVN.
2007/01/11: MLD2P4 has been moved to the new org. Now tackling the
test dirs.
2007/01/09: First try at reorganizing directories. Subdir MLD2P4 still
to be fixed. Documentation still to be updated.
2006/12/11: Documented options in glob_to_loc.
2006/12/06: Fixed raw aggregation.

@ -40,6 +40,12 @@ SLUDEF=-DHave_SLU_ -I$(SLUDIR)
UMFDIR=$(HOME)/LIB/Umfpack_gcc41
UMF=-lumfpack -lamd -L$(UMFDIR)
UMFDEF=-DHave_UMF_ -I$(UMFDIR)
#
# We are using the public domain tool METIS from U. Minnesota. To get it
# check URL http://www.cs.umn.edu:~karypis
#
METIS_LIB = -L$(HOME)/NUMERICAL/metis-4.0 -lmetis
LDLIBS=$(BLACS) $(SLU) $(UMF) $(BLAS) $(METIS_LIB)
# Add -DLargeFptr for 64-bit addresses
CDEFINES=-DAdd_ $(SLUDEF) $(UMFDEF)
@ -53,21 +59,11 @@ RANLIB=ranlib
####################### Section 5 #######################
# Do not edit this #
##########################################################
LIBDIR = lib
LIBNAME = libpsblas.a
TYPEMODS = psb_spmat_type$(.mod) psb_descriptor_type$(.mod) psb_prec_type$(.mod) psb_realloc_mod$(.mod)
CONSTMODS = psb_const_mod$(.mod)
BLASMODS = $(TYPEMODS) psb_psblas_mod$(.mod) psb_comm_mod$(.mod)
METHDMODS = psb_methd_mod$(.mod)
TOOLSMODS = $(CONSTMODS) psi_mod$(.mod) psb_tools_mod$(.mod) psb_serial_mod$(.mod)
PRECMODS = psb_prec_mod$(.mod)
ERRORMODS = psb_error_mod$(.mod)
F90MODS= $(BLASMODS) $(PRECMODS) $(METHDMODS) $(TOOLSMODS) $(ERRORMODS) string$(.mod)
MODS=$(LIBDIR)/psb_const_mod$(.mod) $(LIBDIR)/psb_spmat_type$(.mod) $(LIBDIR)/psb_realloc_mod$(.mod) \
$(LIBDIR)/psb_descriptor_type$(.mod) $(LIBDIR)/psb_prec_type$(.mod) $(LIBDIR)/parts.fh \
$(LIBDIR)/psb_serial_mod$(.mod) $(LIBDIR)/psb_comm_mod$(.mod) $(LIBDIR)/psb_error_mod$(.mod)
LIBDIR=lib
BASELIBNAME=libpsb_base.a
PRECLIBNAME=libpsb_prec.a
METHDLIBNAME=libpsb_krylov.a
UTILLIBNAME=libpsb_util.a
# Under Linux/gmake there is a rule interpreting .mod as Modula source!
$(.mod).o:

@ -10,9 +10,9 @@ F90=/usr/local/gcc42/bin/gfortran
FC=/usr/local/gcc42/bin/gfortran
F77=$(FC)
CC=/usr/local/gcc42/bin/gcc
F90COPT= -O3 -march=pentium4 -msse2 -mfpmath=sse
FCOPT=-O3 -march=pentium4 -msse2 -mfpmath=sse
CCOPT=-O3 -march=pentium4 -msse2 -mfpmath=sse
F90COPT=-O3 -ggdb
FCOPT=-O3 -ggdb
CCOPT=-O3 -ggdb
####################### Section 2 #######################
# Define your linker and linker flags here #
@ -40,6 +40,12 @@ SLUDEF=-DHave_SLU_ -I$(SLUDIR)
UMFDIR=$(HOME)/LIB/Umfpack_gcc41
UMF=-lumfpack -lamd -L$(UMFDIR)
UMFDEF=-DHave_UMF_ -I$(UMFDIR)
#
# We are using the public domain tool METIS from U. Minnesota. To get it
# check URL http://www.cs.umn.edu:~karypis
#
METIS_LIB = -L$(HOME)/NUMERICAL/metis-4.0 -lmetis
LDLIBS=$(BLACS) $(SLU) $(UMF) $(BLAS) $(METIS_LIB)
# Add -DLargeFptr for 64-bit addresses
CDEFINES=-DAdd_ $(SLUDEF) $(UMFDEF)
@ -53,21 +59,11 @@ RANLIB=ranlib
####################### Section 5 #######################
# Do not edit this #
##########################################################
LIBDIR = lib
LIBNAME = libpsblas.a
TYPEMODS = psb_spmat_type$(.mod) psb_descriptor_type$(.mod) psb_prec_type$(.mod) psb_realloc_mod$(.mod)
CONSTMODS = psb_const_mod$(.mod)
BLASMODS = $(TYPEMODS) psb_psblas_mod$(.mod) psb_comm_mod$(.mod)
METHDMODS = psb_methd_mod$(.mod)
TOOLSMODS = $(CONSTMODS) psi_mod$(.mod) psb_tools_mod$(.mod) psb_serial_mod$(.mod)
PRECMODS = psb_prec_mod$(.mod)
ERRORMODS = psb_error_mod$(.mod)
F90MODS= $(BLASMODS) $(PRECMODS) $(METHDMODS) $(TOOLSMODS) $(ERRORMODS) string$(.mod)
MODS=$(LIBDIR)/psb_const_mod$(.mod) $(LIBDIR)/psb_spmat_type$(.mod) $(LIBDIR)/psb_realloc_mod$(.mod) \
$(LIBDIR)/psb_descriptor_type$(.mod) $(LIBDIR)/psb_prec_type$(.mod) $(LIBDIR)/parts.fh \
$(LIBDIR)/psb_serial_mod$(.mod) $(LIBDIR)/psb_comm_mod$(.mod) $(LIBDIR)/psb_error_mod$(.mod)
LIBDIR=lib
BASELIBNAME=libpsb_base.a
PRECLIBNAME=libpsb_prec.a
METHDLIBNAME=libpsb_krylov.a
UTILLIBNAME=libpsb_util.a
# Under Linux/gmake there is a rule interpreting .mod as Modula source!
$(.mod).o:

@ -18,11 +18,11 @@ CCOPT=-O3
####################### Section 2 #######################
# Define your linker and linker flags here #
##########################################################
F90LINK=/usr/local/mpich-ifc91/bin/mpif90 -g -CB -no_cpprt
FLINK=/usr/local/mpich-ifc91/bin/mpif77 -g -CB -no_cpprt
MPF90=/usr/local/mpich-ifc91/bin/mpif90 -g -CB -no_cpprt
MPF77=/usr/local/mpich-ifc91/bin/mpif77 -g -CB -no_cpprt
MPCC=/usr/local/mpich-ifc91/bin/mpicc -g -CB -no_cpprt
F90LINK=/usr/local/mpich-ifc91/bin/mpif90
FLINK=/usr/local/mpich-ifc91/bin/mpif77
MPF90=/usr/local/mpich-ifc91/bin/mpif90
MPF77=/usr/local/mpich-ifc91/bin/mpif77
MPCC=/usr/local/mpich-ifc91/bin/mpicc
####################### Section 3 #######################
# Specify paths to libraries #
@ -42,6 +42,12 @@ SLUDEF=-DHave_SLU_ -I$(SLUDIR)
UMFDIR=$(HOME)/LIB/Umfpack_gcc41
UMF=-lumfpack -lamd -L$(UMFDIR)
UMFDEF=-DHave_UMF_ -I$(UMFDIR)
#
# We are using the public domain tool METIS from U. Minnesota. To get it
# check URL http://www.cs.umn.edu:~karypis
#
METIS_LIB = -L$(HOME)/NUMERICAL/metis-4.0 -lmetis
LDLIBS=$(BLACS) $(SLU) $(UMF) $(BLAS) $(METIS_LIB)
# Add -DLargeFptr for 64-bit addresses
CDEFINES=-DAdd_ $(SLUDEF) $(UMFDEF)
@ -55,21 +61,11 @@ RANLIB=ranlib
####################### Section 5 #######################
# Do not edit this #
##########################################################
LIBDIR = lib
LIBNAME = libpsblas.a
TYPEMODS = psb_spmat_type$(.mod) psb_descriptor_type$(.mod) psb_prec_type$(.mod) psb_realloc_mod$(.mod)
CONSTMODS = psb_const_mod$(.mod)
BLASMODS = $(TYPEMODS) psb_psblas_mod$(.mod) psb_comm_mod$(.mod)
METHDMODS = psb_methd_mod$(.mod)
TOOLSMODS = $(CONSTMODS) psi_mod$(.mod) psb_tools_mod$(.mod) psb_serial_mod$(.mod)
PRECMODS = psb_prec_mod$(.mod)
ERRORMODS = psb_error_mod$(.mod)
F90MODS= $(BLASMODS) $(PRECMODS) $(METHDMODS) $(TOOLSMODS) $(ERRORMODS) string$(.mod)
MODS=$(LIBDIR)/psb_const_mod$(.mod) $(LIBDIR)/psb_spmat_type$(.mod) $(LIBDIR)/psb_realloc_mod$(.mod) \
$(LIBDIR)/psb_descriptor_type$(.mod) $(LIBDIR)/psb_prec_type$(.mod) $(LIBDIR)/parts.fh \
$(LIBDIR)/psb_serial_mod$(.mod) $(LIBDIR)/psb_comm_mod$(.mod) $(LIBDIR)/psb_error_mod$(.mod)
LIBDIR=lib
BASELIBNAME=libpsb_base.a
PRECLIBNAME=libpsb_prec.a
METHDLIBNAME=libpsb_krylov.a
UTILLIBNAME=libpsb_util.a
# Under Linux/gmake there is a rule interpreting .mod as Modula source!
$(.mod).o:

@ -40,6 +40,12 @@ BLACS=-lmpiblacs
#UMFDIR=$(HOME)/LIB/Umfpack
#UMF=-lumfpack -lamd -L$(UMFDIR)
#UMFDEF=-DHave_UMF_ -I$(UMFDIR)
#
# We are using the public domain tool METIS from U. Minnesota. To get it
# check URL http://www.cs.umn.edu:~karypis
#
METIS_LIB = -L$(HOME)/NUMERICAL/metis-4.0 -lmetis
LDLIBS=$(BLACS) $(SLU) $(UMF) $(BLAS) $(METIS_LIB)
# Add -DLargeFptr for 64-bit addresses
CDEFINES=-DNoChange $(SLUDEF) $(UMFDEF) -DLargeFptr
@ -53,21 +59,11 @@ RANLIB=ranlib
####################### Section 5 #######################
# Do not edit this #
##########################################################
LIBDIR = lib
LIBNAME = libpsblas.a
TYPEMODS = psb_spmat_type$(.mod) psb_descriptor_type$(.mod) psb_prec_type$(.mod) psb_realloc_mod$(.mod)
CONSTMODS = psb_const_mod$(.mod)
BLASMODS = $(TYPEMODS) psb_psblas_mod$(.mod) psb_comm_mod$(.mod)
METHDMODS = psb_methd_mod$(.mod)
TOOLSMODS = $(CONSTMODS) psi_mod$(.mod) psb_tools_mod$(.mod) psb_serial_mod$(.mod)
PRECMODS = psb_prec_mod$(.mod)
ERRORMODS = psb_error_mod$(.mod)
F90MODS= $(BLASMODS) $(PRECMODS) $(METHDMODS) $(TOOLSMODS) $(ERRORMODS) string$(.mod)
MODS=$(LIBDIR)/psb_const_mod$(.mod) $(LIBDIR)/psb_spmat_type$(.mod) $(LIBDIR)/psb_realloc_mod$(.mod) \
$(LIBDIR)/psb_descriptor_type$(.mod) $(LIBDIR)/psb_prec_type$(.mod) $(LIBDIR)/parts.fh \
$(LIBDIR)/psb_serial_mod$(.mod) $(LIBDIR)/psb_comm_mod$(.mod) $(LIBDIR)/psb_error_mod$(.mod)
LIBDIR=lib
BASELIBNAME=libpsb_base.a
PRECLIBNAME=libpsb_prec.a
METHDLIBNAME=libpsb_krylov.a
UTILLIBNAME=libpsb_util.a
# Under Linux/gmake there is a rule interpreting .mod as Modula source!
$(.mod).o:

@ -1,18 +1,29 @@
include Make.inc
#PREC=mld2p4
PREC=baseprec
library:
( [ -d lib ] || mkdir lib)
(cd src; make lib)
(cd base; make lib)
(cd $(PREC); make lib )
(cd krylov; make lib)
(cd util; make lib )
@echo "====================================="
@echo "Compilation Successful."
@echo "You can now link to ./lib/libpsblas.a"
clean:
(cd src; make clean)
(cd base; make clean)
(cd $(PREC); make clean )
(cd krylov; make clean)
(cd util; make clean)
cleanlib:
(cd lib; /bin/rm -f *.a *$(.mod) *$(.fh))
veryclean: cleanlib
(cd src; make veryclean)
(cd base; make veryclean)
(cd $(PREC); make veryclean )
(cd krylov; make veryclean)
(cd util; make veryclean)
.PHONY: lib

@ -0,0 +1,27 @@
include ../Make.inc
HERE=.
LIBDIR=../lib
LIBNAME=$(BASELIBNAME)
LIBMOD=psb_base_mod$(.mod)
lib:
(cd modules; make lib LIBNAME=$(BASELIBNAME))
(cd comm; make lib LIBNAME=$(BASELIBNAME))
(cd internals; make lib LIBNAME=$(BASELIBNAME))
(cd tools; make lib LIBNAME=$(BASELIBNAME))
(cd serial; make lib LIBNAME=$(BASELIBNAME))
(cd psblas; make lib LIBNAME=$(BASELIBNAME))
/bin/cp $(HERE)/$(LIBNAME) $(LIBDIR)
/bin/cp $(LIBMOD) $(LIBDIR)
clean:
(cd modules; make clean)
(cd comm; make clean)
(cd internals; make clean)
(cd tools; make clean)
(cd serial; make clean)
(cd psblas; make clean)
veryclean: clean
/bin/rm -f $(HERE)/$(LIBNAME) $(LIBMOD)

@ -4,8 +4,9 @@ OBJS = psb_dgather.o psb_dhalo.o psb_dovrl.o \
psb_ihalo.o psb_zgather.o psb_zhalo.o psb_zovrl.o
MPFOBJS = psb_dscatter.o psb_zscatter.o
INCDIRS = -I ../../lib -I .
LIBDIR = ../../lib
LIBDIR = ..
MODDIR = ../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I .
lib: mpfobjs $(OBJS)
$(AR) $(LIBDIR)/$(LIBNAME) $(MPFOBJS) $(OBJS)

@ -13,8 +13,9 @@ MPFOBJS = psi_dswapdata.o psi_dswaptran.o psi_iswapdata.o \
psi_iswaptran.o psi_desc_index.o \
psi_zswapdata.o psi_zswaptran.o
MPFOBJS2 = psi_extrct_dl.o
INCDIRS = -I ../../lib -I .
LIBDIR = ../../lib
LIBDIR = ..
MODDIR = ../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I .
lib: mpfobjs $(FOBJS) $(FOBJS2) $(COBJS) $(MPFOBJS2)
$(AR) $(LIBDIR)/$(LIBNAME) $(MPFOBJS) $(MPFOBJS2) $(FOBJS) $(FOBJS2) \

@ -34,7 +34,6 @@ subroutine psi_crea_index(desc_a,index_in,index_out,glob_idx,nxch,nsnd,nrcv,info
use psb_descriptor_type
use psb_error_mod
use psb_penv_mod
use psi_mod, only: psi_sort_dl, psi_desc_index, psi_dl_check
implicit none
type(psb_desc_type), intent(in) :: desc_a
@ -51,6 +50,33 @@ subroutine psi_crea_index(desc_a,index_in,index_out,glob_idx,nxch,nsnd,nrcv,info
logical,parameter :: debug=.false.
character(len=20) :: name
interface
subroutine psi_sort_dl(dep_list,l_dep_list,np,info)
integer :: np,dep_list(:,:), l_dep_list(:), info
end subroutine psi_sort_dl
end interface
interface
subroutine psi_dl_check(dep_list,dl_lda,np,length_dl)
integer :: np,dl_lda,length_dl(0:np)
integer :: dep_list(dl_lda,0:np)
end subroutine psi_dl_check
end interface
interface
subroutine psi_desc_index(desc,index_in,dep_list,&
& length_dl,nsnd,nrcv,desc_index,&
& isglob_in,info)
use psb_descriptor_type
type(psb_desc_type) :: desc
integer :: index_in(:),dep_list(:)
integer, allocatable :: desc_index(:)
integer :: length_dl,nsnd,nrcv,info
logical :: isglob_in
end subroutine psi_desc_index
end interface
info = 0
name='psi_crea_index'
call psb_erractionsave(err_act)

@ -122,7 +122,7 @@ contains
implicit none
integer :: n, idx(:)
integer :: n, k, idx(:)
real(kind(1.d0)) :: beta, x(:), y(:)
! Locals

@ -3,18 +3,19 @@ include ../../Make.inc
MODULES = psb_realloc_mod.o psb_string_mod.o psb_spmat_type.o \
psb_desc_type.o psb_spsb_mod.o \
psb_serial_mod.o psb_tools_mod.o \
psb_prec_type.o psb_error_mod.o psb_prec_mod.o \
psb_methd_mod.o psb_const_mod.o \
psb_error_mod.o \
psb_const_mod.o \
psb_comm_mod.o psb_psblas_mod.o psi_mod.o \
psb_check_mod.o blacs_env.o psb_gps_mod.o
# psb_methd_mod.o psb_prec_type.o psb_prec_mod.o
MPFOBJS = psb_penv_mod.o
LIBMOD=psb_base_mod$(.mod)
MPFOBJS=psb_penv_mod.o
OBJS = error.o psb_base_mod.o
LIBDIR = ..
INCDIRS = -I .
OBJS = error.o psb_sparse_mod.o psb_all_mod.o
INCDIRS = -I ../../lib
LIBDIR = ../../lib
psb_realloc_mod.o : psb_error_mod.o
psb_spmat_type.o : psb_realloc_mod.o psb_const_mod.o psb_string_mod.o
@ -26,13 +27,13 @@ psb_check_mod.o: psb_desc_type.o
psb_methd_mod.o: psb_serial_mod.o psb_desc_type.o psb_prec_type.o
psb_tools_mod.o: psb_spmat_type.o psb_desc_type.o psi_mod.o psb_gps_mod.o
psb_gps_mod.o: psb_realloc_mod.o
psb_sparse_mod.o: $(MODULES) $(MPFOBJS)
psb_base_mod.o: $(MODULES) $(MPFOBJS)
lib: mpfobjs $(MODULES) $(OBJS)
$(AR) $(LIBDIR)/$(LIBNAME) $(MODULES) $(OBJS) $(MPFOBJS)
$(RANLIB) $(LIBDIR)/$(LIBNAME)
/bin/cp *$(.mod) ./parts.fh ../../lib
/bin/cp $(LIBMOD) ./parts.fh $(LIBDIR)
mpfobjs:

@ -28,18 +28,15 @@
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
module psb_sparse_mod
module psb_base_mod
use psb_string_mod
use psb_error_mod
use psb_penv_mod
use psb_descriptor_type
use psb_prec_type
use psb_serial_mod
use psb_tools_mod
use psb_psblas_mod
use psb_prec_mod
use psb_methd_mod
use psb_error_mod
use psb_string_mod
end module psb_sparse_mod
use psb_tools_mod
end module psb_base_mod

@ -59,15 +59,13 @@ module psb_const_mod
integer, parameter :: psb_dec_type_=1, psb_m_=2,psb_n_=3
integer, parameter :: psb_n_row_=4, psb_n_col_=5,psb_ctxt_=6
integer, parameter :: psb_loc_to_glob_=7
integer, parameter :: psb_ovl_state_=8
integer, parameter :: psb_mpi_c_=9
integer, parameter :: psb_thal_xch_=11
integer, parameter :: psb_thal_snd_=12
integer, parameter :: psb_thal_rcv_=13
integer, parameter :: psb_tovr_xch_=14
integer, parameter :: psb_tovr_snd_=15
integer, parameter :: psb_tovr_rcv_=16
integer, parameter :: psb_mdata_size_=20
integer, parameter :: psb_mpi_c_=9,psb_mdata_size_=20
integer, parameter :: psb_desc_asb_=3099
integer, parameter :: psb_desc_bld_=psb_desc_asb_+1
integer, parameter :: psb_desc_repl_=3199
@ -75,8 +73,6 @@ module psb_const_mod
integer, parameter :: psb_desc_upd_asb_=psb_desc_upd_+1
integer, parameter :: psb_desc_large_asb_=psb_desc_upd_asb_+1
integer, parameter :: psb_desc_large_bld_=psb_desc_large_asb_+1
integer, parameter :: psb_cd_ovl_bld_=psb_desc_large_bld_+1
integer, parameter :: psb_cd_ovl_asb_=psb_cd_ovl_bld_+1
integer, parameter :: nbits=14
integer, parameter :: hashsize=2**nbits, hashmask=hashsize-1
integer, parameter :: psb_default_large_threshold=4*1024*1024 ! to be reviewed

@ -136,28 +136,6 @@ contains
end function psb_is_asb_desc
logical function psb_is_ovl_bld(desc)
type(psb_desc_type), intent(in) :: desc
psb_is_ovl_bld = (desc%matrix_data(psb_ovl_state_)==psb_cd_ovl_bld_)
end function psb_is_ovl_bld
logical function psb_is_ovl_asb(desc)
type(psb_desc_type), intent(in) :: desc
psb_is_ovl_asb = (desc%matrix_data(psb_ovl_state_)==psb_cd_ovl_asb_)
end function psb_is_ovl_asb
logical function psb_is_ovl_ok(desc)
type(psb_desc_type), intent(in) :: desc
psb_is_ovl_ok = (desc%matrix_data(psb_ovl_state_)==psb_cd_ovl_asb_).or.&
& (desc%matrix_data(psb_ovl_state_)==psb_cd_ovl_bld_)
end function psb_is_ovl_ok
logical function psb_is_ok_dec(dectype)
integer :: dectype
@ -249,69 +227,5 @@ contains
end function psb_is_large_dec
subroutine psb_cd_set_bld(desc,info)
!
! Change state of a descriptor into BUILD.
! If the descriptor is LARGE, check the AVL search tree
! and initialize it if necessary.
!
use psb_const_mod
use psb_error_mod
use psb_penv_mod
implicit none
type(psb_desc_type), intent(inout) :: desc
integer :: info
!locals
integer :: np,me,ictxt, isz, err_act,idx,gidx,lidx
logical, parameter :: debug=.false.,debugprt=.false.
character(len=20) :: name, char_err
if (debug) write(0,*) me,'Entered CDCPY'
if (psb_get_errstatus() /= 0) return
info = 0
call psb_erractionsave(err_act)
name = 'psb_cd_set_bld'
ictxt = psb_cd_get_context(desc)
! check on blacs grid
call psb_info(ictxt, me, np)
if (debug) write(0,*) me,'Entered CDCPY'
if (psb_is_large_desc(desc)) then
if (.not.allocated(desc%ptree)) then
allocate(desc%ptree(2),stat=info)
if (info /= 0) then
info=4000
goto 9999
endif
call InitPairSearchTree(desc%ptree,info)
do idx=1, psb_cd_get_local_cols(desc)
gidx = desc%loc_to_glob(idx)
call SearchInsKeyVal(desc%ptree,gidx,idx,lidx,info)
if (lidx /= idx) then
write(0,*) 'Warning from cdset: mismatch in PTREE ',idx,lidx
endif
enddo
end if
desc%matrix_data(psb_dec_type_) = psb_desc_large_bld_
else
desc%matrix_data(psb_dec_type_) = psb_desc_bld_
end if
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == act_ret) then
return
else
call psb_error(ictxt)
end if
return
end subroutine psb_cd_set_bld
end module psb_descriptor_type

@ -351,8 +351,6 @@ contains
write (0,'("indices in input array are not within problem dimension ",2(i0,2x))')i_e_d(1:2)
case(150)
write (0,'("indices in input array are not belonging to the calling process ",i0)')i_e_d(1)
case(151)
write (0,'("indices in input array are not belonging to the calling process ")')
case(290)
write (0,'("To call this routine you must first call psb_geall on the same matrix")')
case(295)

@ -33,433 +33,394 @@ module psb_serial_mod
use psb_string_mod
interface psb_csdp
subroutine psb_dcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
type(psb_dspmat_type), intent(inout) :: b
integer, intent(out) :: info
integer, intent(in), optional :: ifc,upd,dupl
character, intent(in), optional :: check,trans,unitd
end subroutine psb_dcsdp
subroutine psb_zcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
type(psb_zspmat_type), intent(inout) :: b
integer, intent(out) :: info
integer, intent(in), optional :: ifc,upd,dupl
character, intent(in), optional :: check,trans,unitd
end subroutine psb_zcsdp
subroutine psb_dcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
type(psb_dspmat_type), intent(inout) :: b
integer, intent(out) :: info
integer, intent(in), optional :: ifc,upd,dupl
character, intent(in), optional :: check,trans,unitd
end subroutine psb_dcsdp
subroutine psb_zcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
type(psb_zspmat_type), intent(inout) :: b
integer, intent(out) :: info
integer, intent(in), optional :: ifc,upd,dupl
character, intent(in), optional :: check,trans,unitd
end subroutine psb_zcsdp
end interface
interface psb_csrws
subroutine psb_dcsrws(rw,a,info,trans)
use psb_spmat_type
type(psb_dspmat_type) :: a
real(kind(1.d0)), allocatable :: rw(:)
integer :: info
character, optional :: trans
end subroutine psb_dcsrws
subroutine psb_zcsrws(rw,a,info,trans)
use psb_spmat_type
type(psb_zspmat_type) :: a
complex(kind(1.d0)), allocatable :: rw(:)
integer :: info
character, optional :: trans
end subroutine psb_zcsrws
subroutine psb_dcsrws(rw,a,info,trans)
use psb_spmat_type
type(psb_dspmat_type) :: a
real(kind(1.d0)), allocatable :: rw(:)
integer :: info
character, optional :: trans
end subroutine psb_dcsrws
subroutine psb_zcsrws(rw,a,info,trans)
use psb_spmat_type
type(psb_zspmat_type) :: a
complex(kind(1.d0)), allocatable :: rw(:)
integer :: info
character, optional :: trans
end subroutine psb_zcsrws
end interface
interface psb_cssm
subroutine psb_dcssm(alpha,t,b,beta,c,info,trans,unitd,d)
use psb_spmat_type
type(psb_dspmat_type) :: t
real(kind(1.d0)) :: alpha, beta, b(:,:), c(:,:)
integer :: info
character, optional :: trans, unitd
real(kind(1.d0)), optional, target :: d(:)
end subroutine psb_dcssm
subroutine psb_dcssv(alpha,t,b,beta,c,info,trans,unitd,d)
use psb_spmat_type
type(psb_dspmat_type) :: t
real(kind(1.d0)) :: alpha, beta, b(:), c(:)
integer :: info
character, optional :: trans, unitd
real(kind(1.d0)), optional, target :: d(:)
end subroutine psb_dcssv
subroutine psb_zcssm(alpha,t,b,beta,c,info,trans,unitd,d)
use psb_spmat_type
type(psb_zspmat_type) :: t
complex(kind(1.d0)) :: alpha, beta, b(:,:), c(:,:)
integer :: info
character, optional :: trans, unitd
complex(kind(1.d0)), optional, target :: d(:)
end subroutine psb_zcssm
subroutine psb_zcssv(alpha,t,b,beta,c,info,trans,unitd,d)
use psb_spmat_type
type(psb_zspmat_type) :: t
complex(kind(1.d0)) :: alpha, beta, b(:), c(:)
integer :: info
character, optional :: trans, unitd
complex(kind(1.d0)), optional, target :: d(:)
end subroutine psb_zcssv
subroutine psb_dcssm(alpha,t,b,beta,c,info,trans,unitd,d)
use psb_spmat_type
type(psb_dspmat_type) :: t
real(kind(1.d0)) :: alpha, beta, b(:,:), c(:,:)
integer :: info
character, optional :: trans, unitd
real(kind(1.d0)), optional, target :: d(:)
end subroutine psb_dcssm
subroutine psb_dcssv(alpha,t,b,beta,c,info,trans,unitd,d)
use psb_spmat_type
type(psb_dspmat_type) :: t
real(kind(1.d0)) :: alpha, beta, b(:), c(:)
integer :: info
character, optional :: trans, unitd
real(kind(1.d0)), optional, target :: d(:)
end subroutine psb_dcssv
subroutine psb_zcssm(alpha,t,b,beta,c,info,trans,unitd,d)
use psb_spmat_type
type(psb_zspmat_type) :: t
complex(kind(1.d0)) :: alpha, beta, b(:,:), c(:,:)
integer :: info
character, optional :: trans, unitd
complex(kind(1.d0)), optional, target :: d(:)
end subroutine psb_zcssm
subroutine psb_zcssv(alpha,t,b,beta,c,info,trans,unitd,d)
use psb_spmat_type
type(psb_zspmat_type) :: t
complex(kind(1.d0)) :: alpha, beta, b(:), c(:)
integer :: info
character, optional :: trans, unitd
complex(kind(1.d0)), optional, target :: d(:)
end subroutine psb_zcssv
end interface
interface psb_csmm
subroutine psb_dcsmv(alpha,a,b,beta,c,info,trans)
use psb_spmat_type
type(psb_dspmat_type) :: a
real(kind(1.d0)) :: alpha, beta, b(:), c(:)
integer :: info
character, optional :: trans
end subroutine psb_dcsmv
subroutine psb_dcsmm(alpha,a,b,beta,c,info,trans)
use psb_spmat_type
type(psb_dspmat_type) :: a
real(kind(1.d0)) :: alpha, beta, b(:,:), c(:,:)
integer :: info
character, optional :: trans
end subroutine psb_dcsmm
subroutine psb_zcsmv(alpha,a,b,beta,c,info,trans)
use psb_spmat_type
type(psb_zspmat_type) :: a
complex(kind(1.d0)) :: alpha, beta, b(:), c(:)
integer :: info
character, optional :: trans
end subroutine psb_zcsmv
subroutine psb_zcsmm(alpha,a,b,beta,c,info,trans)
use psb_spmat_type
type(psb_zspmat_type) :: a
complex(kind(1.d0)) :: alpha, beta, b(:,:), c(:,:)
integer :: info
character, optional :: trans
end subroutine psb_zcsmm
subroutine psb_dcsmv(alpha,a,b,beta,c,info,trans)
use psb_spmat_type
type(psb_dspmat_type) :: a
real(kind(1.d0)) :: alpha, beta, b(:), c(:)
integer :: info
character, optional :: trans
end subroutine psb_dcsmv
subroutine psb_dcsmm(alpha,a,b,beta,c,info,trans)
use psb_spmat_type
type(psb_dspmat_type) :: a
real(kind(1.d0)) :: alpha, beta, b(:,:), c(:,:)
integer :: info
character, optional :: trans
end subroutine psb_dcsmm
subroutine psb_zcsmv(alpha,a,b,beta,c,info,trans)
use psb_spmat_type
type(psb_zspmat_type) :: a
complex(kind(1.d0)) :: alpha, beta, b(:), c(:)
integer :: info
character, optional :: trans
end subroutine psb_zcsmv
subroutine psb_zcsmm(alpha,a,b,beta,c,info,trans)
use psb_spmat_type
type(psb_zspmat_type) :: a
complex(kind(1.d0)) :: alpha, beta, b(:,:), c(:,:)
integer :: info
character, optional :: trans
end subroutine psb_zcsmm
end interface
interface psb_fixcoo
subroutine psb_dfixcoo(a,info,idir)
use psb_spmat_type
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
integer, intent(in), optional :: idir
end subroutine psb_dfixcoo
subroutine psb_zfixcoo(a,info,idir)
use psb_spmat_type
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
integer, intent(in), optional :: idir
end subroutine psb_zfixcoo
subroutine psb_dfixcoo(a,info,idir)
use psb_spmat_type
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
integer, intent(in), optional :: idir
end subroutine psb_dfixcoo
subroutine psb_zfixcoo(a,info,idir)
use psb_spmat_type
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
integer, intent(in), optional :: idir
end subroutine psb_zfixcoo
end interface
interface psb_ipcoo2csr
subroutine psb_dipcoo2csr(a,info,rwshr)
use psb_spmat_type
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
logical, optional :: rwshr
end subroutine psb_dipcoo2csr
subroutine psb_zipcoo2csr(a,info,rwshr)
use psb_spmat_type
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
logical, optional :: rwshr
end subroutine psb_zipcoo2csr
subroutine psb_dipcoo2csr(a,info,rwshr)
use psb_spmat_type
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
logical, optional :: rwshr
end subroutine psb_dipcoo2csr
subroutine psb_zipcoo2csr(a,info,rwshr)
use psb_spmat_type
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
logical, optional :: rwshr
end subroutine psb_zipcoo2csr
end interface
interface psb_ipcoo2csc
subroutine psb_dipcoo2csc(a,info,clshr)
use psb_spmat_type
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
logical, optional :: clshr
end subroutine psb_dipcoo2csc
subroutine psb_zipcoo2csc(a,info,clshr)
use psb_spmat_type
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
logical, optional :: clshr
end subroutine psb_zipcoo2csc
subroutine psb_dipcoo2csc(a,info,clshr)
use psb_spmat_type
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
logical, optional :: clshr
end subroutine psb_dipcoo2csc
subroutine psb_zipcoo2csc(a,info,clshr)
use psb_spmat_type
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
logical, optional :: clshr
end subroutine psb_zipcoo2csc
end interface
interface psb_ipcsr2coo
subroutine psb_dipcsr2coo(a,info)
use psb_spmat_type
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
end subroutine psb_dipcsr2coo
subroutine psb_zipcsr2coo(a,info)
use psb_spmat_type
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
end subroutine psb_zipcsr2coo
subroutine psb_dipcsr2coo(a,info)
use psb_spmat_type
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
end subroutine psb_dipcsr2coo
subroutine psb_zipcsr2coo(a,info)
use psb_spmat_type
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
end subroutine psb_zipcsr2coo
end interface
interface psb_csprt
subroutine psb_dcsprt(iout,a,iv,irs,ics,head,ivr,ivc)
use psb_spmat_type
integer, intent(in) :: iout
type(psb_dspmat_type), intent(in) :: a
integer, intent(in), optional :: iv(:)
integer, intent(in), optional :: irs,ics
character(len=*), optional :: head
integer, intent(in), optional :: ivr(:),ivc(:)
end subroutine psb_dcsprt
subroutine psb_zcsprt(iout,a,iv,irs,ics,head,ivr,ivc)
use psb_spmat_type
integer, intent(in) :: iout
type(psb_zspmat_type), intent(in) :: a
integer, intent(in), optional :: iv(:)
integer, intent(in), optional :: irs,ics
character(len=*), optional :: head
integer, intent(in), optional :: ivr(:),ivc(:)
end subroutine psb_zcsprt
subroutine psb_dcsprt(iout,a,iv,irs,ics,head,ivr,ivc)
use psb_spmat_type
integer, intent(in) :: iout
type(psb_dspmat_type), intent(in) :: a
integer, intent(in), optional :: iv(:)
integer, intent(in), optional :: irs,ics
character(len=*), optional :: head
integer, intent(in), optional :: ivr(:),ivc(:)
end subroutine psb_dcsprt
subroutine psb_zcsprt(iout,a,iv,irs,ics,head,ivr,ivc)
use psb_spmat_type
integer, intent(in) :: iout
type(psb_zspmat_type), intent(in) :: a
integer, intent(in), optional :: iv(:)
integer, intent(in), optional :: irs,ics
character(len=*), optional :: head
integer, intent(in), optional :: ivr(:),ivc(:)
end subroutine psb_zcsprt
end interface
interface psb_neigh
subroutine psb_dneigh(a,idx,neigh,n,info,lev)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: idx
integer, intent(out) :: n
integer, allocatable :: neigh(:)
integer, intent(out) :: info
integer, optional, intent(in) :: lev
end subroutine psb_dneigh
subroutine psb_zneigh(a,idx,neigh,n,info,lev)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: idx
integer, intent(out) :: n
integer, allocatable :: neigh(:)
integer, intent(out) :: info
integer, optional, intent(in) :: lev
end subroutine psb_zneigh
subroutine psb_dneigh(a,idx,neigh,n,info,lev)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: idx
integer, intent(out) :: n
integer, allocatable :: neigh(:)
integer, intent(out) :: info
integer, optional, intent(in) :: lev
end subroutine psb_dneigh
subroutine psb_zneigh(a,idx,neigh,n,info,lev)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: idx
integer, intent(out) :: n
integer, allocatable :: neigh(:)
integer, intent(out) :: info
integer, optional, intent(in) :: lev
end subroutine psb_zneigh
end interface
interface psb_coins
subroutine psb_dcoins(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl,rebuild)
use psb_spmat_type
integer, intent(in) :: nz, imin,imax,jmin,jmax
integer, intent(in) :: ia(:),ja(:)
real(kind(1.d0)), intent(in) :: val(:)
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
integer, intent(in), optional :: gtl(:)
logical, optional, intent(in) :: rebuild
end subroutine psb_dcoins
subroutine psb_zcoins(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl,rebuild)
use psb_spmat_type
integer, intent(in) :: nz, imin,imax,jmin,jmax
integer, intent(in) :: ia(:),ja(:)
complex(kind(1.d0)), intent(in) :: val(:)
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
integer, intent(in), optional :: gtl(:)
logical, optional, intent(in) :: rebuild
end subroutine psb_zcoins
subroutine psb_dcoins(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl,rebuild)
use psb_spmat_type
integer, intent(in) :: nz, imin,imax,jmin,jmax
integer, intent(in) :: ia(:),ja(:)
real(kind(1.d0)), intent(in) :: val(:)
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
integer, intent(in), optional :: gtl(:)
logical, optional, intent(in) :: rebuild
end subroutine psb_dcoins
subroutine psb_zcoins(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl,rebuild)
use psb_spmat_type
integer, intent(in) :: nz, imin,imax,jmin,jmax
integer, intent(in) :: ia(:),ja(:)
complex(kind(1.d0)), intent(in) :: val(:)
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
integer, intent(in), optional :: gtl(:)
logical, optional, intent(in) :: rebuild
end subroutine psb_zcoins
end interface
interface psb_symbmm
subroutine psb_dsymbmm(a,b,c,info)
use psb_spmat_type
type(psb_dspmat_type) :: a,b,c
integer :: info
end subroutine psb_dsymbmm
subroutine psb_zsymbmm(a,b,c,info)
use psb_spmat_type
type(psb_zspmat_type) :: a,b,c
integer :: info
end subroutine psb_zsymbmm
subroutine psb_dsymbmm(a,b,c,info)
use psb_spmat_type
type(psb_dspmat_type) :: a,b,c
integer :: info
end subroutine psb_dsymbmm
subroutine psb_zsymbmm(a,b,c,info)
use psb_spmat_type
type(psb_zspmat_type) :: a,b,c
integer :: info
end subroutine psb_zsymbmm
end interface
interface psb_numbmm
subroutine psb_dnumbmm(a,b,c)
use psb_spmat_type
type(psb_dspmat_type) :: a,b,c
end subroutine psb_dnumbmm
subroutine psb_znumbmm(a,b,c)
use psb_spmat_type
type(psb_zspmat_type) :: a,b,c
end subroutine psb_znumbmm
subroutine psb_dnumbmm(a,b,c)
use psb_spmat_type
type(psb_dspmat_type) :: a,b,c
end subroutine psb_dnumbmm
subroutine psb_znumbmm(a,b,c)
use psb_spmat_type
type(psb_zspmat_type) :: a,b,c
end subroutine psb_znumbmm
end interface
interface psb_transp
subroutine psb_dtransp(a,b,c,fmt)
use psb_spmat_type
type(psb_dspmat_type) :: a,b
integer, optional :: c
character(len=*), optional :: fmt
end subroutine psb_dtransp
subroutine psb_ztransp(a,b,c,fmt)
use psb_spmat_type
type(psb_zspmat_type) :: a,b
integer, optional :: c
character(len=*), optional :: fmt
end subroutine psb_ztransp
subroutine psb_dtransp(a,b,c,fmt)
use psb_spmat_type
type(psb_dspmat_type) :: a,b
integer, optional :: c
character(len=*), optional :: fmt
end subroutine psb_dtransp
subroutine psb_ztransp(a,b,c,fmt)
use psb_spmat_type
type(psb_zspmat_type) :: a,b
integer, optional :: c
character(len=*), optional :: fmt
end subroutine psb_ztransp
end interface
interface psb_transc
subroutine psb_ztransc(a,b,c,fmt)
use psb_spmat_type
type(psb_zspmat_type) :: a,b
integer, optional :: c
character(len=*), optional :: fmt
end subroutine psb_ztransc
subroutine psb_ztransc(a,b,c,fmt)
use psb_spmat_type
type(psb_zspmat_type) :: a,b
integer, optional :: c
character(len=*), optional :: fmt
end subroutine psb_ztransc
end interface
interface psb_rwextd
subroutine psb_drwextd(nr,a,info,b)
use psb_spmat_type
integer, intent(in) :: nr
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
type(psb_dspmat_type), intent(in), optional :: b
end subroutine psb_drwextd
subroutine psb_zrwextd(nr,a,info,b)
use psb_spmat_type
integer, intent(in) :: nr
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
type(psb_zspmat_type), intent(in), optional :: b
end subroutine psb_zrwextd
subroutine psb_drwextd(nr,a,info,b)
use psb_spmat_type
integer, intent(in) :: nr
type(psb_dspmat_type), intent(inout) :: a
integer, intent(out) :: info
type(psb_dspmat_type), intent(in), optional :: b
end subroutine psb_drwextd
subroutine psb_zrwextd(nr,a,info,b)
use psb_spmat_type
integer, intent(in) :: nr
type(psb_zspmat_type), intent(inout) :: a
integer, intent(out) :: info
type(psb_zspmat_type), intent(in), optional :: b
end subroutine psb_zrwextd
end interface
interface psb_csnmi
real(kind(1.d0)) function psb_dcsnmi(a,info,trans)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(out) :: info
character, optional :: trans
end function psb_dcsnmi
real(kind(1.d0)) function psb_zcsnmi(a,info,trans)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(out) :: info
character, optional :: trans
end function psb_zcsnmi
real(kind(1.d0)) function psb_dcsnmi(a,info,trans)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(out) :: info
character, optional :: trans
end function psb_dcsnmi
real(kind(1.d0)) function psb_zcsnmi(a,info,trans)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(out) :: info
character, optional :: trans
end function psb_zcsnmi
end interface
interface psb_sp_getdiag
subroutine psb_dspgtdiag(a,d,info)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
real(kind(1.d0)), intent(inout) :: d(:)
integer, intent(out) :: info
end subroutine psb_dspgtdiag
subroutine psb_zspgtdiag(a,d,info)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
complex(kind(1.d0)), intent(inout) :: d(:)
integer, intent(out) :: info
end subroutine psb_zspgtdiag
subroutine psb_dspgtdiag(a,d,info)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
real(kind(1.d0)), intent(inout) :: d(:)
integer, intent(out) :: info
end subroutine psb_dspgtdiag
subroutine psb_zspgtdiag(a,d,info)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
complex(kind(1.d0)), intent(inout) :: d(:)
integer, intent(out) :: info
end subroutine psb_zspgtdiag
end interface
interface psb_sp_scal
subroutine psb_dspscal(a,d,info)
use psb_spmat_type
type(psb_dspmat_type), intent(inout) :: a
real(kind(1.d0)), intent(in) :: d(:)
integer, intent(out) :: info
end subroutine psb_dspscal
subroutine psb_zspscal(a,d,info)
use psb_spmat_type
type(psb_zspmat_type), intent(inout) :: a
complex(kind(1.d0)), intent(in) :: d(:)
integer, intent(out) :: info
end subroutine psb_zspscal
subroutine psb_dspscal(a,d,info)
use psb_spmat_type
type(psb_dspmat_type), intent(inout) :: a
real(kind(1.d0)), intent(in) :: d(:)
integer, intent(out) :: info
end subroutine psb_dspscal
subroutine psb_zspscal(a,d,info)
use psb_spmat_type
type(psb_zspmat_type), intent(inout) :: a
complex(kind(1.d0)), intent(in) :: d(:)
integer, intent(out) :: info
end subroutine psb_zspscal
end interface
interface psb_sp_getblk
subroutine psb_dspgtblk(irw,a,b,info,append,iren,lrw)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: irw
type(psb_dspmat_type), intent(inout) :: b
logical, intent(in), optional :: append
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_dspgtblk
subroutine psb_zspgtblk(irw,a,b,info,append,iren,lrw)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: irw
type(psb_zspmat_type), intent(inout) :: b
logical, intent(in), optional :: append
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_zspgtblk
subroutine psb_dspgtblk(irw,a,b,info,append,iren,lrw)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: irw
type(psb_dspmat_type), intent(inout) :: b
logical, intent(in), optional :: append
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_dspgtblk
subroutine psb_zspgtblk(irw,a,b,info,append,iren,lrw)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: irw
type(psb_zspmat_type), intent(inout) :: b
logical, intent(in), optional :: append
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_zspgtblk
end interface
interface psb_sp_getrow
subroutine psb_dspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: irw
integer, intent(out) :: nz
integer, intent(inout) :: ia(:), ja(:)
real(kind(1.d0)), intent(inout) :: val(:)
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_dspgetrow
subroutine psb_zspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: irw
integer, intent(out) :: nz
integer, intent(inout) :: ia(:), ja(:)
complex(kind(1.d0)), intent(inout) :: val(:)
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_zspgetrow
subroutine psb_dspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: irw
integer, intent(out) :: nz
integer, intent(inout) :: ia(:), ja(:)
real(kind(1.d0)), intent(inout) :: val(:)
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_dspgetrow
subroutine psb_zspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: irw
integer, intent(out) :: nz
integer, intent(inout) :: ia(:), ja(:)
complex(kind(1.d0)), intent(inout) :: val(:)
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_zspgetrow
end interface
interface csrp
subroutine dcsrp(trans,m,n,fida,descra,ia1,ia2,&
& infoa,p,work,lwork,ierror)
integer, intent(in) :: m, n, lwork
integer, intent(out) :: ierror
character, intent(in) :: trans
double precision, intent(inout) :: work(*)
integer, intent(in) :: p(*)
integer, intent(inout) :: ia1(*), ia2(*), infoa(*)
character, intent(in) :: fida*5, descra*11
end subroutine dcsrp
subroutine zcsrp(trans,m,n,fida,descra,ia1,ia2,&
& infoa,p,work,lwork,ierror)
integer, intent(in) :: m, n, lwork
integer, intent(out) :: ierror
character, intent(in) :: trans
complex(kind(1.d0)), intent(inout) :: work(*)
integer, intent(in) :: p(*)
integer, intent(inout) :: ia1(*), ia2(*), infoa(*)
character, intent(in) :: fida*5, descra*11
end subroutine zcsrp
end interface
interface isaperm
logical function isaperm(n,ip)
integer, intent(in) :: n
integer, intent(inout) :: ip(*)
end function isaperm
end interface
interface psb_cest
subroutine psb_cest(afmt, m,n,nnz, lia1, lia2, lar, iup, info)
integer, intent(in) :: m,n,nnz,iup
integer, intent(out) :: lia1, lia2, lar, info
character(len=5) :: afmt
end subroutine psb_cest
end interface
end module psb_serial_mod

@ -557,22 +557,20 @@ Module psb_tools_mod
interface psb_glob_to_loc
subroutine psb_glob_to_loc2(x,y,desc_a,info,iact,owned)
subroutine psb_glob_to_loc2(x,y,desc_a,info,iact)
use psb_descriptor_type
type(psb_desc_type), intent(in) :: desc_a
integer,intent(in) :: x(:)
integer,intent(out) :: y(:)
integer, intent(out) :: info
character, intent(in), optional :: iact
logical, intent(in), optional :: owned
type(psb_desc_type), intent(in) :: desc_a
integer,intent(in) :: x(:)
integer,intent(out) :: y(:)
integer, intent(out) :: info
character, intent(in), optional :: iact
end subroutine psb_glob_to_loc2
subroutine psb_glob_to_loc(x,desc_a,info,iact,owned)
subroutine psb_glob_to_loc(x,desc_a,info,iact)
use psb_descriptor_type
type(psb_desc_type), intent(in) :: desc_a
integer,intent(inout) :: x(:)
integer, intent(out) :: info
character, intent(in), optional :: iact
logical, intent(in), optional :: owned
type(psb_desc_type), intent(in) :: desc_a
integer,intent(inout) :: x(:)
integer, intent(out) :: info
character, intent(in), optional :: iact
end subroutine psb_glob_to_loc
end interface

@ -77,15 +77,15 @@ module psi_mod
end interface
interface
subroutine psi_desc_index(desc,index_in,dep_list,&
& length_dl,nsnd,nrcv,desc_index,isglob_in,info)
use psb_descriptor_type
type(psb_desc_type) :: desc
integer :: index_in(:),dep_list(:)
integer,allocatable :: desc_index(:)
integer :: length_dl,nsnd,nrcv,info
logical :: isglob_in
end subroutine psi_desc_index
subroutine psi_desc_index(desc_data,index_in,dep_list,&
& length_dl,nsnd,nrcv,loc_to_glob,glob_to_loc,desc_index,&
& isglob_in,info)
integer :: desc_data(:),index_in(:),dep_list(:)
integer :: loc_to_glob(:),glob_to_loc(:)
integer,allocatable, intent(inout) :: desc_index(:)
integer :: length_dl,nsnd,nrcv,info
logical :: isglob_in
end subroutine psi_desc_index
end interface
interface
@ -94,13 +94,6 @@ module psi_mod
end subroutine psi_sort_dl
end interface
interface
subroutine psi_dl_check(dep_list,dl_lda,np,length_dl)
integer :: np,dl_lda,length_dl(0:np)
integer :: dep_list(dl_lda,0:np)
end subroutine psi_dl_check
end interface
interface psi_swapdata
subroutine psi_dswapdatam(flag,n,beta,y,desc_a,work,info,data)
use psb_descriptor_type
@ -430,7 +423,6 @@ contains
call psb_errpush(4010,name,a_err='psi_crea_ovr_elem')
goto 9999
end if
cdesc%matrix_data(psb_ovl_state_)=psb_cd_ovl_asb_
! finally bnd_elem
call psi_crea_bnd_elem(idx_out,cdesc,info)

@ -6,10 +6,10 @@ OBJS= psb_ddot.o psb_damax.o psb_dasum.o psb_daxpby.o\
pdtreecomb.o \
psb_zamax.o psb_zasum.o psb_zaxpby.o psb_zdot.o \
psb_znrm2.o psb_znrmi.o psb_zspmm.o psb_zspsm.o
LIBDIR=../../lib
HERE=.
INCDIRS=-I. -I.. -I$(LIBDIR)
LIBDIR = ..
MODDIR = ../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I .
lib: $(OBJS)

@ -14,8 +14,9 @@ FOBJS = psb_cest.o psb_dcoins.o psb_dcsdp.o psb_dcsmm.o psb_dcsmv.o \
psb_zrwextd.o psb_zsymbmm.o psb_znumbmm.o psb_zspscal.o\
psb_getifield.o psb_setifield.o psb_update_mod.o
INCDIRS = -I ../../lib -I .
LIBDIR = ../../lib
LIBDIR = ..
MODDIR = ../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I .
lib: auxd cood csrd jadd f77d dpd lib1
$(AR) $(LIBDIR)/$(LIBNAME) $(FOBJS)
@ -53,11 +54,4 @@ clean:
(cd dp; make clean)
(cd f77; make clean)
veryclean:
/bin/rm -f $(FOBJS) *$(.mod)
(cd aux; make veryclean)
(cd coo; make veryclean)
(cd csr; make veryclean)
(cd jad; make veryclean)
(cd dp; make veryclean)
(cd f77; make veryclean)
veryclean: clean

@ -12,9 +12,10 @@ OBJS=$(FOBJS)
#
# Where the library should go, and how it is called.
# Note that we are regenerating most of libsparker.a on the fly.
LIBDIR=../../../lib
#LIBNAME=libsparker.a
INCDIRS=-I. -I$(LIBDIR)
SPARKERDIR=..
LIBDIR = ../..
MODDIR = ../../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I . -I$(SPARKERDIR)
#
# No change should be needed below

@ -12,11 +12,12 @@ OBJS=$(FOBJS)
#
# Where the library should go, and how it is called.
# Note that we are regenerating most of libsparker.a on the fly.
LIBDIR=../../../lib
#LIBNAME=libsparker.a
LIBFILE=$(LIBDIR)/$(LIBNAME)
SPARKERDIR=..
INCDIRS=-I. -I$(SPARKERDIR) -I$(LIBDIR)
LIBDIR = ../..
MODDIR = ../../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I . -I$(SPARKERDIR)
LIBFILE=$(LIBDIR)/$(LIBNAME)
#
# No change should be needed below

@ -13,10 +13,11 @@ OBJS=$(FOBJS)
#
# Where the library should go, and how it is called.
# Note that we are regenerating most of libsparker.a on the fly.
LIBDIR=../../../lib
#LIBNAME=libsparker.a
SPARKERDIR=..
LIBDIR = ../..
MODDIR = ../../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I . -I$(SPARKERDIR)
LIBFILE=$(LIBDIR)/$(LIBNAME)
INCDIRS=-I. -I$(LIBDIR)
#
# No change should be needed below

@ -19,10 +19,12 @@ OBJS=$(FOBJS)
#
# Where the library should go, and how it is called.
# Note that we are regenerating most of libsparker.a on the fly.
LIBDIR=../../../lib
SPARKERDIR=..
LIBDIR = ../..
MODDIR = ../../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I . -I$(SPARKERDIR)
#LIBNAME=libsparker.a
LIBFILE=$(LIBDIR)/$(LIBNAME)
INCDIRS=-I. -I$(LIBDIR)
#
# No change should be needed below

@ -14,10 +14,12 @@ OBJS=$(FOBJS)
#
# Where the library should go, and how it is called.
# Note that we are regenerating most of libsparker.a on the fly.
LIBDIR=../../../lib
SPARKERDIR=..
LIBDIR = ../..
MODDIR = ../../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I . -I$(SPARKERDIR)
#LIBNAME=libsparker.a
LIBFILE=$(LIBDIR)/$(LIBNAME)
INCDIRS=-I. -I$(LIBDIR)
#
# No change should be needed below

@ -11,10 +11,12 @@ OBJS=$(FOBJS)
#
# Where the library should go, and how it is called.
# Note that we are regenerating most of libsparker.a on the fly.
LIBDIR=../../../lib
SPARKERDIR=..
LIBDIR = ../..
MODDIR = ../../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I . -I$(SPARKERDIR)
#LIBNAME=libsparker.a
LIBFILE=$(LIBDIR)/$(LIBNAME)
INCDIRS=-I. -I$(LIBDIR)
#
# No change should be needed below

@ -48,8 +48,6 @@ subroutine psb_dcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
use psb_error_mod
use psb_spmat_type
use psb_string_mod
use psb_serial_mod, only : psb_cest
implicit none
!....Parameters...
@ -72,6 +70,14 @@ subroutine psb_dcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
logical, parameter :: debug=.false.
character(len=20) :: name, ch_err
interface psb_cest
subroutine psb_cest(afmt, m,n,nnz, lia1, lia2, lar, iup, info)
integer, intent(in) :: m,n,nnz,iup
integer, intent(out) :: lia1, lia2, lar, info
character, intent(inout) :: afmt*5
end subroutine psb_cest
end interface
name='psb_csdp'
info = 0
call psb_erractionsave(err_act)

@ -57,7 +57,7 @@ subroutine psb_dcsprt(iout,a,iv,eirs,eics,head,ivr,ivc)
character(len=*), optional :: head
integer, intent(in), optional :: ivr(:), ivc(:)
character(len=*), parameter :: frmtr='(2(i16,1x),e16.8,2(i16,1x))'
character(len=*), parameter :: frmtr='(2(i6,1x),e16.8,2(i6,1x))'
integer :: irs,ics,i,j
if (present(eirs)) then

@ -40,15 +40,27 @@
!
subroutine psb_dnumbmm(a,b,c)
use psb_realloc_mod
use psb_spmat_type
use psb_serial_mod, only : psb_sp_getrow
implicit none
type(psb_dspmat_type) :: a,b,c
real(kind(1.d0)), allocatable :: temp(:)
integer :: info
interface psb_sp_getrow
subroutine psb_dspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: irw
integer, intent(out) :: nz
integer, intent(inout) :: ia(:), ja(:)
real(kind(1.d0)), intent(inout) :: val(:)
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_dspgetrow
end interface
allocate(temp(max(a%m,a%k,b%m,b%k)),stat=info)
if (info /= 0) then

@ -42,8 +42,6 @@
subroutine psb_dspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
use psb_spmat_type
use psb_string_mod
use psb_serial_mod, only: psb_sp_getblk
implicit none
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: irw
integer, intent(out) :: nz
@ -52,6 +50,23 @@ subroutine psb_dspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
interface psb_spgtblk
subroutine psb_dspgtblk(irw,a,b,info,append,iren,lrw)
! Output is always in COO format into B, irrespective of
! the input format
use psb_spmat_type
use psb_const_mod
implicit none
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: irw
type(psb_dspmat_type), intent(inout) :: b
integer,intent(out) :: info
logical, intent(in), optional :: append
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
end subroutine psb_dspgtblk
end interface
integer :: lrw_, ierr(5), err_act
type(psb_dspmat_type) :: b
@ -78,9 +93,9 @@ subroutine psb_dspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
call psb_sp_all(lrw_-irw+1,lrw_-irw+1,b,info)
if (present(iren)) then
call psb_sp_getblk(irw,a,b,info,iren=iren,lrw=lrw_)
call psb_spgtblk(irw,a,b,info,iren=iren,lrw=lrw_)
else
call psb_sp_getblk(irw,a,b,info,lrw=lrw_)
call psb_spgtblk(irw,a,b,info,lrw=lrw_)
end if
if (info /= 0) then
info=136

@ -45,13 +45,25 @@ subroutine psb_dspgtdiag(a,d,info)
use psb_spmat_type
use psb_error_mod
use psb_const_mod
use psb_serial_mod, only : psb_sp_getblk
implicit none
type(psb_dspmat_type), intent(in) :: a
real(kind(1.d0)), intent(inout) :: d(:)
integer, intent(out) :: info
interface psb_spgtblk
subroutine psb_dspgtblk(irw,a,b,info,append,iren,lrw)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: irw
type(psb_dspmat_type), intent(inout) :: b
logical, intent(in), optional :: append
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_dspgtblk
end interface
type(psb_dspmat_type) :: tmpa
integer :: i,j,k,nr, nz, err_act, ii, rng, irb, nrb
character(len=20) :: name, ch_err
@ -90,10 +102,10 @@ subroutine psb_dspgtdiag(a,d,info)
write(0,*)'in spgtdiag'
do i=1, rng, nrb
irb=min(i+nrb-1,rng)
call psb_sp_getblk(i,a,tmpa,info,lrw=irb)
call psb_spgtblk(i,a,tmpa,info,lrw=irb)
if(info.ne.0) then
info=4010
ch_err='psb_sp_getblk'
ch_err='psb_spgtblk'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if

@ -41,7 +41,6 @@
subroutine psb_dsymbmm(a,b,c,info)
use psb_spmat_type
use psb_string_mod
use psb_serial_mod, only : psb_sp_getrow
implicit none
type(psb_dspmat_type) :: a,b,c
@ -56,6 +55,19 @@ subroutine psb_dsymbmm(a,b,c,info)
integer, allocatable :: ic(:),jc(:)
end subroutine symbmm
end interface
interface psb_sp_getrow
subroutine psb_dspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
use psb_spmat_type
type(psb_dspmat_type), intent(in) :: a
integer, intent(in) :: irw
integer, intent(out) :: nz
integer, intent(inout) :: ia(:), ja(:)
real(kind(1.d0)), intent(inout) :: val(:)
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_dspgetrow
end interface
character(len=20) :: name, ch_err
integer :: err_act
@ -88,7 +100,6 @@ subroutine psb_dsymbmm(a,b,c,info)
endif
nze = max(a%m+1,2*a%m)
call psb_sp_reall(c,nze,info)
!
! Note: we need to test whether there is a performance impact
! in not using the original Douglas & Bank code.

@ -48,8 +48,6 @@ subroutine psb_zcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
use psb_error_mod
use psb_spmat_type
use psb_string_mod
use psb_serial_mod, only : psb_cest
implicit none
!....Parameters...
@ -72,6 +70,14 @@ subroutine psb_zcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
logical, parameter :: debug=.false.
character(len=20) :: name, ch_err
interface psb_cest
subroutine psb_cest(afmt, m,n,nnz, lia1, lia2, lar, iup, info)
integer, intent(in) :: m,n,nnz,iup
integer, intent(out) :: lia1, lia2, lar, info
character, intent(inout) :: afmt*5
end subroutine psb_cest
end interface
name='psb_csdp'
info = 0
call psb_erractionsave(err_act)

@ -40,15 +40,27 @@
!
subroutine psb_znumbmm(a,b,c)
use psb_realloc_mod
use psb_spmat_type
use psb_serial_mod, only : psb_sp_getrow
implicit none
type(psb_zspmat_type) :: a,b,c
complex(kind(1.d0)), allocatable :: temp(:)
integer :: info
interface psb_sp_getrow
subroutine psb_zspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: irw
integer, intent(out) :: nz
integer, intent(inout) :: ia(:), ja(:)
complex(kind(1.d0)), intent(inout) :: val(:)
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_zspgetrow
end interface
allocate(temp(max(a%m,a%k,b%m,b%k)),stat=info)
if (info /= 0) then

@ -42,8 +42,6 @@
subroutine psb_zspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
use psb_spmat_type
use psb_string_mod
use psb_serial_mod, only: psb_sp_getblk
implicit none
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: irw
integer, intent(out) :: nz
@ -52,6 +50,23 @@ subroutine psb_zspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
interface psb_spgtblk
subroutine psb_zspgtblk(irw,a,b,info,append,iren,lrw)
! Output is always in COO format into B, irrespective of
! the input format
use psb_spmat_type
use psb_const_mod
implicit none
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: irw
type(psb_zspmat_type), intent(inout) :: b
integer,intent(out) :: info
logical, intent(in), optional :: append
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
end subroutine psb_zspgtblk
end interface
integer :: lrw_, ierr(5), err_act
type(psb_zspmat_type) :: b
@ -78,9 +93,9 @@ subroutine psb_zspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
call psb_sp_all(lrw_-irw+1,lrw_-irw+1,b,info)
if (present(iren)) then
call psb_sp_getblk(irw,a,b,info,iren=iren,lrw=lrw_)
call psb_spgtblk(irw,a,b,info,iren=iren,lrw=lrw_)
else
call psb_sp_getblk(irw,a,b,info,lrw=lrw_)
call psb_spgtblk(irw,a,b,info,lrw=lrw_)
end if
if (info /= 0) then
info=136

@ -45,13 +45,25 @@ subroutine psb_zspgtdiag(a,d,info)
use psb_spmat_type
use psb_error_mod
use psb_const_mod
use psb_serial_mod, only : psb_sp_getblk
implicit none
type(psb_zspmat_type), intent(in) :: a
complex(kind(1.d0)), intent(inout) :: d(:)
integer, intent(out) :: info
interface psb_spgtblk
subroutine psb_zspgtblk(irw,a,b,info,append,iren,lrw)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: irw
type(psb_zspmat_type), intent(inout) :: b
logical, intent(in), optional :: append
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_zspgtblk
end interface
type(psb_zspmat_type) :: tmpa
integer :: i,j,k,nr, nz, err_act, ii, rng, irb, nrb
character(len=20) :: name, ch_err
@ -90,10 +102,10 @@ subroutine psb_zspgtdiag(a,d,info)
write(0,*)'in spgtdiag'
do i=1, rng, nrb
irb=min(i+nrb-1,rng)
call psb_sp_getblk(i,a,tmpa,info,lrw=irb)
call psb_spgtblk(i,a,tmpa,info,lrw=irb)
if(info.ne.0) then
info=4010
ch_err='psb_sp_getblk'
ch_err='psb_spgtblk'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if

@ -41,7 +41,6 @@
subroutine psb_zsymbmm(a,b,c,info)
use psb_spmat_type
use psb_string_mod
use psb_serial_mod, only : psb_sp_getrow
implicit none
type(psb_zspmat_type) :: a,b,c
@ -57,6 +56,19 @@ subroutine psb_zsymbmm(a,b,c,info)
end subroutine symbmm
end interface
interface psb_sp_getrow
subroutine psb_zspgetrow(irw,a,nz,ia,ja,val,info,iren,lrw)
use psb_spmat_type
type(psb_zspmat_type), intent(in) :: a
integer, intent(in) :: irw
integer, intent(out) :: nz
integer, intent(inout) :: ia(:), ja(:)
complex(kind(1.d0)), intent(inout) :: val(:)
integer, intent(in), target, optional :: iren(:)
integer, intent(in), optional :: lrw
integer, intent(out) :: info
end subroutine psb_zspgetrow
end interface
character(len=20) :: name, ch_err
integer :: err_act
name='psb_symbmm'
@ -101,7 +113,6 @@ subroutine psb_zsymbmm(a,b,c,info)
call inner_symbmm(a,b,c,itemp,info)
endif
call psb_realloc(size(c%ia1),c%aspk,info)
c%pl(1) = 0
c%pr(1) = 0
c%m=a%m

@ -15,8 +15,9 @@ FOBJS = psb_dallc.o psb_dasb.o psb_dcsrp.o psb_cdprt.o \
MPFOBJS = psb_dsphalo.o psb_zsphalo.o psb_cdasb.o psb_dcdovr.o psb_zcdovr.o
INCDIRS = -I ../../lib -I .
LIBDIR = ../../lib
LIBDIR = ..
MODDIR = ../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I .
lib: mpfobjs $(FOBJS)
$(AR) $(LIBDIR)/$(LIBNAME) $(MPFOBJS) $(FOBJS)
@ -26,8 +27,6 @@ lib: mpfobjs $(FOBJS)
mpfobjs:
(make $(MPFOBJS) F90="$(MPF90)" FC="$(MPF90)" FCOPT="$(F90COPT)")
psb_glob_to_loc.o: psb_glob_to_loc.f90
$(F90) $(F90COPT) $(INCDIRS) -c $< $(F90INLINE)
clean:
/bin/rm -f $(MPFOBJS) $(FOBJS)

@ -57,7 +57,7 @@ subroutine psb_cdall(m, n, parts, ictxt, desc_a, info)
!locals
Integer :: counter,i,j,np,me,loc_row,err,loc_col,nprocs,&
& l_ov_ix,l_ov_el,idx, err_act, itmpov, k, ns, glx
& l_ov_ix,l_ov_el,idx, err_act, itmpov, k, ns
integer :: int_err(5),exch(2)
integer, allocatable :: prc_v(:), temp_ovrlap(:), ov_idx(:),ov_el(:)
logical, parameter :: debug=.false.
@ -116,17 +116,13 @@ subroutine psb_cdall(m, n, parts, ictxt, desc_a, info)
!count local rows number
! allocate work vector
if (m > psb_cd_get_large_threshold()) then
allocate(desc_a%matrix_data(psb_mdata_size_),&
& temp_ovrlap(m),prc_v(np),stat=info)
else
allocate(desc_a%glob_to_loc(m),desc_a%matrix_data(psb_mdata_size_),&
& temp_ovrlap(m),prc_v(np),stat=info)
end if
if (info /= 0) then
allocate(prc_v(np),desc_a%glob_to_loc(m),&
&desc_a%matrix_data(psb_mdata_size_),temp_ovrlap(m),stat=info)
if (info /= no_err) then
info=2025
err=info
int_err(1)=m
call psb_errpush(info,name,i_err=int_err)
call psb_errpush(err,name,int_err)
goto 9999
endif
@ -135,189 +131,76 @@ subroutine psb_cdall(m, n, parts, ictxt, desc_a, info)
counter = 0
itmpov = 0
temp_ovrlap(:) = -1
if ( m >psb_cd_get_large_threshold()) then
desc_a%matrix_data(psb_dec_type_) = psb_desc_large_bld_
loc_col = (m+np-1)/np
allocate(desc_a%loc_to_glob(loc_col), desc_a%lprm(1),&
& desc_a%ptree(2),stat=info)
if (info == 0) call InitPairSearchTree(desc_a%ptree,info)
if (info /= 0) then
info=2025
int_err(1)=loc_col
call psb_errpush(info,name,i_err=int_err)
goto 9999
end if
! set LOC_TO_GLOB array to all "-1" values
desc_a%lprm(1) = 0
desc_a%loc_to_glob(:) = -1
k = 0
do i=1,m
if (info == 0) then
call parts(i,m,np,prc_v,nprocs)
if (nprocs > np) then
info=570
int_err(1)=3
int_err(2)=np
int_err(3)=nprocs
int_err(4)=i
err=info
call psb_errpush(err,name,int_err)
goto 9999
else if (nprocs <= 0) then
info=575
int_err(1)=3
int_err(2)=nprocs
int_err(3)=i
err=info
call psb_errpush(err,name,int_err)
goto 9999
else
do j=1,nprocs
if ((prc_v(j) > np-1).or.(prc_v(j) < 0)) then
info=580
int_err(1)=3
int_err(2)=prc_v(j)
int_err(3)=i
err=info
call psb_errpush(err,name,int_err)
goto 9999
end if
end do
endif
j=1
do
if (j > nprocs) exit
if (prc_v(j) == me) exit
j=j+1
enddo
if (j <= nprocs) then
if (prc_v(j) == me) then
! this point belongs to me
k = k + 1
call psb_check_size((k+1),desc_a%loc_to_glob,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
desc_a%loc_to_glob(k) = i
call SearchInsKeyVal(desc_a%ptree,i,k,glx,info)
if (nprocs > 1) then
call psb_check_size((itmpov+3+nprocs),temp_ovrlap,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
itmpov = itmpov + 1
temp_ovrlap(itmpov) = i
itmpov = itmpov + 1
temp_ovrlap(itmpov) = nprocs
temp_ovrlap(itmpov+1:itmpov+nprocs) = prc_v(1:nprocs)
itmpov = itmpov + nprocs
endif
do i=1,m
if (info == 0) then
call parts(i,m,np,prc_v,nprocs)
if (nprocs > np) then
info=570
int_err(1)=3
int_err(2)=np
int_err(3)=nprocs
int_err(4)=i
err=info
call psb_errpush(err,name,int_err)
goto 9999
else if (nprocs <= 0) then
info=575
int_err(1)=3
int_err(2)=nprocs
int_err(3)=i
err=info
call psb_errpush(err,name,int_err)
goto 9999
else
do j=1,nprocs
if ((prc_v(j) > np-1).or.(prc_v(j) < 0)) then
info=580
int_err(1)=3
int_err(2)=prc_v(j)
int_err(3)=i
err=info
call psb_errpush(err,name,int_err)
goto 9999
end if
end if
end if
enddo
if (info /= 0) then
info=4000
call psb_errpush(info,name)
goto 9999
endif
loc_row = k
else
desc_a%matrix_data(psb_dec_type_) = psb_desc_bld_
do i=1,m
if (info == 0) then
call parts(i,m,np,prc_v,nprocs)
if (nprocs > np) then
info=570
int_err(1)=3
int_err(2)=np
int_err(3)=nprocs
int_err(4)=i
err=info
call psb_errpush(err,name,int_err)
goto 9999
else if (nprocs <= 0) then
info=575
int_err(1)=3
int_err(2)=nprocs
int_err(3)=i
err=info
call psb_errpush(err,name,int_err)
goto 9999
else
do j=1,nprocs
if ((prc_v(j) > np-1).or.(prc_v(j) < 0)) then
info=580
int_err(1)=3
int_err(2)=prc_v(j)
int_err(3)=i
err=info
call psb_errpush(err,name,int_err)
goto 9999
end if
end do
endif
desc_a%glob_to_loc(i) = -(np+prc_v(1)+1)
j=1
do
if (j > nprocs) exit
if (prc_v(j) == me) exit
j=j+1
enddo
if (j <= nprocs) then
if (prc_v(j) == me) then
! this point belongs to me
counter=counter+1
desc_a%glob_to_loc(i) = counter
if (nprocs > 1) then
call psb_check_size((itmpov+3+nprocs),temp_ovrlap,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
end do
endif
desc_a%glob_to_loc(i) = -(np+prc_v(1)+1)
j=1
do
if (j > nprocs) exit
if (prc_v(j) == me) exit
j=j+1
enddo
if (j <= nprocs) then
if (prc_v(j) == me) then
! this point belongs to me
counter=counter+1
desc_a%glob_to_loc(i) = counter
if (nprocs > 1) then
if ((itmpov+2+nprocs) > size(temp_ovrlap)) then
ns = max(itmpov+2+nprocs,int(1.25*size(temp_ovrlap)))
call psb_realloc(ns,temp_ovrlap,info,pad=-1)
if (info /= 0) then
info=2025
int_err(1)=m
err=info
call psb_errpush(err,name,int_err)
goto 9999
end if
itmpov = itmpov + 1
temp_ovrlap(itmpov) = i
itmpov = itmpov + 1
temp_ovrlap(itmpov) = nprocs
temp_ovrlap(itmpov+1:itmpov+nprocs) = prc_v(1:nprocs)
itmpov = itmpov + nprocs
endif
endif
end if
itmpov = itmpov + 1
temp_ovrlap(itmpov) = i
itmpov = itmpov + 1
temp_ovrlap(itmpov) = nprocs
temp_ovrlap(itmpov+1:itmpov+nprocs) = prc_v(1:nprocs)
itmpov = itmpov + nprocs
endif
end if
endif
enddo
! estimate local cols number
loc_row=counter
loc_col=min(2*loc_row,m)
allocate(desc_a%loc_to_glob(loc_col),&
&desc_a%lprm(1),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
! set LOC_TO_GLOB array to all "-1" values
desc_a%lprm(1) = 0
desc_a%loc_to_glob(:) = -1
do i=1,m
k = desc_a%glob_to_loc(i)
if (k > 0) then
desc_a%loc_to_glob(k) = i
endif
enddo
end if
end if
endif
enddo
loc_row=counter
! check on parts function
if (debug) write(*,*) 'PSB_CDALL: End main loop:' ,loc_row,itmpov,info
@ -377,20 +260,36 @@ subroutine psb_cdall(m, n, parts, ictxt, desc_a, info)
call psb_transfer(ov_idx,desc_a%ovrlap_index,info)
if (info == 0) call psb_transfer(ov_el,desc_a%ovrlap_elem,info)
if (info == 0) deallocate(prc_v,temp_ovrlap,stat=info)
deallocate(prc_v,temp_ovrlap,stat=info)
if (info /= no_err) then
info=4000
err=info
call psb_errpush(err,name)
Goto 9999
endif
! At this point overlap_elem is OK.
desc_a%matrix_data(psb_ovl_state_) = psb_cd_ovl_asb_
! estimate local cols number
loc_col=min(2*loc_row,m)
allocate(desc_a%loc_to_glob(loc_col),&
&desc_a%lprm(1),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
! set LOC_TO_GLOB array to all "-1" values
desc_a%lprm(1) = 0
desc_a%loc_to_glob(:) = -1
do i=1,m
k = desc_a%glob_to_loc(i)
if (k > 0) then
desc_a%loc_to_glob(k) = i
endif
enddo
! set fields in desc_a%MATRIX_DATA....
desc_a%matrix_data(psb_n_row_) = loc_row
desc_a%matrix_data(psb_n_col_) = loc_row
call psb_cd_set_bld(desc_a,info)
call psb_realloc(1,desc_a%halo_index, info)
if (info /= no_err) then
@ -402,8 +301,10 @@ subroutine psb_cdall(m, n, parts, ictxt, desc_a, info)
desc_a%halo_index(:) = -1
desc_a%matrix_data(psb_m_) = m
desc_a%matrix_data(psb_n_) = n
desc_a%matrix_data(psb_dec_type_) = psb_desc_bld_
desc_a%matrix_data(psb_ctxt_) = ictxt
call psb_get_mpicomm(ictxt,desc_a%matrix_data(psb_mpi_c_))

@ -325,13 +325,12 @@ subroutine psb_cdalv(m, v, ictxt, desc_a, info, flag)
call psb_transfer(ov_idx,desc_a%ovrlap_index,info)
if (info == 0) call psb_transfer(ov_el,desc_a%ovrlap_elem,info)
if (info == 0) deallocate(temp_ovrlap,stat=info)
deallocate(temp_ovrlap,stat=info)
if (info /= 0) then
info=4000
call psb_errpush(info,name)
goto 9999
endif
desc_a%matrix_data(psb_ovl_state_) = psb_cd_ovl_asb_
! set fields in desc_a%MATRIX_DATA....
desc_a%matrix_data(psb_n_row_) = loc_row

@ -130,16 +130,6 @@ subroutine psb_cdasb(desc_a,info)
end if
if (psb_is_large_dec(dectype) ) then
if (allocated(desc_a%ptree)) then
call FreePairSearchTree(desc_a%ptree)
deallocate(desc_a%ptree,stat=info)
if (info /= 0) then
info=2059
call psb_errpush(info,name)
goto 9999
end if
end if
desc_a%matrix_data(psb_dec_type_) = psb_desc_large_asb_
!!$ write(0,*) 'Done large dec asmbly',desc_a%matrix_data(psb_dec_type_),&
!!$ & psb_desc_large_asb_,psb_is_asb_dec(desc_a%matrix_data(psb_dec_type_))

@ -117,7 +117,7 @@ subroutine psb_cddec(nloc, ictxt, desc_a, info)
!locals
Integer :: i,j,np,me,err,n,itmpov, k,&
& l_ov_ix,l_ov_el,idx, err_act,m, ip,glx
& l_ov_ix,l_ov_el,idx, err_act,m, ip
Integer :: INT_ERR(5), thalo(1), tovr(1)
integer, allocatable :: nlv(:)
logical, parameter :: debug=.false.
@ -164,107 +164,57 @@ subroutine psb_cddec(nloc, ictxt, desc_a, info)
!count local rows number
if ( m >psb_cd_get_large_threshold()) then
allocate(desc_a%loc_to_glob(nloc), desc_a%lprm(1),&
& desc_a%ptree(2),desc_a%matrix_data(psb_mdata_size_),stat=info)
if (info == 0) call InitPairSearchTree(desc_a%ptree,info)
if (info /= 0) then
info=2025
int_err(1)=nloc
call psb_errpush(info,name,i_err=int_err)
goto 9999
end if
! set LOC_TO_GLOB array to all "-1" values
desc_a%lprm(1) = 0
desc_a%loc_to_glob(:) = -1
desc_a%matrix_data(psb_n_row_) = nloc
desc_a%matrix_data(psb_n_col_) = nloc
desc_a%matrix_data(psb_m_) = m
desc_a%matrix_data(psb_n_) = m
desc_a%matrix_data(psb_dec_type_) = psb_desc_large_bld_
desc_a%matrix_data(psb_ctxt_) = ictxt
call psb_get_mpicomm(ictxt,desc_a%matrix_data(psb_mpi_c_))
do ip=0, np-1
if (ip==me) then
do i=1, nlv(ip)
call SearchInsKeyVal(desc_a%ptree,j,i,glx,info)
desc_a%loc_to_glob(i) = j
j = j + 1
enddo
else
do i=1, nlv(ip)
j = j + 1
enddo
endif
enddo
tovr = -1
thalo = -1
desc_a%lprm(:) = 0
call psi_cnv_dsc(thalo,tovr,desc_a,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psi_bld_cdesc')
goto 9999
end if
! allocate work vector
!!$ allocate(desc_a%glob_to_loc(m),desc_a%matrix_data(psb_mdata_size_),&
!!$ & desc_a%loc_to_glob(nloc),desc_a%lprm(1),&
!!$ & desc_a%ovrlap_index(1),desc_a%ovrlap_elem(1),&
!!$ & desc_a%halo_index(1),desc_a%bnd_elem(1),stat=info)
allocate(desc_a%glob_to_loc(m),desc_a%matrix_data(psb_mdata_size_),&
& desc_a%loc_to_glob(m),desc_a%lprm(1),stat=info)
if (info /= 0) then
info=2025
int_err(1)=m
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
desc_a%matrix_data(psb_dec_type_) = psb_desc_large_asb_
else
allocate(desc_a%glob_to_loc(m),desc_a%matrix_data(psb_mdata_size_),&
& desc_a%loc_to_glob(m),desc_a%lprm(1),stat=info)
if (info /= 0) then
info=2025
int_err(1)=m
call psb_errpush(info,name,i_err=int_err)
goto 9999
desc_a%matrix_data(psb_n_row_) = nloc
desc_a%matrix_data(psb_n_col_) = nloc
desc_a%matrix_data(psb_m_) = m
desc_a%matrix_data(psb_n_) = m
desc_a%matrix_data(psb_dec_type_) = psb_desc_bld_
desc_a%matrix_data(psb_ctxt_) = ictxt
call psb_get_mpicomm(ictxt,desc_a%matrix_data(psb_mpi_c_))
j = 1
do ip=0, np-1
if (ip==me) then
do i=1, nlv(ip)
desc_a%glob_to_loc(j) = i
desc_a%loc_to_glob(i) = j
j = j + 1
enddo
else
do i=1, nlv(ip)
desc_a%glob_to_loc(j) = -(np+ip+1)
j = j + 1
enddo
endif
enddo
tovr = -1
thalo = -1
desc_a%lprm(:) = 0
desc_a%matrix_data(psb_n_row_) = nloc
desc_a%matrix_data(psb_n_col_) = nloc
desc_a%matrix_data(psb_m_) = m
desc_a%matrix_data(psb_n_) = m
desc_a%matrix_data(psb_dec_type_) = psb_desc_bld_
desc_a%matrix_data(psb_ctxt_) = ictxt
call psb_get_mpicomm(ictxt,desc_a%matrix_data(psb_mpi_c_))
j = 1
do ip=0, np-1
if (ip==me) then
do i=1, nlv(ip)
desc_a%glob_to_loc(j) = i
desc_a%loc_to_glob(i) = j
j = j + 1
enddo
else
do i=1, nlv(ip)
desc_a%glob_to_loc(j) = -(np+ip+1)
j = j + 1
enddo
endif
enddo
tovr = -1
thalo = -1
desc_a%lprm(:) = 0
desc_a%matrix_data(psb_ovl_state_) = psb_cd_ovl_bld_
call psi_cnv_dsc(thalo,tovr,desc_a,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psi_bld_cdesc')
goto 9999
end if
desc_a%matrix_data(psb_dec_type_) = psb_desc_asb_
call psi_cnv_dsc(thalo,tovr,desc_a,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psi_bld_cdesc')
goto 9999
end if
endif
desc_a%matrix_data(psb_dec_type_) = psb_desc_asb_
call psb_erractionrestore(err_act)
return

@ -45,9 +45,15 @@ subroutine psb_cdren(trans,iperm,desc_a,info)
use psb_error_mod
use psb_penv_mod
use psb_string_mod
use psb_serial_mod
implicit none
interface isaperm
logical function isaperm(n,ip)
integer, intent(in) :: n
integer, intent(inout) :: ip(*)
end function isaperm
end interface
!...parameters....
type(psb_desc_type), intent(inout) :: desc_a
integer, intent(inout) :: iperm(:)
@ -57,6 +63,8 @@ subroutine psb_cdren(trans,iperm,desc_a,info)
integer :: i,j,np,me, n_col, kh, nh
integer :: dectype
integer :: ictxt,n_row, int_err(5), err_act
real(kind(1.d0)) :: time(10), mpi_wtime, real_err(6)
external mpi_wtime
logical, parameter :: debug=.false.
character(len=20) :: name
@ -65,6 +73,8 @@ subroutine psb_cdren(trans,iperm,desc_a,info)
call psb_erractionsave(err_act)
name = 'psb_dcren'
time(1) = mpi_wtime()
ictxt = psb_cd_get_context(desc_a)
dectype = psb_cd_get_dectype(desc_a)
n_row = psb_cd_get_local_rows(desc_a)
@ -197,6 +207,13 @@ subroutine psb_cdren(trans,iperm,desc_a,info)
endif
time(4) = mpi_wtime()
time(4) = time(4) - time(3)
if (debug) then
call psb_amx(ictxt, time(4))
write (*, *) ' comm structs assembly: ', time(4)*1.d-3
end if
call psb_erractionrestore(err_act)
return

@ -48,8 +48,6 @@ Subroutine psb_dcdovr(a,desc_a,novr,desc_ov,info)
use psb_serial_mod
use psb_descriptor_type
Use psb_prec_type
Use psb_prec_mod
use psb_error_mod
use psb_penv_mod
use psb_tools_mod
@ -157,165 +155,87 @@ Subroutine psb_dcdovr(a,desc_a,novr,desc_ov,info)
l_tmp_ovr_idx = novr*(3*Max(2*index_dim,1)+1)
l_tmp_halo = novr*(3*Size(desc_a%halo_index))
call psb_cd_set_bld(desc_ov,info)
desc_ov%matrix_data(psb_ovl_state_)=psb_cd_ovl_bld_
if (psb_is_large_desc(desc_a)) then
desc_ov%matrix_data(psb_dec_type_) = psb_desc_large_bld_
else
desc_ov%matrix_data(psb_dec_type_) = psb_desc_bld_
end if
If(debug) then
Write(0,*)'Start cdovrbld',me,lworks,lworkr
call psb_barrier(ictxt)
endif
if (.false.) then
!
! The real work goes on in here....
!
Call psb_cdovrbld(novr,desc_ov,desc_a,a,&
& l_tmp_halo,l_tmp_ovr_idx,lworks,lworkr,info)
if (info /= 0) then
info=4010
ch_err='psb_cdovrbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
If(debug) then
Write(0,*)'Done cdovrbld',me,lworks,lworkr
call psb_barrier(ictxt)
endif
Allocate(brvindx(np+1),rvsz(np),sdsz(np),bsdindx(np+1),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
Allocate(works(lworks),workr(lworkr),t_halo_in(l_tmp_halo),&
& t_halo_out(l_tmp_halo), temp(lworkr),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
call psb_sp_all(blk,max(lworks,lworkr),info)
if (info /= 0) then
info=4010
ch_err='psb_sp_all'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
blk%fida='COO'
Allocate(tmp_ovr_idx(l_tmp_ovr_idx),tmp_halo(l_tmp_halo),&
& halo(size(desc_a%halo_index)),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
halo(:) = desc_a%halo_index(:)
desc_ov%ovrlap_elem(:) = -1
tmp_ovr_idx(:) = -1
tmp_halo(:) = -1
counter_e = 1
tot_recv = 0
counter_h = 1
counter_o = 1
! Init overlap with desc_a%ovrlap (if any)
counter = 1
Do While (desc_a%ovrlap_index(counter) /= -1)
proc = desc_a%ovrlap_index(counter+psb_proc_id_)
n_elem_recv = desc_a%ovrlap_index(counter+psb_n_elem_recv_)
n_elem_send = desc_a%ovrlap_index(counter+n_elem_recv+psb_n_elem_send_)
Do j=0,n_elem_recv-1
idx = desc_a%ovrlap_index(counter+psb_elem_recv_+j)
If(idx > Size(desc_ov%loc_to_glob)) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
gidx = desc_ov%loc_to_glob(idx)
call psb_check_size((counter_o+3),tmp_ovr_idx,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
else
tmp_ovr_idx(counter_o)=proc
tmp_ovr_idx(counter_o+1)=1
tmp_ovr_idx(counter_o+2)=gidx
tmp_ovr_idx(counter_o+3)=-1
counter_o=counter_o+3
end Do
counter=counter+n_elem_recv+n_elem_send+2
end Do
Allocate(brvindx(np+1),rvsz(np),sdsz(np),bsdindx(np+1),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
Allocate(works(lworks),workr(lworkr),t_halo_in(l_tmp_halo),&
& t_halo_out(l_tmp_halo), temp(lworkr),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
!
! A picture is in order to understand what goes on here.
! I is the internal part; H is halo, R row, C column. The final
! matrix with N levels of overlap looks like this
!
! I | Hc1 | 0 | 0 |
! -------|-----|-----|-----|
! Hr1 | Hd1 | Hc2 | 0 |
! -------|-----|-----|-----|
! 0 | Hr2 | Hd2 | Hc2 |
! -------|-----|-----|-----|
! 0 | 0 | Hr3 | Hd3 | Hc3
!
! At the start we already have I and Hc1, so we know the row
! indices that will make up Hr1, and also who owns them. As we
! actually get those rows, we receive the column indices in Hc2;
! these define the row indices for Hr2, and so on. When we have
! reached the desired level HrN, we may ignore HcN.
!
!
Do i_ovr = 1, novr
call psb_sp_all(blk,max(lworks,lworkr),info)
if (info /= 0) then
info=4010
ch_err='psb_sp_all'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (debug) write(0,*) me,'Running on overlap level ',i_ovr,' of ',novr
blk%fida='COO'
!
! At this point, halo contains a valid halo corresponding to the
! matrix enlarged with the elements in the frontier for I_OVR-1.
! At the start, this is just the halo for A; the rows for indices in
! the first halo will contain column indices defining the second halo
! level and so on.
!
bsdindx(:) = 0
sdsz(:) = 0
brvindx(:) = 0
rvsz(:) = 0
idxr = 0
idxs = 0
counter = 1
counter_t = 1
Do While (halo(counter) /= -1)
tot_elem=0
proc=halo(counter+psb_proc_id_)
n_elem_recv=halo(counter+psb_n_elem_recv_)
n_elem_send=halo(counter+n_elem_recv+psb_n_elem_send_)
If ((counter+n_elem_recv+n_elem_send) > Size(halo)) then
info = -1
call psb_errpush(info,name)
goto 9999
end If
tot_recv=tot_recv+n_elem_recv
if (debug) write(0,*) me,' CDOVRBLD tot_recv:',proc,n_elem_recv,tot_recv
!
!
! The format of the halo vector exists in two forms: 1. Temporary
! 2. Assembled. In this loop we are using the (assembled) halo_in and
! copying it into (temporary) halo_out; this is because tmp_halo will
! be enlarged with the new column indices received, and will reassemble
! everything for the next iteration.
!
Allocate(tmp_ovr_idx(l_tmp_ovr_idx),tmp_halo(l_tmp_halo),&
& halo(size(desc_a%halo_index)),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
halo(:) = desc_a%halo_index(:)
desc_ov%ovrlap_elem(:) = -1
tmp_ovr_idx(:) = -1
tmp_halo(:) = -1
counter_e = 1
tot_recv = 0
counter_h = 1
counter_o = 1
! Init overlap with desc_a%ovrlap (if any)
counter = 1
Do While (desc_a%ovrlap_index(counter) /= -1)
proc = desc_a%ovrlap_index(counter+psb_proc_id_)
n_elem_recv = desc_a%ovrlap_index(counter+psb_n_elem_recv_)
n_elem_send = desc_a%ovrlap_index(counter+n_elem_recv+psb_n_elem_send_)
!
! add recv elements in halo_index into ovrlap_index
!
Do j=0,n_elem_recv-1
If((counter+psb_elem_recv_+j)>Size(halo)) then
info=-2
call psb_errpush(info,name)
goto 9999
end If
idx = halo(counter+psb_elem_recv_+j)
idx = desc_a%ovrlap_index(counter+psb_elem_recv_+j)
If(idx > Size(desc_ov%loc_to_glob)) then
info=-3
call psb_errpush(info,name)
@ -336,343 +256,443 @@ Subroutine psb_dcdovr(a,desc_a,novr,desc_ov,info)
tmp_ovr_idx(counter_o+2)=gidx
tmp_ovr_idx(counter_o+3)=-1
counter_o=counter_o+3
end Do
counter=counter+n_elem_recv+n_elem_send+2
end Do
if (.not.psb_is_large_desc(desc_ov)) then
call psb_check_size((counter_h+3),tmp_halo,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
tmp_halo(counter_h)=proc
tmp_halo(counter_h+1)=1
tmp_halo(counter_h+2)=idx
tmp_halo(counter_h+3)=-1
counter_h=counter_h+3
end if
Enddo
if (debug) write(0,*) me,'Checktmp_o_i Loop Mid1',tmp_ovr_idx(1:10)
counter = counter+n_elem_recv
!
! A picture is in order to understand what goes on here.
! I is the internal part; H is halo, R row, C column. The final
! matrix with N levels of overlap looks like this
!
! I | Hc1 | 0 | 0 |
! -------|-----|-----|-----|
! Hr1 | Hd1 | Hc2 | 0 |
! -------|-----|-----|-----|
! 0 | Hr2 | Hd2 | Hc2 |
! -------|-----|-----|-----|
! 0 | 0 | Hr3 | Hd3 | Hc3
!
! At the start we already have I and Hc1, so we know the row
! indices that will make up Hr1, and also who owns them. As we
! actually get those rows, we receive the column indices in Hc2;
! these define the row indices for Hr2, and so on. When we have
! reached the desired level HrN, we may ignore HcN.
!
!
Do i_ovr = 1, novr
if (debug) write(0,*) me,'Running on overlap level ',i_ovr,' of ',novr
!
! add send elements in halo_index into ovrlap_index
! At this point, halo contains a valid halo corresponding to the
! matrix enlarged with the elements in the frontier for I_OVR-1.
! At the start, this is just the halo for A; the rows for indices in
! the first halo will contain column indices defining the second halo
! level and so on.
!
Do j=0,n_elem_send-1
idx = halo(counter+psb_elem_send_+j)
gidx = desc_ov%loc_to_glob(idx)
if (idx > psb_cd_get_local_rows(Desc_a)) &
& write(0,*) me,i_ovr,'Out of local rows ',&
& idx,psb_cd_get_local_rows(Desc_a)
call psb_check_size((counter_o+3),tmp_ovr_idx,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
bsdindx(:) = 0
sdsz(:) = 0
brvindx(:) = 0
rvsz(:) = 0
idxr = 0
idxs = 0
counter = 1
counter_t = 1
Do While (halo(counter) /= -1)
tot_elem=0
proc=halo(counter+psb_proc_id_)
n_elem_recv=halo(counter+psb_n_elem_recv_)
n_elem_send=halo(counter+n_elem_recv+psb_n_elem_send_)
If ((counter+n_elem_recv+n_elem_send) > Size(halo)) then
info = -1
call psb_errpush(info,name)
goto 9999
end if
tmp_ovr_idx(counter_o)=proc
tmp_ovr_idx(counter_o+1)=1
tmp_ovr_idx(counter_o+2)=gidx
tmp_ovr_idx(counter_o+3)=-1
counter_o=counter_o+3
end If
tot_recv=tot_recv+n_elem_recv
if (debug) write(0,*) me,' CDOVRBLD tot_recv:',proc,n_elem_recv,tot_recv
!
!
! The format of the halo vector exists in two forms: 1. Temporary
! 2. Assembled. In this loop we are using the (assembled) halo_in and
! copying it into (temporary) halo_out; this is because tmp_halo will
! be enlarged with the new column indices received, and will reassemble
! everything for the next iteration.
!
!
! Prepare to exchange the halo rows with the other proc.
! add recv elements in halo_index into ovrlap_index
!
If (i_ovr < (novr)) Then
n_elem = psb_sp_get_nnz_row(idx,a)
Do j=0,n_elem_recv-1
If((counter+psb_elem_recv_+j)>Size(halo)) then
info=-2
call psb_errpush(info,name)
goto 9999
end If
call psb_check_size((idxs+tot_elem+n_elem),works,info)
idx = halo(counter+psb_elem_recv_+j)
If(idx > Size(desc_ov%loc_to_glob)) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
gidx = desc_ov%loc_to_glob(idx)
call psb_check_size((counter_o+3),tmp_ovr_idx,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
If((n_elem) > size(blk%ia2)) Then
isz = max((3*size(blk%ia2))/2,(n_elem))
if (debug) write(0,*) me,'Realloc blk',isz
call psb_sp_reall(blk,isz,info)
tmp_ovr_idx(counter_o)=proc
tmp_ovr_idx(counter_o+1)=1
tmp_ovr_idx(counter_o+2)=gidx
tmp_ovr_idx(counter_o+3)=-1
counter_o=counter_o+3
if (.not.psb_is_large_desc(desc_ov)) then
call psb_check_size((counter_h+3),tmp_halo,info,pad=-1)
if (info /= 0) then
info=4010
ch_err='psb_sp_reall'
call psb_errpush(info,name,a_err=ch_err)
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
End If
call psb_sp_getblk(idx,a,blk,info)
tmp_halo(counter_h)=proc
tmp_halo(counter_h+1)=1
tmp_halo(counter_h+2)=idx
tmp_halo(counter_h+3)=-1
counter_h=counter_h+3
end if
Enddo
if (debug) write(0,*) me,'Checktmp_o_i Loop Mid1',tmp_ovr_idx(1:10)
counter = counter+n_elem_recv
!
! add send elements in halo_index into ovrlap_index
!
Do j=0,n_elem_send-1
idx = halo(counter+psb_elem_send_+j)
gidx = desc_ov%loc_to_glob(idx)
if (idx > psb_cd_get_local_rows(Desc_a)) &
& write(0,*) me,i_ovr,'Out of local rows ',&
& idx,psb_cd_get_local_rows(Desc_a)
call psb_check_size((counter_o+3),tmp_ovr_idx,info,pad=-1)
if (info /= 0) then
info=4010
ch_err='psb_sp_getblk'
call psb_errpush(info,name,a_err=ch_err)
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
!!$ write(0,*) me,'Iteration: ',j,i_ovr
Do jj=1,n_elem
works(idxs+tot_elem+jj)=desc_ov%loc_to_glob(blk%ia2(jj))
End Do
tot_elem=tot_elem+n_elem
End If
Enddo
tmp_ovr_idx(counter_o)=proc
tmp_ovr_idx(counter_o+1)=1
tmp_ovr_idx(counter_o+2)=gidx
tmp_ovr_idx(counter_o+3)=-1
counter_o=counter_o+3
!
! Prepare to exchange the halo rows with the other proc.
!
If (i_ovr < (novr)) Then
n_elem = psb_sp_get_nnz_row(idx,a)
if (i_ovr < novr) then
if (tot_elem > 1) then
call imsr(tot_elem,works(idxs+1))
lx = works(idxs+1)
i = 1
j = 1
do
j = j + 1
if (j > tot_elem) exit
if (works(idxs+j) /= lx) then
i = i + 1
works(idxs+i) = works(idxs+j)
lx = works(idxs+i)
call psb_check_size((idxs+tot_elem+n_elem),works,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
end do
tot_elem = i
endif
if (debug) write(0,*) me,'Checktmp_o_i Loop Mid2',tmp_ovr_idx(1:10)
sdsz(proc+1) = tot_elem
idxs = idxs + tot_elem
end if
counter = counter+n_elem_send+3
if (debug) write(0,*) me,'Checktmp_o_i Loop End',tmp_ovr_idx(1:10)
Enddo
if (debug) write(0,*)me,'End phase 1 CDOVRBLD', m, n_col, tot_recv
if (i_ovr < novr) then
!
! Exchange data requests with everybody else: so far we have
! accumulated RECV requests, we have an all-to-all to build
! matchings SENDs.
!
call mpi_alltoall(sdsz,1,mpi_integer,rvsz,1,mpi_integer,icomm,info)
if (info /= 0) then
info=4010
ch_err='mpi_alltoall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
idxs = 0
idxr = 0
counter = 1
Do
proc=halo(counter)
if (proc == -1) exit
n_elem_recv = halo(counter+psb_n_elem_recv_)
counter = counter+n_elem_recv
n_elem_send = halo(counter+psb_n_elem_send_)
bsdindx(proc+1) = idxs
idxs = idxs + sdsz(proc+1)
brvindx(proc+1) = idxr
idxr = idxr + rvsz(proc+1)
If((n_elem) > size(blk%ia2)) Then
isz = max((3*size(blk%ia2))/2,(n_elem))
if (debug) write(0,*) me,'Realloc blk',isz
call psb_sp_reall(blk,isz,info)
if (info /= 0) then
info=4010
ch_err='psb_sp_reall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
End If
call psb_sp_getblk(idx,a,blk,info)
if (info /= 0) then
info=4010
ch_err='psb_sp_getblk'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
!!$ write(0,*) me,'Iteration: ',j,i_ovr
Do jj=1,n_elem
works(idxs+tot_elem+jj)=desc_ov%loc_to_glob(blk%ia2(jj))
End Do
tot_elem=tot_elem+n_elem
End If
Enddo
if (i_ovr < novr) then
if (tot_elem > 1) then
call imsr(tot_elem,works(idxs+1))
lx = works(idxs+1)
i = 1
j = 1
do
j = j + 1
if (j > tot_elem) exit
if (works(idxs+j) /= lx) then
i = i + 1
works(idxs+i) = works(idxs+j)
lx = works(idxs+i)
end if
end do
tot_elem = i
endif
if (debug) write(0,*) me,'Checktmp_o_i Loop Mid2',tmp_ovr_idx(1:10)
sdsz(proc+1) = tot_elem
idxs = idxs + tot_elem
end if
counter = counter+n_elem_send+3
if (debug) write(0,*) me,'Checktmp_o_i Loop End',tmp_ovr_idx(1:10)
Enddo
iszr=sum(rvsz)
if (max(iszr,1) > lworkr) then
call psb_realloc(max(iszr,1),workr,info)
if (debug) write(0,*)me,'End phase 1 CDOVRBLD', m, n_col, tot_recv
if (i_ovr < novr) then
!
! Exchange data requests with everybody else: so far we have
! accumulated RECV requests, we have an all-to-all to build
! matchings SENDs.
!
call mpi_alltoall(sdsz,1,mpi_integer,rvsz,1,mpi_integer,icomm,info)
if (info /= 0) then
info=4010
ch_err='psb_realloc'
ch_err='mpi_alltoall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
lworkr=max(iszr,1)
end if
call mpi_alltoallv(works,sdsz,bsdindx,mpi_integer,&
& workr,rvsz,brvindx,mpi_integer,icomm,info)
if (info /= 0) then
info=4010
ch_err='mpi_alltoallv'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (debug) write(0,*) 'ISZR :',iszr
idxs = 0
idxr = 0
counter = 1
Do
proc=halo(counter)
if (proc == -1) exit
n_elem_recv = halo(counter+psb_n_elem_recv_)
counter = counter+n_elem_recv
n_elem_send = halo(counter+psb_n_elem_send_)
bsdindx(proc+1) = idxs
idxs = idxs + sdsz(proc+1)
brvindx(proc+1) = idxr
idxr = idxr + rvsz(proc+1)
counter = counter+n_elem_send+3
Enddo
iszr=sum(rvsz)
if (max(iszr,1) > lworkr) then
call psb_realloc(max(iszr,1),workr,info)
if (info /= 0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
lworkr=max(iszr,1)
end if
if (psb_is_large_desc(desc_a)) then
call psb_check_size(iszr,maskr,info)
call mpi_alltoallv(works,sdsz,bsdindx,mpi_integer,&
& workr,rvsz,brvindx,mpi_integer,icomm,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
ch_err='mpi_alltoallv'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psi_idx_cnv(iszr,workr,maskr,desc_ov,info)
iszs = count(maskr(1:iszr)<=0)
if (iszs > size(works)) call psb_realloc(iszs,works,info)
j = 0
do i=1,iszr
if (maskr(i) < 0) then
j=j+1
works(j) = workr(i)
end if
end do
!
! fnd_owner on desc_a because we want the procs who
! owned the rows from the beginning!
!
call psi_fnd_owner(iszs,works,temp,desc_a,info)
n_col=psb_cd_get_local_cols(desc_ov)
do i=1,iszs
idx = works(i)
n_col = psb_cd_get_local_cols(desc_ov)
call psi_idx_ins_cnv(idx,lidx,desc_ov,info)
if (psb_cd_get_local_cols(desc_ov) > n_col ) then
!
! This is a new index. Assigning a local index as
! we receive them guarantees that all indices for HALO(I)
! will be less than those for HALO(J) whenever I<J
!
proc_id = temp(i)
call psb_check_size((counter_t+3),t_halo_in,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
t_halo_in(counter_t)=proc_id
t_halo_in(counter_t+1)=1
t_halo_in(counter_t+2)=lidx
counter_t=counter_t+3
if (.false.) write(0,*) me,' CDOVRBLD: Added t_halo_in ',&
&proc_id,lidx,idx
endif
end Do
else
Do i=1,iszr
idx=workr(i)
if (idx <1) then
write(0,*) me,'Error in CDOVRBLD level',i_ovr,' : ',idx,i,iszr
else If (desc_ov%glob_to_loc(idx) < -np) Then
!
! This is a new index. Assigning a local index as
! we receive them guarantees that all indices for HALO(I)
! will be less than those for HALO(J) whenever I<J
!
n_col=n_col+1
proc_id=-desc_ov%glob_to_loc(idx)-np-1
call psb_check_size(n_col,desc_ov%loc_to_glob,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
desc_ov%glob_to_loc(idx)=n_col
desc_ov%loc_to_glob(n_col)=idx
if (debug) write(0,*) 'ISZR :',iszr
call psb_check_size((counter_t+3),t_halo_in,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
if (psb_is_large_desc(desc_a)) then
call psb_check_size(iszr,maskr,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
call psi_idx_cnv(iszr,workr,maskr,desc_ov,info)
iszs = count(maskr(1:iszr)<=0)
if (iszs > size(works)) call psb_realloc(iszs,works,info)
j = 0
do i=1,iszr
if (maskr(i) < 0) then
j=j+1
works(j) = workr(i)
end if
end do
!
! fnd_owner on desc_a because we want the procs who
! owned the rows from the beginning!
!
call psi_fnd_owner(iszs,works,temp,desc_a,info)
n_col=psb_cd_get_local_cols(desc_ov)
do i=1,iszs
idx = works(i)
n_col = psb_cd_get_local_cols(desc_ov)
call psi_idx_ins_cnv(idx,lidx,desc_ov,info)
if (psb_cd_get_local_cols(desc_ov) > n_col ) then
!
! This is a new index. Assigning a local index as
! we receive them guarantees that all indices for HALO(I)
! will be less than those for HALO(J) whenever I<J
!
proc_id = temp(i)
call psb_check_size((counter_t+3),t_halo_in,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
t_halo_in(counter_t)=proc_id
t_halo_in(counter_t+1)=1
t_halo_in(counter_t+2)=lidx
counter_t=counter_t+3
if (.false.) write(0,*) me,' CDOVRBLD: Added t_halo_in ',&
&proc_id,lidx,idx
endif
end Do
else
Do i=1,iszr
idx=workr(i)
if (idx <1) then
write(0,*) me,'Error in CDOVRBLD level',i_ovr,' : ',idx,i,iszr
else If (desc_ov%glob_to_loc(idx) < -np) Then
!
! This is a new index. Assigning a local index as
! we receive them guarantees that all indices for HALO(I)
! will be less than those for HALO(J) whenever I<J
!
n_col=n_col+1
proc_id=-desc_ov%glob_to_loc(idx)-np-1
call psb_check_size(n_col,desc_ov%loc_to_glob,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
desc_ov%glob_to_loc(idx)=n_col
desc_ov%loc_to_glob(n_col)=idx
call psb_check_size((counter_t+3),t_halo_in,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
t_halo_in(counter_t)=proc_id
t_halo_in(counter_t+1)=1
t_halo_in(counter_t+2)=n_col
counter_t=counter_t+3
if (debug) write(0,*) me,' CDOVRBLD: Added into t_halo_in from recv',&
&proc_id,n_col,idx
else if (desc_ov%glob_to_loc(idx) < 0) Then
if (debug) write(0,*) me,'Wrong input to cdovrbld??',&
&idx,desc_ov%glob_to_loc(idx)
End If
End Do
desc_ov%matrix_data(psb_n_col_)=n_col
t_halo_in(counter_t)=proc_id
t_halo_in(counter_t+1)=1
t_halo_in(counter_t+2)=n_col
counter_t=counter_t+3
if (debug) write(0,*) me,' CDOVRBLD: Added into t_halo_in from recv',&
&proc_id,n_col,idx
else if (desc_ov%glob_to_loc(idx) < 0) Then
if (debug) write(0,*) me,'Wrong input to cdovrbld??',&
&idx,desc_ov%glob_to_loc(idx)
End If
End Do
desc_ov%matrix_data(psb_n_col_)=n_col
end if
end if
end if
!!$ desc_ov%matrix_data(psb_n_row_)=desc_ov%matrix_data(psb_n_col_)
!
! Ok, now we have a temporary halo with all the info for the
! next round. If we need to keep going, convert the halo format
! from temporary to final, so that we can work out the next iteration.
! This uses one of the convert_comm internals, i.e. we are doing
! the equivalent of a partial call to convert_comm
!
!
! Ok, now we have a temporary halo with all the info for the
! next round. If we need to keep going, convert the halo format
! from temporary to final, so that we can work out the next iteration.
! This uses one of the convert_comm internals, i.e. we are doing
! the equivalent of a partial call to convert_comm
!
If (i_ovr < (novr)) Then
If (i_ovr < (novr)) Then
t_halo_in(counter_t)=-1
t_halo_in(counter_t)=-1
if (debug) write(0,*) me,'Checktmp_o_i 1',tmp_ovr_idx(1:10)
if (debug) write(0,*) me,'Calling Crea_Halo'
if (debug) write(0,*) me,'Checktmp_o_i 1',tmp_ovr_idx(1:10)
if (debug) write(0,*) me,'Calling Crea_Halo'
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,.false.,&
& nxch,nsnd,nrcv,info)
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,.false.,&
& nxch,nsnd,nrcv,info)
if (debug) then
write(0,*) me,'Done Crea_Index'
call psb_barrier(ictxt)
end if
if (debug) write(0,*) me,'Checktmp_o_i 2',tmp_ovr_idx(1:10)
if (debug) write(0,*) me,'Done Crea_Halo'
call psb_transfer(t_halo_out,halo,info)
!
! At this point we have built the halo necessary for I_OVR+1.
!
End If
if (debug) write(0,*) me,'Checktmp_o_i ',tmp_ovr_idx(1:10)
if (debug) then
write(0,*) me,'Done Crea_Index'
call psb_barrier(ictxt)
end if
if (debug) write(0,*) me,'Checktmp_o_i 2',tmp_ovr_idx(1:10)
if (debug) write(0,*) me,'Done Crea_Halo'
call psb_transfer(t_halo_out,halo,info)
!
! At this point we have built the halo necessary for I_OVR+1.
!
End If
if (debug) write(0,*) me,'Checktmp_o_i ',tmp_ovr_idx(1:10)
End Do
End Do
desc_ov%matrix_data(psb_m_)=psb_cd_get_global_rows(desc_a)
desc_ov%matrix_data(psb_n_)=psb_cd_get_global_cols(desc_a)
desc_ov%matrix_data(psb_m_)=psb_cd_get_global_rows(desc_a)
desc_ov%matrix_data(psb_n_)=psb_cd_get_global_cols(desc_a)
tmp_halo(counter_h:)=-1
tmp_ovr_idx(counter_o:)=-1
tmp_halo(counter_h:)=-1
tmp_ovr_idx(counter_o:)=-1
!
! At this point we have gathered all the indices in the halo at
! N levels of overlap. Just call cnv_dsc. This is
! the same routine as gets called inside CDASB.
!
!
! At this point we have gathered all the indices in the halo at
! N levels of overlap. Just call cnv_dsc. This is
! the same routine as gets called inside CDASB.
!
if (debug) then
write(0,*) 'psb_cdovrbld: converting indexes'
call psb_barrier(ictxt)
end if
!.... convert comunication stuctures....
! Note that we have to keep local_rows until the very end,
! because otherwise the halo build mechanism of cdasb
! will not work.
call psb_transfer(tmp_ovr_idx,desc_ov%ovrlap_index,info)
call psb_transfer(tmp_halo,desc_ov%halo_index,info)
call psb_cdasb(desc_ov,info)
desc_ov%matrix_data(psb_n_row_)=desc_ov%matrix_data(psb_n_col_)
if (debug) then
write(0,*) me,'Done CDASB'
call psb_barrier(ictxt)
end if
if (debug) then
write(0,*) 'psb_cdovrbld: converting indexes'
call psb_barrier(ictxt)
end if
!.... convert comunication stuctures....
! Note that we have to keep local_rows until the very end,
! because otherwise the halo build mechanism of cdasb
! will not work.
call psb_transfer(tmp_ovr_idx,desc_ov%ovrlap_index,info)
call psb_transfer(tmp_halo,desc_ov%halo_index,info)
call psb_cdasb(desc_ov,info)
desc_ov%matrix_data(psb_n_row_)=desc_ov%matrix_data(psb_n_col_)
if (debug) then
write(0,*) me,'Done CDASB'
call psb_barrier(ictxt)
end if
if (info == 0) call psb_sp_free(blk,info)
if (info /= 0) then
ch_err='sp_free'
call psb_errpush(4013,name,a_err=ch_err,i_err=(/info,0,0,0,0/))
goto 9999
end if
if (info == 0) call psb_sp_free(blk,info)
if (info /= 0) then
ch_err='sp_free'
call psb_errpush(4013,name,a_err=ch_err,i_err=(/info,0,0,0,0/))
goto 9999
end if
end if
call psb_erractionrestore(err_act)
return

@ -47,6 +47,28 @@ subroutine psb_dcsrp(trans,iperm,a, desc_a, info)
use psb_penv_mod
! implicit none
interface dcsrp
subroutine dcsrp(trans,m,n,fida,descra,ia1,ia2,&
& infoa,p,work,lwork,ierror)
integer, intent(in) :: m, n, lwork
integer, intent(out) :: ierror
character, intent(in) :: trans
double precision, intent(inout) :: work(*)
integer, intent(in) :: p(*)
integer, intent(inout) :: ia1(*), ia2(*), infoa(*)
character, intent(in) :: fida*5, descra*11
end subroutine dcsrp
end interface
interface isaperm
logical function isaperm(n,ip)
integer, intent(in) :: n
integer, intent(inout) :: ip(*)
end function isaperm
end interface
!...parameters....
type(psb_dspmat_type), intent(inout) :: a
@ -61,8 +83,12 @@ subroutine psb_dcsrp(trans,iperm,a, desc_a, info)
integer :: ictxt,n_row,err_act, int_err(5)
character(len=20) :: name, char_err
real(kind(1.d0)) :: time(10), mpi_wtime
external mpi_wtime
logical, parameter :: debug=.false.
time(1) = mpi_wtime()
ictxt = psb_cd_get_context(desc_a)
dectype = psb_cd_get_dectype(desc_a)
n_row = psb_cd_get_local_rows(desc_a)
@ -133,7 +159,7 @@ subroutine psb_dcsrp(trans,iperm,a, desc_a, info)
! hmm, maybe we should just move all of this onto a different level,
! have a specialized subroutine, and do it in the solver context????
if (debug) write(0,*) 'spasb: calling dcsrp',size(work_dcsdp)
call csrp(trans,n_row,n_col,a%fida,a%descra,a%ia1,a%ia2,a%infoa,&
call dcsrp(trans,n_row,n_col,a%fida,a%descra,a%ia1,a%ia2,a%infoa,&
& ipt,work_dcsdp,size(work_dcsdp),info)
if(info /= no_err) then
info=4010
@ -144,6 +170,13 @@ subroutine psb_dcsrp(trans,iperm,a, desc_a, info)
deallocate(ipt,work_dcsdp)
time(4) = mpi_wtime()
time(4) = time(4) - time(3)
if (debug) then
call psb_amx(ictxt, time(4))
write (*, *) ' comm structs assembly: ', time(4)*1.d-3
end if
call psb_erractionrestore(err_act)
return

@ -69,6 +69,15 @@ subroutine psb_dgelp(trans,iperm,x,desc_a,info)
integer, intent(in) :: p(*)
end subroutine dgelp
end interface
interface isaperm
logical function isaperm(n,ip)
integer, intent(in) :: n
integer, intent(inout) :: ip(*)
end function isaperm
end interface
character(len=20) :: name, ch_err
name = 'psb_dgelp'
@ -205,6 +214,14 @@ subroutine psb_dgelpv(trans,iperm,x,desc_a,info)
end subroutine dgelp
end interface
interface isaperm
logical function isaperm(n,ip)
integer, intent(in) :: n
integer, intent(inout) :: ip(*)
end function isaperm
end interface
character(len=20) :: name, ch_err
name = 'psb_dgelpv'

@ -73,6 +73,11 @@ subroutine psb_dinsvi(m, irw, val, x, desc_a, info, dupl)
call psb_erractionsave(err_act)
name = 'psb_dinsvi'
!!$ if (.not.allocated(desc_a%glob_to_loc)) then
!!$ info=3110
!!$ call psb_errpush(info,name)
!!$ return
!!$ end if
if ((.not.allocated(desc_a%matrix_data))) then
int_err(1)=3110
call psb_errpush(info,name)
@ -252,6 +257,11 @@ subroutine psb_dinsi(m, irw, val, x, desc_a, info, dupl)
call psb_erractionsave(err_act)
name = 'psb_dinsi'
!!$ if (.not.allocated(desc_a%glob_to_loc)) then
!!$ info=3110
!!$ call psb_errpush(info,name)
!!$ return
!!$ end if
if ((.not.allocated(desc_a%matrix_data))) then
int_err(1)=3110
call psb_errpush(info,name)

@ -49,6 +49,49 @@ subroutine psb_dspcnv(a,b,desc_a,info)
use psb_error_mod
use psb_penv_mod
implicit none
interface dcsdp
subroutine dcsdp(check,trans,m,n,unitd,d,&
& fida,descra,a,ia1,ia2,infoa,&
& pl,fidh,descrh,h,ih1,ih2,infoh,pr,lh,lh1,lh2,&
& work,lwork,ierror)
integer, intent(in) :: lh, lwork, lh1, lh2, m, n
integer, intent(out) :: ierror
character, intent(in) :: check, trans, unitd
real(kind(1.d0)), intent(in) :: d(*), a(*)
real(kind(1.d0)), intent(out) :: h(*)
real(kind(1.d0)), intent(inout) :: work(*)
integer, intent(in) :: ia1(*), ia2(*), infoa(*)
integer, intent(out) :: ih1(*), ih2(*), pl(*),pr(*), infoh(*)
character, intent(in) :: fida*5, descra*11
character, intent(out) :: fidh*5, descrh*11
end subroutine dcsdp
end interface
interface dcsrp
subroutine dcsrp(trans,m,n,fida,descra,ia1,ia2,&
& infoa,p,work,lwork,ierror)
integer, intent(in) :: m, n, lwork
integer, intent(out) :: ierror
character, intent(in) :: trans
real(kind(1.d0)), intent(inout) :: work(*)
integer, intent(in) :: p(*)
integer, intent(inout) :: ia1(*), ia2(*), infoa(*)
character, intent(in) :: fida*5, descra*11
end subroutine dcsrp
end interface
interface dcsprt
subroutine dcsprt(m,n,fida,descra,a,ia1,ia2,infoa ,iout,ierror)
integer, intent(in) :: iout,m, n
integer, intent(out) :: ierror
real(kind(1.d0)), intent(in) :: a(*)
integer, intent(in) :: ia1(*), ia2(*), infoa(*)
character, intent(in) :: fida*5, descra*11
end subroutine dcsprt
end interface
!...parameters....
type(psb_dspmat_type), intent(in) :: a
@ -57,11 +100,17 @@ subroutine psb_dspcnv(a,b,desc_a,info)
integer, intent(out) :: info
!....locals....
integer :: int_err(5)
integer :: ia1_size,ia2_size,aspk_size,err_act&
& ,i,err,np,me,n_col
integer, allocatable :: i_temp(:)
integer :: dectype
real(kind(1.d0)) :: d(1)
integer,allocatable :: i_temp(:)
real(kind(1.d0)),allocatable :: work_dcsdp(:)
integer :: ia1_size,ia2_size,aspk_size,&
& err_act,i,np,me,n_col,l_dcsdp
integer :: lwork_dcsdp,dectype
integer :: ictxt,n_row
character :: check*1, trans*1, unitd*1
real(kind(1.d0)) :: time(10), mpi_wtime
external mpi_wtime
logical, parameter :: debug=.false.
character(len=20) :: name, ch_err
@ -70,6 +119,7 @@ subroutine psb_dspcnv(a,b,desc_a,info)
name = 'psb_dspcnv'
call psb_erractionsave(err_act)
time(1) = mpi_wtime()
ictxt = psb_cd_get_context(desc_a)
dectype = psb_cd_get_dectype(desc_a)
@ -99,15 +149,46 @@ subroutine psb_dspcnv(a,b,desc_a,info)
if (debug) write (0, *) name,' sizes',ia1_size,ia2_size,aspk_size
! convert only without check
check='N'
trans='N'
unitd='U'
! l_dcsdp is the size requested for dcsdp procedure
l_dcsdp=(ia1_size+100)
b%m=n_row
b%k=n_col
call psb_sp_all(b,ia1_size,ia2_size,aspk_size,info)
allocate(work_dcsdp(l_dcsdp),stat=info)
if (info /= 0) then
info=2025
int_err(1)=l_dcsdp
call psb_errpush(info, name, i_err=int_err)
goto 9999
endif
call psb_csdp(a,b,info)
lwork_dcsdp=size(work_dcsdp)
! set infoa(1) to nnzero
b%pl(:) = 0
b%pr(:) = 0
if (debug) write (0, *) name,' calling dcsdp',lwork_dcsdp,&
&size(work_dcsdp)
! convert aspk,ia1,ia2 in requested representation mode
if (debug) then
endif
! result is put in b
call dcsdp(check,trans,n_row,n_col,unitd,d,a%fida,a%descra,&
& a%aspk,a%ia1,a%ia2,a%infoa,&
& b%pl,b%fida,b%descra,b%aspk,b%ia1,b%ia2,b%infoa,b%pr,&
& size(b%aspk),size(b%ia1),size(b%ia2),&
& work_dcsdp,size(work_dcsdp),info)
if(info /= no_err) then
info=4010
ch_err='psb_csdp'
ch_err='dcsdp'
call psb_errpush(info, name, a_err=ch_err)
goto 9999
end if
@ -147,6 +228,9 @@ subroutine psb_dspcnv(a,b,desc_a,info)
endif
if (debug) write (0, *) me,name,' from dcsdp ',&
&b%fida,' pl ', b%pl(:),'pr',b%pr(:)
call psb_erractionrestore(err_act)
return

@ -49,7 +49,6 @@ Subroutine psb_dsphalo(a,desc_a,blk,info,rwcnv,clcnv,outfmt)
use psb_serial_mod
use psb_descriptor_type
Use psb_prec_type
use psb_realloc_mod
use psb_tools_mod, only : psb_glob_to_loc, psb_loc_to_glob
use psb_error_mod

@ -52,7 +52,6 @@ subroutine psb_dspins(nz,ia,ja,val,a,desc_a,info,rebuild)
use psb_const_mod
use psb_error_mod
use psb_penv_mod
use psb_tools_mod
implicit none
!....parameters...
@ -71,6 +70,27 @@ subroutine psb_dspins(nz,ia,ja,val,a,desc_a,info,rebuild)
logical :: rebuild_
integer, allocatable :: ila(:),jla(:)
interface psb_cdins
subroutine psb_cdins(nz,ia,ja,desc_a,info,ila,jla)
use psb_descriptor_type
implicit none
type(psb_desc_type), intent(inout) :: desc_a
integer, intent(in) :: nz,ia(:),ja(:)
integer, intent(out) :: info
integer, optional, intent(out) :: ila(:), jla(:)
end subroutine psb_cdins
end interface
interface psb_glob_to_loc
subroutine psb_glob_to_loc(x,desc_a,info,iact)
use psb_descriptor_type
implicit none
type(psb_desc_type), intent(in) :: desc_a
integer, intent(inout) :: x(:)
integer, intent(out) :: info
character, intent(in), optional :: iact
end subroutine psb_glob_to_loc
end interface
character(len=20) :: name, ch_err
info = 0

@ -11,24 +11,25 @@ subroutine psb_get_ovrlap(ovrel,desc,info)
character(len=20) :: name
info = 0
name='psb_get_overlap'
name='psi_get_overlap'
call psb_erractionsave(err_act)
if (psb_is_ovl_asb(desc)) then
i=0
j=1
do while(desc%ovrlap_elem(j) /= -1)
i = i +1
j = j + 2
enddo
i=0
j=1
do while(desc%ovrlap_elem(j) /= -1)
i = i +1
j = j + 2
enddo
call psb_realloc(i,ovrel,info)
if (i > 0) then
allocate(ovrel(i),stat=info)
if (info /= 0 ) then
info = 4000
call psb_errpush(info,name)
goto 9999
end if
i=0
j=1
do while(desc%ovrlap_elem(j) /= -1)
@ -38,10 +39,17 @@ subroutine psb_get_ovrlap(ovrel,desc,info)
enddo
else
info = 1122
call psb_errpush(info,name)
goto 9999
if (allocated(ovrel)) then
deallocate(ovrel,stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Deallocate')
goto 9999
end if
end if
end if
call psb_erractionrestore(err_act)
return

@ -40,7 +40,7 @@
! info - integer. Eventually returns an error code.
! iact - integer(optional). A character defining the behaviour of this subroutine when is found an index not belonging to the calling process
!
subroutine psb_glob_to_loc2(x,y,desc_a,info,iact,owned)
subroutine psb_glob_to_loc2(x,y,desc_a,info,iact)
use psb_descriptor_type
use psb_const_mod
@ -52,63 +52,51 @@ subroutine psb_glob_to_loc2(x,y,desc_a,info,iact,owned)
!...parameters....
type(psb_desc_type), intent(in) :: desc_a
integer, intent(in) :: x(:)
integer, intent(out) :: y(:), info
character, intent(in), optional :: iact
logical, intent(in), optional :: owned
integer, intent(in) :: x(:)
integer, intent(out) :: y(:), info
character, intent(in), optional :: iact
!....locals....
integer :: n, i, tmp
character :: act
integer :: int_err(5), err_act
real(kind(1.d0)) :: real_val
integer, parameter :: zero=0
logical :: owned_
integer :: n, i, tmp
character :: act
integer :: int_err(5), err_act
real(kind(1.d0)) :: real_val
integer, parameter :: zero=0
character(len=20) :: name
integer :: ictxt, iam, np
if(psb_get_errstatus() /= 0) return
info=0
name = 'glob_to_loc'
ictxt = psb_cd_get_context(desc_a)
call psb_info(ictxt,iam,np)
call psb_erractionsave(err_act)
if (present(iact)) then
act=iact
else
act='I'
act='A'
endif
act = toupper(act)
if (present(owned)) then
owned_ = owned
else
owned_ = .false.
end if
int_err=0
real_val = 0.d0
n = size(x)
call psi_idx_cnv(n,x,y,desc_a,info,owned=owned_)
call psi_idx_cnv(n,x,y,desc_a,info)
select case(act)
case('I')
case('E','I')
call psb_erractionrestore(err_act)
return
case('W')
if (count(y(1:n)<0) >0) then
write(0,'("Out of bounds input in subroutine glob_to_loc")')
if ((info /= 0).or.(count(y(1:n)<0) >0)) then
write(0,'("Error ",i5," in subroutine glob_to_loc")') info
end if
case('E','A')
if (count(y(1:n)<0) >0) then
info = 151
case('A')
if ((info /= 0).or.(count(y(1:n)<0) >0)) then
call psb_errpush(info,name)
goto 9999
end if
end select
if (info /= 0) then
call psb_errpush(info,name)
goto 9999
end if
call psb_erractionrestore(err_act)
return
@ -165,7 +153,7 @@ end subroutine psb_glob_to_loc2
! info - integer. Eventually returns an error code.
! iact - integer(optional). A character defining the behaviour of this subroutine when is found an index not belonging to the calling process
!
subroutine psb_glob_to_loc(x,desc_a,info,iact,owned)
subroutine psb_glob_to_loc(x,desc_a,info,iact)
use psb_penv_mod
use psb_descriptor_type
@ -180,58 +168,51 @@ subroutine psb_glob_to_loc(x,desc_a,info,iact,owned)
integer, intent(inout) :: x(:)
integer, intent(out) :: info
character, intent(in), optional :: iact
logical, intent(in), optional :: owned
!....locals....
integer :: n, i, tmp, nk, key, idx, ih, nh, lb, ub, lm
character :: act
integer :: int_err(5), err_act
integer :: int_err(5), err_act, dectype
real(kind(1.d0)) :: real_val, t0, t1,t2
integer, parameter :: zero=0
logical :: owned_
character(len=20) :: name
integer :: ictxt, iam, np
if(psb_get_errstatus() /= 0) return
info=0
name = 'glob_to_loc'
ictxt = psb_cd_get_context(desc_a)
ictxt = desc_a%matrix_data(psb_ctxt_)
call psb_info(ictxt,iam,np)
call psb_erractionsave(err_act)
dectype = desc_a%matrix_data(psb_dec_type_)
if (present(iact)) then
act=iact
else
act='I'
act='A'
endif
act = toupper(act)
if (present(owned)) then
owned_ = owned
else
owned_ = .false.
end if
n = size(x)
call psi_idx_cnv(n,x,desc_a,info,owned=owned_)
call psi_idx_cnv(n,x,desc_a,info)
select case(act)
case('I')
case('E','I')
call psb_erractionrestore(err_act)
return
case('W')
if (count(x(1:n)<0) >0) then
write(0,'("Out of bounds input in subroutine glob_to_loc")')
if ((info /= 0).or.(count(x(1:n)<0) >0)) then
write(0,'("Error ",i5," in subroutine glob_to_loc")') info
end if
case('E','A')
if (count(x(1:n)<0) >0) then
info = 151
case('A')
if ((info /= 0).or.(count(x(1:n)<0) >0)) then
write(0,*) count(x(1:n)<0)
call psb_errpush(info,name)
goto 9999
end if
end select
if (info /= 0) then
call psb_errpush(info,name)
goto 9999
end if
call psb_erractionrestore(err_act)
return
@ -245,69 +226,69 @@ subroutine psb_glob_to_loc(x,desc_a,info,iact,owned)
end if
return
!!$contains
!!$
!!$ subroutine inlbsrch(ipos,key,n,v)
!!$ implicit none
!!$ integer ipos, key, n
!!$ integer v(n)
!!$
!!$ integer lb, ub, m
!!$
!!$
!!$ lb = 1
!!$ ub = n
!!$ ipos = -1
!!$
!!$ do
!!$ if (lb > ub) return
!!$ m = (lb+ub)/2
!!$ if (key.eq.v(m)) then
!!$ ipos = m
!!$ return
!!$ else if (key.lt.v(m)) then
!!$ ub = m-1
!!$ else
!!$ lb = m + 1
!!$ end if
!!$ enddo
!!$ return
!!$ end subroutine inlbsrch
!!$
!!$ subroutine inner_cnv(n,x,hashsize,hashmask,hashv,glb_lc)
!!$ integer :: n, hashsize,hashmask,x(:), hashv(0:),glb_lc(:,:)
!!$ integer :: i, ih, key, idx,nh,tmp,lb,ub,lm
!!$ do i=1, n
!!$ key = x(i)
!!$ ih = iand(key,hashmask)
!!$ idx = hashv(ih)
!!$ nh = hashv(ih+1) - hashv(ih)
!!$ if (nh > 0) then
!!$ tmp = -1
!!$ lb = idx
!!$ ub = idx+nh-1
!!$ do
!!$ if (lb>ub) exit
!!$ lm = (lb+ub)/2
!!$ if (key==glb_lc(lm,1)) then
!!$ tmp = lm
!!$ exit
!!$ else if (key<glb_lc(lm,1)) then
!!$ ub = lm - 1
!!$ else
!!$ lb = lm + 1
!!$ end if
!!$ end do
!!$ else
!!$ tmp = -1
!!$ end if
!!$ if (tmp > 0) then
!!$ x(i) = glb_lc(tmp,2)
!!$ else
!!$ x(i) = tmp
!!$ end if
!!$ end do
!!$ end subroutine inner_cnv
contains
subroutine inlbsrch(ipos,key,n,v)
implicit none
integer ipos, key, n
integer v(n)
integer lb, ub, m
lb = 1
ub = n
ipos = -1
do
if (lb > ub) return
m = (lb+ub)/2
if (key.eq.v(m)) then
ipos = m
return
else if (key.lt.v(m)) then
ub = m-1
else
lb = m + 1
end if
enddo
return
end subroutine inlbsrch
subroutine inner_cnv(n,x,hashsize,hashmask,hashv,glb_lc)
integer :: n, hashsize,hashmask,x(:), hashv(0:),glb_lc(:,:)
integer :: i, ih, key, idx,nh,tmp,lb,ub,lm
do i=1, n
key = x(i)
ih = iand(key,hashmask)
idx = hashv(ih)
nh = hashv(ih+1) - hashv(ih)
if (nh > 0) then
tmp = -1
lb = idx
ub = idx+nh-1
do
if (lb>ub) exit
lm = (lb+ub)/2
if (key==glb_lc(lm,1)) then
tmp = lm
exit
else if (key<glb_lc(lm,1)) then
ub = lm - 1
else
lb = lm + 1
end if
end do
else
tmp = -1
end if
if (tmp > 0) then
x(i) = glb_lc(tmp,2)
else
x(i) = tmp
end if
end do
end subroutine inner_cnv
end subroutine psb_glob_to_loc

@ -73,6 +73,11 @@ subroutine psb_iinsvi(m, irw, val, x, desc_a, info, dupl)
call psb_erractionsave(err_act)
name = 'psb_insvi'
!!$ if (.not.allocated(desc_a%glob_to_loc)) then
!!$ info=3110
!!$ call psb_errpush(info,name)
!!$ return
!!$ end if
if ((.not.allocated(desc_a%matrix_data))) then
int_err(1)=3110
call psb_errpush(info,name)
@ -251,6 +256,11 @@ subroutine psb_iinsi(m,irw, val, x, desc_a, info, dupl)
call psb_erractionsave(err_act)
name = 'psb_iinsi'
!!$ if (.not.allocated(desc_a%glob_to_loc)) then
!!$ info=3110
!!$ call psb_errpush(info,name)
!!$ return
!!$ end if
if ((.not.allocated(desc_a%matrix_data))) then
int_err(1)=3110
call psb_errpush(info,name)

@ -100,12 +100,12 @@ subroutine psb_loc_to_glob2(x,y,desc_a,info,iact)
if (info /= 0) then
select case(act)
case('I')
case('E')
call psb_erractionrestore(err_act)
return
case('W')
write(0,'("Error ",i5," in subroutine glob_to_loc")') info
case('E','A')
case('A')
call psb_errpush(info,name)
goto 9999
end select
@ -223,12 +223,12 @@ subroutine psb_loc_to_glob(x,desc_a,info,iact)
if (info /= 0) then
select case(act)
case('I')
case('E')
call psb_erractionrestore(err_act)
return
case('W')
write(0,'("Error ",i5," in subroutine glob_to_loc")') info
case('A','E')
case('A')
call psb_errpush(info,name)
goto 9999
end select

@ -47,8 +47,6 @@ Subroutine psb_zcdovr(a,desc_a,novr,desc_ov,info)
use psb_serial_mod
use psb_descriptor_type
Use psb_prec_type
Use psb_prec_mod
use psb_error_mod
use psb_penv_mod
use psb_tools_mod
@ -102,11 +100,11 @@ Subroutine psb_zcdovr(a,desc_a,novr,desc_ov,info)
If(debug) Write(0,*)'IN CDOVR1',novr ,m,nnzero,n_col
if (novr<0) then
info=10
int_err(1)=1
int_err(2)=novr
call psb_errpush(info,name,i_err=int_err)
goto 9999
info=10
int_err(1)=1
int_err(2)=novr
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
if (debug) write(0,*) 'Calling desccpy'
@ -139,9 +137,9 @@ Subroutine psb_zcdovr(a,desc_a,novr,desc_ov,info)
!
nztot = psb_sp_get_nnzeros(a)
if (nztot>0) then
lovr = ((nztot+m-1)/m)*nhalo*novr
lworks = ((nztot+m-1)/m)*nhalo
lworkr = ((nztot+m-1)/m)*nhalo
lovr = ((nztot+m-1)/m)*nhalo*novr
lworks = ((nztot+m-1)/m)*nhalo
lworkr = ((nztot+m-1)/m)*nhalo
else
info=-1
call psb_errpush(info,name)
@ -156,164 +154,87 @@ Subroutine psb_zcdovr(a,desc_a,novr,desc_ov,info)
l_tmp_ovr_idx = novr*(3*Max(2*index_dim,1)+1)
l_tmp_halo = novr*(3*Size(desc_a%halo_index))
call psb_cd_set_bld(desc_ov,info)
desc_ov%matrix_data(psb_ovl_state_)=psb_cd_ovl_bld_
if (psb_is_large_desc(desc_a)) then
desc_ov%matrix_data(psb_dec_type_) = psb_desc_large_bld_
else
desc_ov%matrix_data(psb_dec_type_) = psb_desc_bld_
end if
If(debug) then
Write(0,*)'Start cdovrbld',me,lworks,lworkr
call psb_barrier(ictxt)
endif
if (.false.) then
!
! The real work goes on in here....
!
Call psb_cdovrbld(novr,desc_ov,desc_a,a,&
& l_tmp_halo,l_tmp_ovr_idx,lworks,lworkr,info)
if (info /= 0) then
info=4010
ch_err='psb_cdovrbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
Allocate(brvindx(np+1),rvsz(np),sdsz(np),bsdindx(np+1),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
Allocate(works(lworks),workr(lworkr),t_halo_in(l_tmp_halo),&
& t_halo_out(l_tmp_halo), temp(lworkr),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
call psb_sp_all(blk,max(lworks,lworkr),info)
if (info /= 0) then
info=4010
ch_err='psb_sp_all'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
blk%fida='COO'
If(debug) then
Write(0,*)'Done cdovrbld',me,lworks,lworkr
call psb_barrier(ictxt)
endif
Allocate(tmp_ovr_idx(l_tmp_ovr_idx),tmp_halo(l_tmp_halo),&
& halo(size(desc_a%halo_index)),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
halo(:) = desc_a%halo_index(:)
desc_ov%ovrlap_elem(:) = -1
tmp_ovr_idx(:) = -1
tmp_halo(:) = -1
counter_e = 1
tot_recv = 0
counter_h = 1
counter_o = 1
! Init overlap with desc_a%ovrlap (if any)
counter = 1
Do While (desc_a%ovrlap_index(counter) /= -1)
proc = desc_a%ovrlap_index(counter+psb_proc_id_)
n_elem_recv = desc_a%ovrlap_index(counter+psb_n_elem_recv_)
n_elem_send = desc_a%ovrlap_index(counter+n_elem_recv+psb_n_elem_send_)
Do j=0,n_elem_recv-1
idx = desc_a%ovrlap_index(counter+psb_elem_recv_+j)
If(idx > Size(desc_ov%loc_to_glob)) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
gidx = desc_ov%loc_to_glob(idx)
call psb_check_size((counter_o+3),tmp_ovr_idx,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
else
tmp_ovr_idx(counter_o)=proc
tmp_ovr_idx(counter_o+1)=1
tmp_ovr_idx(counter_o+2)=gidx
tmp_ovr_idx(counter_o+3)=-1
counter_o=counter_o+3
end Do
counter=counter+n_elem_recv+n_elem_send+2
end Do
Allocate(brvindx(np+1),rvsz(np),sdsz(np),bsdindx(np+1),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
Allocate(works(lworks),workr(lworkr),t_halo_in(l_tmp_halo),&
& t_halo_out(l_tmp_halo), temp(lworkr),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
!
! A picture is in order to understand what goes on here.
! I is the internal part; H is halo, R row, C column. The final
! matrix with N levels of overlap looks like this
!
! I | Hc1 | 0 | 0 |
! -------|-----|-----|-----|
! Hr1 | Hd1 | Hc2 | 0 |
! -------|-----|-----|-----|
! 0 | Hr2 | Hd2 | Hc2 |
! -------|-----|-----|-----|
! 0 | 0 | Hr3 | Hd3 | Hc3
!
! At the start we already have I and Hc1, so we know the row
! indices that will make up Hr1, and also who owns them. As we
! actually get those rows, we receive the column indices in Hc2;
! these define the row indices for Hr2, and so on. When we have
! reached the desired level HrN, we may ignore HcN.
!
!
Do i_ovr = 1, novr
call psb_sp_all(blk,max(lworks,lworkr),info)
if (info /= 0) then
info=4010
ch_err='psb_sp_all'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (debug) write(0,*) me,'Running on overlap level ',i_ovr,' of ',novr
blk%fida='COO'
!
! At this point, halo contains a valid halo corresponding to the
! matrix enlarged with the elements in the frontier for I_OVR-1.
! At the start, this is just the halo for A; the rows for indices in
! the first halo will contain column indices defining the second halo
! level and so on.
!
bsdindx(:) = 0
sdsz(:) = 0
brvindx(:) = 0
rvsz(:) = 0
idxr = 0
idxs = 0
counter = 1
counter_t = 1
Do While (halo(counter) /= -1)
tot_elem=0
proc=halo(counter+psb_proc_id_)
n_elem_recv=halo(counter+psb_n_elem_recv_)
n_elem_send=halo(counter+n_elem_recv+psb_n_elem_send_)
If ((counter+n_elem_recv+n_elem_send) > Size(halo)) then
info = -1
call psb_errpush(info,name)
goto 9999
end If
tot_recv=tot_recv+n_elem_recv
if (debug) write(0,*) me,' CDOVRBLD tot_recv:',proc,n_elem_recv,tot_recv
!
!
! The format of the halo vector exists in two forms: 1. Temporary
! 2. Assembled. In this loop we are using the (assembled) halo_in and
! copying it into (temporary) halo_out; this is because tmp_halo will
! be enlarged with the new column indices received, and will reassemble
! everything for the next iteration.
!
Allocate(tmp_ovr_idx(l_tmp_ovr_idx),tmp_halo(l_tmp_halo),&
& halo(size(desc_a%halo_index)),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
halo(:) = desc_a%halo_index(:)
desc_ov%ovrlap_elem(:) = -1
tmp_ovr_idx(:) = -1
tmp_halo(:) = -1
counter_e = 1
tot_recv = 0
counter_h = 1
counter_o = 1
! Init overlap with desc_a%ovrlap (if any)
counter = 1
Do While (desc_a%ovrlap_index(counter) /= -1)
proc = desc_a%ovrlap_index(counter+psb_proc_id_)
n_elem_recv = desc_a%ovrlap_index(counter+psb_n_elem_recv_)
n_elem_send = desc_a%ovrlap_index(counter+n_elem_recv+psb_n_elem_send_)
!
! add recv elements in halo_index into ovrlap_index
!
Do j=0,n_elem_recv-1
If((counter+psb_elem_recv_+j)>Size(halo)) then
info=-2
call psb_errpush(info,name)
goto 9999
end If
idx = halo(counter+psb_elem_recv_+j)
idx = desc_a%ovrlap_index(counter+psb_elem_recv_+j)
If(idx > Size(desc_ov%loc_to_glob)) then
info=-3
call psb_errpush(info,name)
@ -334,342 +255,443 @@ Subroutine psb_zcdovr(a,desc_a,novr,desc_ov,info)
tmp_ovr_idx(counter_o+2)=gidx
tmp_ovr_idx(counter_o+3)=-1
counter_o=counter_o+3
if (.not.psb_is_large_desc(desc_ov)) then
call psb_check_size((counter_h+3),tmp_halo,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
end Do
counter=counter+n_elem_recv+n_elem_send+2
end Do
tmp_halo(counter_h)=proc
tmp_halo(counter_h+1)=1
tmp_halo(counter_h+2)=idx
tmp_halo(counter_h+3)=-1
counter_h=counter_h+3
end if
Enddo
if (debug) write(0,*) me,'Checktmp_o_i Loop Mid1',tmp_ovr_idx(1:10)
counter = counter+n_elem_recv
!
! add send elements in halo_index into ovrlap_index
!
Do j=0,n_elem_send-1
!
! A picture is in order to understand what goes on here.
! I is the internal part; H is halo, R row, C column. The final
! matrix with N levels of overlap looks like this
!
! I | Hc1 | 0 | 0 |
! -------|-----|-----|-----|
! Hr1 | Hd1 | Hc2 | 0 |
! -------|-----|-----|-----|
! 0 | Hr2 | Hd2 | Hc2 |
! -------|-----|-----|-----|
! 0 | 0 | Hr3 | Hd3 | Hc3
!
! At the start we already have I and Hc1, so we know the row
! indices that will make up Hr1, and also who owns them. As we
! actually get those rows, we receive the column indices in Hc2;
! these define the row indices for Hr2, and so on. When we have
! reached the desired level HrN, we may ignore HcN.
!
!
Do i_ovr = 1, novr
idx = halo(counter+psb_elem_send_+j)
gidx = desc_ov%loc_to_glob(idx)
if (idx > psb_cd_get_local_rows(Desc_a)) &
& write(0,*) me,i_ovr,'Out of local rows ',&
& idx,psb_cd_get_local_rows(Desc_a)
if (debug) write(0,*) me,'Running on overlap level ',i_ovr,' of ',novr
call psb_check_size((counter_o+3),tmp_ovr_idx,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
!
! At this point, halo contains a valid halo corresponding to the
! matrix enlarged with the elements in the frontier for I_OVR-1.
! At the start, this is just the halo for A; the rows for indices in
! the first halo will contain column indices defining the second halo
! level and so on.
!
bsdindx(:) = 0
sdsz(:) = 0
brvindx(:) = 0
rvsz(:) = 0
idxr = 0
idxs = 0
counter = 1
counter_t = 1
Do While (halo(counter) /= -1)
tot_elem=0
proc=halo(counter+psb_proc_id_)
n_elem_recv=halo(counter+psb_n_elem_recv_)
n_elem_send=halo(counter+n_elem_recv+psb_n_elem_send_)
If ((counter+n_elem_recv+n_elem_send) > Size(halo)) then
info = -1
call psb_errpush(info,name)
goto 9999
end if
tmp_ovr_idx(counter_o)=proc
tmp_ovr_idx(counter_o+1)=1
tmp_ovr_idx(counter_o+2)=gidx
tmp_ovr_idx(counter_o+3)=-1
counter_o=counter_o+3
end If
tot_recv=tot_recv+n_elem_recv
if (debug) write(0,*) me,' CDOVRBLD tot_recv:',proc,n_elem_recv,tot_recv
!
!
! The format of the halo vector exists in two forms: 1. Temporary
! 2. Assembled. In this loop we are using the (assembled) halo_in and
! copying it into (temporary) halo_out; this is because tmp_halo will
! be enlarged with the new column indices received, and will reassemble
! everything for the next iteration.
!
!
! Prepare to exchange the halo rows with the other proc.
! add recv elements in halo_index into ovrlap_index
!
If (i_ovr < (novr)) Then
n_elem = psb_sp_get_nnz_row(idx,a)
Do j=0,n_elem_recv-1
If((counter+psb_elem_recv_+j)>Size(halo)) then
info=-2
call psb_errpush(info,name)
goto 9999
end If
call psb_check_size((idxs+tot_elem+n_elem),works,info)
idx = halo(counter+psb_elem_recv_+j)
If(idx > Size(desc_ov%loc_to_glob)) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
gidx = desc_ov%loc_to_glob(idx)
call psb_check_size((counter_o+3),tmp_ovr_idx,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
If((n_elem) > size(blk%ia2)) Then
isz = max((3*size(blk%ia2))/2,(n_elem))
if (debug) write(0,*) me,'Realloc blk',isz
call psb_sp_reall(blk,isz,info)
tmp_ovr_idx(counter_o)=proc
tmp_ovr_idx(counter_o+1)=1
tmp_ovr_idx(counter_o+2)=gidx
tmp_ovr_idx(counter_o+3)=-1
counter_o=counter_o+3
if (.not.psb_is_large_desc(desc_ov)) then
call psb_check_size((counter_h+3),tmp_halo,info,pad=-1)
if (info /= 0) then
info=4010
ch_err='psb_sp_reall'
call psb_errpush(info,name,a_err=ch_err)
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
End If
call psb_sp_getblk(idx,a,blk,info)
tmp_halo(counter_h)=proc
tmp_halo(counter_h+1)=1
tmp_halo(counter_h+2)=idx
tmp_halo(counter_h+3)=-1
counter_h=counter_h+3
end if
Enddo
if (debug) write(0,*) me,'Checktmp_o_i Loop Mid1',tmp_ovr_idx(1:10)
counter = counter+n_elem_recv
!
! add send elements in halo_index into ovrlap_index
!
Do j=0,n_elem_send-1
idx = halo(counter+psb_elem_send_+j)
gidx = desc_ov%loc_to_glob(idx)
if (idx > psb_cd_get_local_rows(Desc_a)) &
& write(0,*) me,i_ovr,'Out of local rows ',&
& idx,psb_cd_get_local_rows(Desc_a)
call psb_check_size((counter_o+3),tmp_ovr_idx,info,pad=-1)
if (info /= 0) then
info=4010
ch_err='psb_sp_getblk'
call psb_errpush(info,name,a_err=ch_err)
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
!!$ write(0,*) me,'Iteration: ',j,i_ovr
Do jj=1,n_elem
works(idxs+tot_elem+jj)=desc_ov%loc_to_glob(blk%ia2(jj))
End Do
tot_elem=tot_elem+n_elem
End If
Enddo
tmp_ovr_idx(counter_o)=proc
tmp_ovr_idx(counter_o+1)=1
tmp_ovr_idx(counter_o+2)=gidx
tmp_ovr_idx(counter_o+3)=-1
counter_o=counter_o+3
!
! Prepare to exchange the halo rows with the other proc.
!
If (i_ovr < (novr)) Then
n_elem = psb_sp_get_nnz_row(idx,a)
if (i_ovr < novr) then
if (tot_elem > 1) then
call imsr(tot_elem,works(idxs+1))
lx = works(idxs+1)
i = 1
j = 1
do
j = j + 1
if (j > tot_elem) exit
if (works(idxs+j) /= lx) then
i = i + 1
works(idxs+i) = works(idxs+j)
lx = works(idxs+i)
call psb_check_size((idxs+tot_elem+n_elem),works,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
end do
tot_elem = i
endif
if (debug) write(0,*) me,'Checktmp_o_i Loop Mid2',tmp_ovr_idx(1:10)
sdsz(proc+1) = tot_elem
idxs = idxs + tot_elem
end if
counter = counter+n_elem_send+3
if (debug) write(0,*) me,'Checktmp_o_i Loop End',tmp_ovr_idx(1:10)
Enddo
if (debug) write(0,*)me,'End phase 1 CDOVRBLD', m, n_col, tot_recv
if (i_ovr < novr) then
!
! Exchange data requests with everybody else: so far we have
! accumulated RECV requests, we have an all-to-all to build
! matchings SENDs.
!
call mpi_alltoall(sdsz,1,mpi_integer,rvsz,1,mpi_integer,icomm,info)
if (info /= 0) then
info=4010
ch_err='mpi_alltoall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
idxs = 0
idxr = 0
counter = 1
Do
proc=halo(counter)
if (proc == -1) exit
n_elem_recv = halo(counter+psb_n_elem_recv_)
counter = counter+n_elem_recv
n_elem_send = halo(counter+psb_n_elem_send_)
bsdindx(proc+1) = idxs
idxs = idxs + sdsz(proc+1)
brvindx(proc+1) = idxr
idxr = idxr + rvsz(proc+1)
If((n_elem) > size(blk%ia2)) Then
isz = max((3*size(blk%ia2))/2,(n_elem))
if (debug) write(0,*) me,'Realloc blk',isz
call psb_sp_reall(blk,isz,info)
if (info /= 0) then
info=4010
ch_err='psb_sp_reall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
End If
call psb_sp_getblk(idx,a,blk,info)
if (info /= 0) then
info=4010
ch_err='psb_sp_getblk'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
!!$ write(0,*) me,'Iteration: ',j,i_ovr
Do jj=1,n_elem
works(idxs+tot_elem+jj)=desc_ov%loc_to_glob(blk%ia2(jj))
End Do
tot_elem=tot_elem+n_elem
End If
Enddo
if (i_ovr < novr) then
if (tot_elem > 1) then
call imsr(tot_elem,works(idxs+1))
lx = works(idxs+1)
i = 1
j = 1
do
j = j + 1
if (j > tot_elem) exit
if (works(idxs+j) /= lx) then
i = i + 1
works(idxs+i) = works(idxs+j)
lx = works(idxs+i)
end if
end do
tot_elem = i
endif
if (debug) write(0,*) me,'Checktmp_o_i Loop Mid2',tmp_ovr_idx(1:10)
sdsz(proc+1) = tot_elem
idxs = idxs + tot_elem
end if
counter = counter+n_elem_send+3
if (debug) write(0,*) me,'Checktmp_o_i Loop End',tmp_ovr_idx(1:10)
Enddo
iszr=sum(rvsz)
if (max(iszr,1) > lworkr) then
call psb_realloc(max(iszr,1),workr,info)
if (debug) write(0,*)me,'End phase 1 CDOVRBLD', m, n_col, tot_recv
if (i_ovr < novr) then
!
! Exchange data requests with everybody else: so far we have
! accumulated RECV requests, we have an all-to-all to build
! matchings SENDs.
!
call mpi_alltoall(sdsz,1,mpi_integer,rvsz,1,mpi_integer,icomm,info)
if (info /= 0) then
info=4010
ch_err='psb_realloc'
ch_err='mpi_alltoall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
lworkr=max(iszr,1)
end if
call mpi_alltoallv(works,sdsz,bsdindx,mpi_integer,&
& workr,rvsz,brvindx,mpi_integer,icomm,info)
if (info /= 0) then
info=4010
ch_err='mpi_alltoallv'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (debug) write(0,*) 'ISZR :',iszr
idxs = 0
idxr = 0
counter = 1
Do
proc=halo(counter)
if (proc == -1) exit
n_elem_recv = halo(counter+psb_n_elem_recv_)
counter = counter+n_elem_recv
n_elem_send = halo(counter+psb_n_elem_send_)
bsdindx(proc+1) = idxs
idxs = idxs + sdsz(proc+1)
brvindx(proc+1) = idxr
idxr = idxr + rvsz(proc+1)
counter = counter+n_elem_send+3
Enddo
iszr=sum(rvsz)
if (max(iszr,1) > lworkr) then
call psb_realloc(max(iszr,1),workr,info)
if (info /= 0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
lworkr=max(iszr,1)
end if
if (psb_is_large_desc(desc_a)) then
call psb_check_size(iszr,maskr,info)
call mpi_alltoallv(works,sdsz,bsdindx,mpi_integer,&
& workr,rvsz,brvindx,mpi_integer,icomm,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
ch_err='mpi_alltoallv'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psi_idx_cnv(iszr,workr,maskr,desc_ov,info)
iszs = count(maskr(1:iszr)<=0)
if (iszs > size(works)) call psb_realloc(iszs,works,info)
j = 0
do i=1,iszr
if (maskr(i) < 0) then
j=j+1
works(j) = workr(i)
end if
end do
!
! fnd_owner on desc_a because we want the procs who
! owned the rows from the beginning!
!
call psi_fnd_owner(iszs,works,temp,desc_a,info)
n_col=psb_cd_get_local_cols(desc_ov)
do i=1,iszs
idx = works(i)
n_col = psb_cd_get_local_cols(desc_ov)
call psi_idx_ins_cnv(idx,lidx,desc_ov,info)
if (psb_cd_get_local_cols(desc_ov) > n_col ) then
!
! This is a new index. Assigning a local index as
! we receive them guarantees that all indices for HALO(I)
! will be less than those for HALO(J) whenever I<J
!
proc_id = temp(i)
call psb_check_size((counter_t+3),t_halo_in,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
t_halo_in(counter_t)=proc_id
t_halo_in(counter_t+1)=1
t_halo_in(counter_t+2)=lidx
counter_t=counter_t+3
if (.false.) write(0,*) me,' CDOVRBLD: Added t_halo_in ',&
&proc_id,lidx,idx
endif
end Do
else
Do i=1,iszr
idx=workr(i)
if (idx <1) then
write(0,*) me,'Error in CDOVRBLD level',i_ovr,' : ',idx,i,iszr
else If (desc_ov%glob_to_loc(idx) < -np) Then
!
! This is a new index. Assigning a local index as
! we receive them guarantees that all indices for HALO(I)
! will be less than those for HALO(J) whenever I<J
!
n_col=n_col+1
proc_id=-desc_ov%glob_to_loc(idx)-np-1
call psb_check_size(n_col,desc_ov%loc_to_glob,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
desc_ov%glob_to_loc(idx)=n_col
desc_ov%loc_to_glob(n_col)=idx
if (debug) write(0,*) 'ISZR :',iszr
call psb_check_size((counter_t+3),t_halo_in,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
if (psb_is_large_desc(desc_a)) then
call psb_check_size(iszr,maskr,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
call psi_idx_cnv(iszr,workr,maskr,desc_ov,info)
iszs = count(maskr(1:iszr)<=0)
if (iszs > size(works)) call psb_realloc(iszs,works,info)
j = 0
do i=1,iszr
if (maskr(i) < 0) then
j=j+1
works(j) = workr(i)
end if
end do
!
! fnd_owner on desc_a because we want the procs who
! owned the rows from the beginning!
!
call psi_fnd_owner(iszs,works,temp,desc_a,info)
n_col=psb_cd_get_local_cols(desc_ov)
do i=1,iszs
idx = works(i)
n_col = psb_cd_get_local_cols(desc_ov)
call psi_idx_ins_cnv(idx,lidx,desc_ov,info)
if (psb_cd_get_local_cols(desc_ov) > n_col ) then
!
! This is a new index. Assigning a local index as
! we receive them guarantees that all indices for HALO(I)
! will be less than those for HALO(J) whenever I<J
!
proc_id = temp(i)
call psb_check_size((counter_t+3),t_halo_in,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
t_halo_in(counter_t)=proc_id
t_halo_in(counter_t+1)=1
t_halo_in(counter_t+2)=lidx
counter_t=counter_t+3
if (.false.) write(0,*) me,' CDOVRBLD: Added t_halo_in ',&
&proc_id,lidx,idx
endif
end Do
else
Do i=1,iszr
idx=workr(i)
if (idx <1) then
write(0,*) me,'Error in CDOVRBLD level',i_ovr,' : ',idx,i,iszr
else If (desc_ov%glob_to_loc(idx) < -np) Then
!
! This is a new index. Assigning a local index as
! we receive them guarantees that all indices for HALO(I)
! will be less than those for HALO(J) whenever I<J
!
n_col=n_col+1
proc_id=-desc_ov%glob_to_loc(idx)-np-1
call psb_check_size(n_col,desc_ov%loc_to_glob,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
desc_ov%glob_to_loc(idx)=n_col
desc_ov%loc_to_glob(n_col)=idx
call psb_check_size((counter_t+3),t_halo_in,info,pad=-1)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_check_size')
goto 9999
end if
t_halo_in(counter_t)=proc_id
t_halo_in(counter_t+1)=1
t_halo_in(counter_t+2)=n_col
counter_t=counter_t+3
if (debug) write(0,*) me,' CDOVRBLD: Added into t_halo_in from recv',&
&proc_id,n_col,idx
else if (desc_ov%glob_to_loc(idx) < 0) Then
if (debug) write(0,*) me,'Wrong input to cdovrbld??',&
&idx,desc_ov%glob_to_loc(idx)
End If
End Do
desc_ov%matrix_data(psb_n_col_)=n_col
t_halo_in(counter_t)=proc_id
t_halo_in(counter_t+1)=1
t_halo_in(counter_t+2)=n_col
counter_t=counter_t+3
if (debug) write(0,*) me,' CDOVRBLD: Added into t_halo_in from recv',&
&proc_id,n_col,idx
else if (desc_ov%glob_to_loc(idx) < 0) Then
if (debug) write(0,*) me,'Wrong input to cdovrbld??',&
&idx,desc_ov%glob_to_loc(idx)
End If
End Do
desc_ov%matrix_data(psb_n_col_)=n_col
end if
end if
end if
!!$ desc_ov%matrix_data(psb_n_row_)=desc_ov%matrix_data(psb_n_col_)
!
! Ok, now we have a temporary halo with all the info for the
! next round. If we need to keep going, convert the halo format
! from temporary to final, so that we can work out the next iteration.
! This uses one of the convert_comm internals, i.e. we are doing
! the equivalent of a partial call to convert_comm
!
!
! Ok, now we have a temporary halo with all the info for the
! next round. If we need to keep going, convert the halo format
! from temporary to final, so that we can work out the next iteration.
! This uses one of the convert_comm internals, i.e. we are doing
! the equivalent of a partial call to convert_comm
!
If (i_ovr < (novr)) Then
If (i_ovr < (novr)) Then
t_halo_in(counter_t)=-1
t_halo_in(counter_t)=-1
if (debug) write(0,*) me,'Checktmp_o_i 1',tmp_ovr_idx(1:10)
if (debug) write(0,*) me,'Calling Crea_Halo'
if (debug) write(0,*) me,'Checktmp_o_i 1',tmp_ovr_idx(1:10)
if (debug) write(0,*) me,'Calling Crea_Halo'
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,.false.,&
& nxch,nsnd,nrcv,info)
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,.false.,&
& nxch,nsnd,nrcv,info)
if (debug) then
write(0,*) me,'Done Crea_Index'
call psb_barrier(ictxt)
end if
if (debug) write(0,*) me,'Checktmp_o_i 2',tmp_ovr_idx(1:10)
if (debug) write(0,*) me,'Done Crea_Halo'
call psb_transfer(t_halo_out,halo,info)
!
! At this point we have built the halo necessary for I_OVR+1.
!
End If
if (debug) write(0,*) me,'Checktmp_o_i ',tmp_ovr_idx(1:10)
if (debug) then
write(0,*) me,'Done Crea_Index'
call psb_barrier(ictxt)
end if
if (debug) write(0,*) me,'Checktmp_o_i 2',tmp_ovr_idx(1:10)
if (debug) write(0,*) me,'Done Crea_Halo'
call psb_transfer(t_halo_out,halo,info)
!
! At this point we have built the halo necessary for I_OVR+1.
!
End If
if (debug) write(0,*) me,'Checktmp_o_i ',tmp_ovr_idx(1:10)
End Do
End Do
desc_ov%matrix_data(psb_m_)=psb_cd_get_global_rows(desc_a)
desc_ov%matrix_data(psb_n_)=psb_cd_get_global_cols(desc_a)
desc_ov%matrix_data(psb_m_)=psb_cd_get_global_rows(desc_a)
desc_ov%matrix_data(psb_n_)=psb_cd_get_global_cols(desc_a)
tmp_halo(counter_h:)=-1
tmp_ovr_idx(counter_o:)=-1
tmp_halo(counter_h:)=-1
tmp_ovr_idx(counter_o:)=-1
!
! At this point we have gathered all the indices in the halo at
! N levels of overlap. Just call cnv_dsc. This is
! the same routine as gets called inside CDASB.
!
!
! At this point we have gathered all the indices in the halo at
! N levels of overlap. Just call cnv_dsc. This is
! the same routine as gets called inside CDASB.
!
if (debug) then
write(0,*) 'psb_cdovrbld: converting indexes'
call psb_barrier(ictxt)
end if
!.... convert comunication stuctures....
! Note that we have to keep local_rows until the very end,
! because otherwise the halo build mechanism of cdasb
! will not work.
call psb_transfer(tmp_ovr_idx,desc_ov%ovrlap_index,info)
call psb_transfer(tmp_halo,desc_ov%halo_index,info)
call psb_cdasb(desc_ov,info)
desc_ov%matrix_data(psb_n_row_)=desc_ov%matrix_data(psb_n_col_)
if (debug) then
write(0,*) me,'Done CDASB'
call psb_barrier(ictxt)
end if
if (debug) then
write(0,*) 'psb_cdovrbld: converting indexes'
call psb_barrier(ictxt)
end if
!.... convert comunication stuctures....
! Note that we have to keep local_rows until the very end,
! because otherwise the halo build mechanism of cdasb
! will not work.
call psb_transfer(tmp_ovr_idx,desc_ov%ovrlap_index,info)
call psb_transfer(tmp_halo,desc_ov%halo_index,info)
call psb_cdasb(desc_ov,info)
desc_ov%matrix_data(psb_n_row_)=desc_ov%matrix_data(psb_n_col_)
if (debug) then
write(0,*) me,'Done CDASB'
call psb_barrier(ictxt)
end if
if (info == 0) call psb_sp_free(blk,info)
if (info /= 0) then
ch_err='sp_free'
call psb_errpush(4013,name,a_err=ch_err,i_err=(/info,0,0,0,0/))
goto 9999
end if
if (info == 0) call psb_sp_free(blk,info)
if (info /= 0) then
ch_err='sp_free'
call psb_errpush(4013,name,a_err=ch_err,i_err=(/info,0,0,0,0/))
goto 9999
end if
end if
call psb_erractionrestore(err_act)
return

@ -45,9 +45,29 @@ subroutine psb_zcsrp(trans,iperm,a, desc_a, info)
use psb_serial_mod
use psb_const_mod
use psb_penv_mod
use psb_serial_mod
implicit none
! implicit none
interface
subroutine zcsrp(trans,m,n,fida,descra,ia1,ia2,&
& infoa,p,work,lwork,ierror)
integer, intent(in) :: m, n, lwork
integer, intent(out) :: ierror
character, intent(in) :: trans
complex(kind(1.d0)), intent(inout) :: work(*)
integer, intent(in) :: p(*)
integer, intent(inout) :: ia1(*), ia2(*), infoa(*)
character, intent(in) :: fida*5, descra*11
end subroutine zcsrp
end interface
interface isaperm
logical function isaperm(n,ip)
integer, intent(in) :: n
integer, intent(inout) :: ip(*)
end function isaperm
end interface
!...parameters....
type(psb_zspmat_type), intent(inout) :: a
@ -138,7 +158,7 @@ subroutine psb_zcsrp(trans,iperm,a, desc_a, info)
! hmm, maybe we should just move all of this onto a different level,
! have a specialized subroutine, and do it in the solver context????
if (debug) write(0,*) 'spasb: calling dcsrp',size(work_dcsdp)
call csrp(trans,n_row,n_col,a%fida,a%descra,a%ia1,a%ia2,a%infoa,&
call zcsrp(trans,n_row,n_col,a%fida,a%descra,a%ia1,a%ia2,a%infoa,&
& ipt,work_dcsdp,size(work_dcsdp),info)
if(info /= no_err) then
info=4010

@ -71,6 +71,14 @@ subroutine psb_zgelp(trans,iperm,x,desc_a,info)
end subroutine zgelp
end interface
interface isaperm
logical function isaperm(n,ip)
integer, intent(in) :: n
integer, intent(inout) :: ip(*)
end function isaperm
end interface
character(len=20) :: name, ch_err
name = 'psb_zgelp'
@ -205,8 +213,16 @@ subroutine psb_zgelpv(trans,iperm,x,desc_a,info)
integer, intent(in) :: p(*)
end subroutine zgelp
end interface
character(len=20) :: name, ch_err
interface isaperm
logical function isaperm(n,ip)
integer, intent(in) :: n
integer, intent(inout) :: ip(*)
end function isaperm
end interface
character(len=20) :: name, ch_err
name = 'psb_zgelpv'
if(psb_get_errstatus() /= 0) return

@ -74,6 +74,11 @@ subroutine psb_zinsvi(m, irw, val, x, desc_a, info, dupl)
call psb_erractionsave(err_act)
name = 'psb_zinsvi'
!!$ if (.not.allocated(desc_a%glob_to_loc)) then
!!$ info=3110
!!$ call psb_errpush(info,name)
!!$ return
!!$ end if
if ((.not.allocated(desc_a%matrix_data))) then
int_err(1)=3110
call psb_errpush(info,name)
@ -252,6 +257,11 @@ subroutine psb_zinsi(m,irw, val, x, desc_a, info, dupl)
call psb_erractionsave(err_act)
name = 'psb_zinsi'
!!$ if (.not.allocated(desc_a%glob_to_loc)) then
!!$ info=3110
!!$ call psb_errpush(info,name)
!!$ return
!!$ end if
if ((.not.allocated(desc_a%matrix_data))) then
int_err(1)=3110
call psb_errpush(info,name)

@ -1,44 +1,44 @@
!!$
!!$ Parallel Sparse BLAS v2.0
!!$ (C) Copyright 2006 Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ 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_zspcnv.f90
!
! Subroutine: psb_zspcnv
! converts sparse matrix a into b
!
! Parameters:
! a - type(<psb_zspmat_type>). The sparse input matrix.
! b - type(<psb_zspmat_type>). The sparse output matrix.
! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code.
!
!!$
!!$ Parallel Sparse BLAS v2.0
!!$ (C) Copyright 2006 Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ 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_zspcnv.f90
!
! Subroutine: psb_zspcnv
! converts sparse matrix a into b
!
! Parameters:
! a - type(<psb_zspmat_type>). The sparse input matrix.
! b - type(<psb_zspmat_type>). The sparse output matrix.
! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code.
!
subroutine psb_zspcnv(a,b,desc_a,info)
use psb_descriptor_type
@ -50,6 +50,49 @@ subroutine psb_zspcnv(a,b,desc_a,info)
use psb_penv_mod
implicit none
interface zcsdp
subroutine zcsdp(check,trans,m,n,unitd,d,&
& fida,descra,a,ia1,ia2,infoa,&
& pl,fidh,descrh,h,ih1,ih2,infoh,pr,lh,lh1,lh2,&
& work,lwork,ierror)
integer, intent(in) :: lh, lwork, lh1, lh2, m, n
integer, intent(out) :: ierror
character, intent(in) :: check, trans, unitd
complex(kind(1.d0)), intent(in) :: d(*), a(*)
complex(kind(1.d0)), intent(out) :: h(*)
complex(kind(1.d0)), intent(inout) :: work(*)
integer, intent(in) :: ia1(*), ia2(*), infoa(*)
integer, intent(out) :: ih1(*), ih2(*), pl(*),pr(*), infoh(*)
character, intent(in) :: fida*5, descra*11
character, intent(out) :: fidh*5, descrh*11
end subroutine zcsdp
end interface
interface zcsrp
subroutine zcsrp(trans,m,n,fida,descra,ia1,ia2,&
& infoa,p,work,lwork,ierror)
integer, intent(in) :: m, n, lwork
integer, intent(out) :: ierror
character, intent(in) :: trans
complex(kind(1.d0)), intent(inout) :: work(*)
integer, intent(in) :: p(*)
integer, intent(inout) :: ia1(*), ia2(*), infoa(*)
character, intent(in) :: fida*5, descra*11
end subroutine zcsrp
end interface
interface zcsprt
subroutine zcsprt(m,n,fida,descra,a,ia1,ia2,infoa ,iout,ierror)
integer, intent(in) :: iout,m, n
integer, intent(out) :: ierror
complex(kind(1.d0)), intent(in) :: a(*)
integer, intent(in) :: ia1(*), ia2(*), infoa(*)
character, intent(in) :: fida*5, descra*11
end subroutine zcsprt
end interface
!...parameters....
type(psb_zspmat_type), intent(in) :: a
@ -58,11 +101,17 @@ subroutine psb_zspcnv(a,b,desc_a,info)
integer, intent(out) :: info
!....locals....
integer :: int_err(5)
complex(kind(1.d0)) :: d(1)
integer,allocatable :: i_temp(:)
complex(kind(1.d0)),allocatable :: work_dcsdp(:)
integer :: ia1_size,ia2_size,aspk_size,err_act&
& ,i,err,np,me,n_col
integer, allocatable :: i_temp(:)
integer :: dectype
& ,i,err,np,me,n_col,l_dcsdp
integer :: lwork_dcsdp,dectype
integer :: ictxt,n_row
character :: check*1, trans*1, unitd*1
real(kind(1.d0)) :: time(10), mpi_wtime
external mpi_wtime
logical, parameter :: debug=.false.
character(len=20) :: name, ch_err
@ -71,6 +120,7 @@ subroutine psb_zspcnv(a,b,desc_a,info)
name = 'psb_zspcnv'
call psb_erractionsave(err_act)
time(1) = mpi_wtime()
ictxt = psb_cd_get_context(desc_a)
dectype = psb_cd_get_dectype(desc_a)
@ -100,15 +150,46 @@ subroutine psb_zspcnv(a,b,desc_a,info)
if (debug) write (0, *) name,' sizes',ia1_size,ia2_size,aspk_size
! convert only without check
check='N'
trans='N'
unitd='U'
! l_dcsdp is the size requested for dcsdp procedure
l_dcsdp=(ia1_size+100)
b%m=n_row
b%k=n_col
call psb_sp_all(b,ia1_size,ia2_size,aspk_size,info)
allocate(work_dcsdp(l_dcsdp),stat=info)
if (info /= 0) then
info=2025
int_err(1)=l_dcsdp
call psb_errpush(info, name, i_err=int_err)
goto 9999
endif
call psb_csdp(a,b,info)
lwork_dcsdp=size(work_dcsdp)
! set infoa(1) to nnzero
b%pl(:) = 0
b%pr(:) = 0
if (debug) write (0, *) name,' calling dcsdp',lwork_dcsdp,&
&size(work_dcsdp)
! convert aspk,ia1,ia2 in requested representation mode
if (debug) then
endif
! result is put in b
call zcsdp(check,trans,n_row,n_col,unitd,d,a%fida,a%descra,&
& a%aspk,a%ia1,a%ia2,a%infoa,&
& b%pl,b%fida,b%descra,b%aspk,b%ia1,b%ia2,b%infoa,b%pr,&
& size(b%aspk),size(b%ia1),size(b%ia2),&
& work_dcsdp,size(work_dcsdp),info)
if(info /= no_err) then
info=4010
ch_err='psb_csdp'
ch_err='zcsdp'
call psb_errpush(info, name, a_err=ch_err)
goto 9999
end if
@ -148,6 +229,9 @@ subroutine psb_zspcnv(a,b,desc_a,info)
endif
if (debug) write (0, *) me,name,' from zcsdp ',&
&b%fida,' pl ', b%pl(:),'pr',b%pr(:)
call psb_erractionrestore(err_act)
return

@ -49,7 +49,6 @@ Subroutine psb_zsphalo(a,desc_a,blk,info,rwcnv,clcnv,outfmt)
use psb_serial_mod
use psb_descriptor_type
Use psb_prec_type
use psb_realloc_mod
use psb_tools_mod, only : psb_glob_to_loc, psb_loc_to_glob
use psb_error_mod

@ -71,6 +71,27 @@ subroutine psb_zspins(nz,ia,ja,val,a,desc_a,info,rebuild)
logical :: rebuild_
integer, allocatable :: ila(:),jla(:)
!!$ interface psb_cdins
!!$ subroutine psb_cdins(nz,ia,ja,desc_a,info,ila,jla)
!!$ use psb_descriptor_type
!!$ implicit none
!!$ type(psb_desc_type), intent(inout) :: desc_a
!!$ integer, intent(in) :: nz,ia(:),ja(:)
!!$ integer, intent(out) :: info
!!$ integer, optional, intent(out) :: ila(:), jla(:)
!!$ end subroutine psb_cdins
!!$ end interface
!!$
!!$ interface psb_glob_to_loc
!!$ subroutine psb_glob_to_loc(x,desc_a,info,iact)
!!$ use psb_descriptor_type
!!$ implicit none
!!$ type(psb_desc_type), intent(in) :: desc_a
!!$ integer, intent(inout) :: x(:)
!!$ integer, intent(out) :: info
!!$ character, intent(in), optional :: iact
!!$ end subroutine psb_glob_to_loc
!!$ end interface
character(len=20) :: name, ch_err
info = 0

@ -1,24 +1,28 @@
include ../../Make.inc
include ../Make.inc
LIBDIR=../../lib
HERE=.
OBJS= psb_dcgstab.o psb_dcg.o psb_dcgs.o \
LIBDIR=../lib
OBJS=psb_krylov_mod.o psb_dcgstab.o psb_dcg.o psb_dcgs.o \
psb_dbicg.o psb_dcgstabl.o psb_dgmresr.o\
psb_zcgstab.o psb_zcgs.o
INCDIRS=-I. -I.. -I$(LIBDIR)
lib: $(OBJS)
$(AR) $(LIBDIR)/$(LIBNAME) $(OBJS)
$(RANLIB) $(LIBDIR)/$(LIBNAME)
LIBMOD=psb_krylov_mod$(.mod)
LOCAL_MODS=$(LIBMOD)
LIBNAME=$(METHDLIBNAME)
#$(F90OBJS): $(MODS)
INCDIRS=-I. -I$(LIBDIR)
lib: $(OBJS)
$(AR) $(HERE)/$(LIBNAME) $(OBJS)
$(RANLIB) $(HERE)/$(LIBNAME)
/bin/cp $(HERE)/$(LIBNAME) $(LIBDIR)
/bin/cp $(LIBMOD) $(LIBDIR)
veryclean: clean
/bin/rm -f $(LIBNAME)
/bin/rm -f $(HERE)/$(LIBNAME)
clean:
/bin/rm -f $(OBJS) $(LOCAL_MODS)

@ -75,15 +75,8 @@
!
subroutine psb_dbicg(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err, itrace,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_psblas_mod
use psb_tools_mod
use psb_const_mod
use psb_base_mod
use psb_prec_mod
use psb_error_mod
use psb_penv_mod
implicit none
!!$ parameters

@ -75,15 +75,8 @@
!
Subroutine psb_dcg(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err, itrace, istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_psblas_mod
use psb_tools_mod
use psb_const_mod
use psb_base_mod
use psb_prec_mod
use psb_error_mod
use psb_penv_mod
implicit none
!!$ Parameters

@ -73,15 +73,8 @@
!
Subroutine psb_dcgs(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err,itrace,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_psblas_mod
use psb_tools_mod
use psb_const_mod
use psb_base_mod
use psb_prec_mod
use psb_error_mod
use psb_penv_mod
implicit none
!!$ parameters

@ -74,15 +74,8 @@
!
Subroutine psb_dcgstab(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err,itrace, istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_psblas_mod
use psb_tools_mod
use psb_const_mod
use psb_base_mod
use psb_prec_mod
use psb_error_mod
use psb_penv_mod
Implicit None
!!$ parameters
Type(psb_dspmat_type), Intent(in) :: a

@ -80,15 +80,8 @@
!
Subroutine psb_dcgstabl(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err,itrace,irst,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_psblas_mod
use psb_tools_mod
use psb_const_mod
use psb_base_mod
use psb_prec_mod
use psb_error_mod
use psb_penv_mod
implicit none
!!$ parameters

@ -82,15 +82,8 @@
!
Subroutine psb_dgmresr(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err,itrace,irst,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_psblas_mod
use psb_tools_mod
use psb_const_mod
use psb_base_mod
use psb_prec_mod
use psb_error_mod
use psb_penv_mod
implicit none
!!$ Parameters

@ -28,20 +28,21 @@
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
Module psb_methd_mod
Module psb_krylov_mod
use psb_base_mod
use psb_prec_mod
interface psb_krylov
module procedure psb_dkrylov, psb_zkrylov
end interface
interface psb_cg
subroutine psb_dcg(a,prec,b,x,eps,&
& desc_a,info,itmax,iter,err,itrace,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_base_mod
use psb_prec_mod
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(kind(1.d0)), intent(in) :: b(:)
@ -58,9 +59,8 @@ Module psb_methd_mod
interface psb_bicg
subroutine psb_dbicg(a,prec,b,x,eps,&
& desc_a,info,itmax,iter,err,itrace,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_base_mod
use psb_prec_mod
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(kind(1.d0)), intent(in) :: b(:)
@ -77,9 +77,8 @@ Module psb_methd_mod
interface psb_bicgstab
subroutine psb_dcgstab(a,prec,b,x,eps,&
& desc_a,info,itmax,iter,err,itrace,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_base_mod
use psb_prec_mod
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(kind(1.d0)), intent(in) :: b(:)
@ -93,9 +92,8 @@ Module psb_methd_mod
end subroutine psb_dcgstab
subroutine psb_zcgstab(a,prec,b,x,eps,&
& desc_a,info,itmax,iter,err,itrace,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_base_mod
use psb_prec_mod
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
complex(kind(1.d0)), intent(in) :: b(:)
@ -112,9 +110,8 @@ Module psb_methd_mod
interface psb_bicgstabl
Subroutine psb_dcgstabl(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err, itrace,irst,istop)
use psb_serial_mod
use psb_descriptor_type
Use psb_prec_type
use psb_base_mod
use psb_prec_mod
Type(psb_dspmat_type), Intent(in) :: a
Type(psb_desc_type), Intent(in) :: desc_a
type(psb_dprec_type), intent(in) :: prec
@ -131,9 +128,8 @@ Module psb_methd_mod
interface psb_rgmres
Subroutine psb_dgmresr(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err,itrace,irst,istop)
use psb_serial_mod
use psb_descriptor_type
Use psb_prec_type
use psb_base_mod
use psb_prec_mod
!!$ parameters
Type(psb_dspmat_type), Intent(in) :: a
Type(psb_desc_type), Intent(in) :: desc_a
@ -151,9 +147,8 @@ Module psb_methd_mod
interface psb_cgs
subroutine psb_dcgs(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err,itrace,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_base_mod
use psb_prec_mod
!!$ parameters
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
@ -168,9 +163,8 @@ Module psb_methd_mod
end subroutine psb_dcgs
subroutine psb_zcgs(a,prec,b,x,eps,&
& desc_a,info,itmax,iter,err,itrace,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_base_mod
use psb_prec_mod
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
complex(kind(1.d0)), intent(in) :: b(:)
@ -189,12 +183,6 @@ contains
Subroutine psb_dkrylov(method,a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err,itrace,irst,istop)
use psb_serial_mod
use psb_descriptor_type
Use psb_prec_type
use psb_string_mod
use psb_penv_mod
!!$ parameters
character(len=*) :: method
Type(psb_dspmat_type), Intent(in) :: a
Type(psb_desc_type), Intent(in) :: desc_a
@ -207,9 +195,8 @@ contains
Integer, Optional, Intent(out) :: iter
Real(Kind(1.d0)), Optional, Intent(out) :: err
integer :: ictxt, me, np
integer :: itmax_, itrace_, irst_, istop_, iter_
real(kind(1.d0)) :: err_
integer :: itmax_, itrace_, irst_, istop_, iter_
real(kind(1.d0)) :: err_
if (present(itmax)) then
itmax_ = itmax
@ -235,8 +222,6 @@ contains
istop_ = 1
end if
ictxt = psb_cd_get_context(desc_a)
call psb_info(ictxt,me,np)
select case(toupper(method))
case('CG')
@ -258,8 +243,6 @@ contains
call psb_bicgstabl(a,prec,b,x,eps,desc_a,info,&
&itmax_,iter_,err_,itrace_,irst_,istop_)
case default
if (me==0) write(0,*) &
& 'psb_krylov: unknown method, defaulting to BiCGSTAB'
call psb_bicgstab(a,prec,b,x,eps,desc_a,info,&
&itmax_,iter_,err_,itrace_,istop_)
end select
@ -277,12 +260,6 @@ contains
Subroutine psb_zkrylov(method,a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err,itrace,irst,istop)
use psb_serial_mod
use psb_descriptor_type
Use psb_prec_type
use psb_string_mod
use psb_penv_mod
!!$ parameters
character(len=*) :: method
Type(psb_zspmat_type), Intent(in) :: a
Type(psb_desc_type), Intent(in) :: desc_a
@ -295,8 +272,6 @@ contains
Integer, Optional, Intent(out) :: iter
Real(Kind(1.d0)), Optional, Intent(out) :: err
integer :: ictxt, me, np
integer :: itmax_, itrace_, irst_, istop_, iter_
real(kind(1.d0)) :: err_
@ -325,9 +300,6 @@ contains
end if
ictxt = psb_cd_get_context(desc_a)
call psb_info(ictxt,me,np)
select case(toupper(method))
!!$ case('CG')
!!$ call psb_cg(a,prec,b,x,eps,desc_a,info,&
@ -348,8 +320,6 @@ contains
!!$ call psb_bicgstabl(a,prec,b,x,eps,desc_a,info,&
!!$ &itmax_,iter_,err_,itrace_,irst_,istop_)
case default
if (me==0) write(0,*) &
& 'psb_krylov: unknown method, defaulting to BiCGSTAB'
call psb_bicgstab(a,prec,b,x,eps,desc_a,info,&
&itmax_,iter_,err_,itrace_,istop_)
end select
@ -366,7 +336,7 @@ contains
end module psb_methd_mod
end module psb_krylov_mod

@ -73,15 +73,8 @@
!
Subroutine psb_zcgs(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err,itrace,istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_psblas_mod
use psb_tools_mod
use psb_const_mod
use psb_base_mod
use psb_prec_mod
use psb_error_mod
use psb_penv_mod
implicit none
!!$ parameters

@ -60,7 +60,8 @@
!
! Parameters:
! a - type(<psb_zspmat_type>). The sparse matrix containing A.
! prec - type(<psb_prec_type>). The data structure containing the preconditioner.
! prec - type(<psb_prec_type>). The data structure containing the
! preconditioner.
! b - real,dimension(:). The right hand side.
! x - real,dimension(:). The vector of unknowns.
! eps - real. The error tolerance.
@ -74,15 +75,8 @@
!
Subroutine psb_zcgstab(a,prec,b,x,eps,desc_a,info,&
&itmax,iter,err,itrace, istop)
use psb_serial_mod
use psb_descriptor_type
use psb_prec_type
use psb_psblas_mod
use psb_tools_mod
use psb_const_mod
use psb_base_mod
use psb_prec_mod
use psb_error_mod
use psb_penv_mod
Implicit None
!!$ parameters
Type(psb_zspmat_type), Intent(in) :: a

@ -0,0 +1,33 @@
MD2P4
Multilevel Domain Decomposition Parallel Preconditioner Package for PSBLAS
for
Parallel Sparse BLAS v2.0
(C) Copyright 2006 Salvatore Filippone University of Rome Tor Vergata
Alfredo Buttari University of Rome Tor Vergata
Daniela di Serafino Second University of Naples
Pasqua D'Ambra ICAR-CNR
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 MD2P4 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 MD2P4 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.

@ -1,31 +1,40 @@
include ../../Make.inc
include ../Make.inc
LIBDIR=../../lib/
LIBDIR=../lib
HERE=.
INCDIRS=-I. -I$(LIBDIR)
MPFOBJS=psb_dilu_bld.o psb_dbldaggrmat.o psb_zilu_bld.o psb_zbldaggrmat.o
MODOBJS= psb_prec_type.o psb_prec_mod.o
MPFOBJS=psb_dbldaggrmat.o psb_zbldaggrmat.o
F90OBJS=psb_dasmatbld.o psb_dslu_bld.o psb_dumf_bld.o psb_dilu_fct.o\
psb_dmlprc_bld.o psb_dsp_renum.o\
psb_dprecbld.o psb_dprecfree.o psb_dprecset.o \
psb_dmlprc_bld.o psb_dsp_renum.o psb_dilu_bld.o \
psb_dprecbld.o psb_dprecfree.o psb_dprecset.o \
psb_dbaseprc_bld.o psb_ddiagsc_bld.o psb_dgenaggrmap.o \
psb_dprc_aply.o psb_dmlprc_aply.o \
psb_dbaseprc_aply.o psb_dbjac_aply.o\
psb_zasmatbld.o psb_zslu_bld.o psb_zumf_bld.o psb_zilu_fct.o\
psb_zmlprc_bld.o psb_zsp_renum.o\
psb_zprecbld.o psb_zprecfree.o psb_zprecset.o \
psb_zmlprc_bld.o psb_zsp_renum.o psb_zilu_bld.o \
psb_zprecbld.o psb_zprecfree.o psb_zprecset.o \
psb_zbaseprc_bld.o psb_zdiagsc_bld.o psb_zgenaggrmap.o \
psb_zprc_aply.o psb_zmlprc_aply.o \
psb_zbaseprc_aply.o psb_zbjac_aply.o\
$(MPFOBJS)
COBJS=psb_slu_impl.o psb_umf_impl.o psb_zslu_impl.o psb_zumf_impl.o
INCDIRS=-I. -I.. -I$(LIBDIR)
OBJS=$(F90OBJS) $(COBJS) $(MPFOBJS) $(MODOBJS)
LIBMOD=psb_prec_mod$(.mod)
LOCAL_MODS=$(LIBMOD) psb_prec_type$(.mod)
LIBNAME=$(PRECLIBNAME)
OBJS=$(F90OBJS) $(COBJS) $(MPFOBJS)
lib: mpobjs $(OBJS)
$(AR) $(HERE)/$(LIBNAME) $(OBJS)
$(RANLIB) $(HERE)/$(LIBNAME)
/bin/cp $(HERE)/$(LIBNAME) $(LIBDIR)
/bin/cp $(LIBMOD) $(LIBDIR)
lib: mpobjs $(OBJS)
$(AR) $(LIBDIR)/$(LIBNAME) $(MPFOBJS) $(OBJS)
$(RANLIB) $(LIBDIR)/$(LIBNAME)
$(F90OBJS) $(MPFOBJS): $(MODOBJS)
psb_prec_mod.o: psb_prec_type.o
mpobjs:
(make $(MPFOBJS) F90="$(MPF90)" F90COPT="$(F90COPT)")
@ -35,5 +44,3 @@ veryclean: clean
clean:
/bin/rm -f $(OBJS) $(LOCAL_MODS)
veryclean: clean

@ -53,13 +53,8 @@
!*****************************************************************************
Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_serial_mod
use psb_descriptor_type
Use psb_prec_type
use psb_tools_mod
use psb_const_mod
use psb_error_mod
use psb_penv_mod
use psb_base_mod
use psb_prec_type
Implicit None
! .. Array Arguments ..

@ -40,14 +40,8 @@ subroutine psb_dbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
! where K is a a basic preconditioner stored in prec
!
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_psblas_mod
use psb_const_mod
use psb_error_mod
use psb_penv_mod
use psb_prec_mod, only : psb_bjac_aply
implicit none
type(psb_desc_type),intent(in) :: desc_data
@ -68,6 +62,19 @@ subroutine psb_dbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
external mpi_wtime
character(len=20) :: name, ch_err
interface psb_bjac_aply
subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type), intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: prec
real(kind(0.d0)),intent(inout) :: x(:), y(:)
real(kind(0.d0)),intent(in) :: alpha,beta
character(len=1) :: trans
real(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_dbjac_aply
end interface
name='psb_dbaseprc_aply'
info = 0

@ -36,19 +36,8 @@
!!$
subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
use psb_serial_mod
Use psb_spmat_type
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_tools_mod
use psb_comm_mod
use psb_const_mod
use psb_psblas_mod
use psb_error_mod
use psb_penv_mod
use psb_prec_mod, only: psb_diagsc_bld, psb_ilu_bld, &
& psb_slu_bld,psb_umf_bld
Implicit None
type(psb_dspmat_type), target :: a
@ -57,6 +46,55 @@ subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
integer, intent(out) :: info
character, intent(in), optional :: upd
interface psb_diagsc_bld
subroutine psb_ddiagsc_bld(a,desc_data,p,upd,info)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_dspmat_type), intent(in), target :: a
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
end subroutine psb_ddiagsc_bld
end interface
interface psb_ilu_bld
subroutine psb_dilu_bld(a,desc_data,p,upd,info)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_dspmat_type), intent(in), target :: a
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
end subroutine psb_dilu_bld
end interface
interface psb_slu_bld
subroutine psb_dslu_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
implicit none
type(psb_dspmat_type), intent(inout) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_dslu_bld
end interface
interface psb_umf_bld
subroutine psb_dumf_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
implicit none
type(psb_dspmat_type), intent(inout) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_dumf_bld
end interface
! Local scalars
Integer :: err, nnzero, n_row, n_col,I,j,k,ictxt,&

@ -42,13 +42,8 @@ subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
! but since both are INTENT(IN) this should be legal.
!
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_psblas_mod
use psb_const_mod
use psb_error_mod
use psb_penv_mod
implicit none
type(psb_desc_type), intent(in) :: desc_data
@ -190,12 +185,6 @@ subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
tx = dzero
ty = dzero
!!$ open(50+me)
!!$ call psb_csprt(50+me,prec%av(ap_nd_))
!!$ call flush(50+me)
!!$ close(50+me)
!!$ call psb_barrier(ictxt)
select case(prec%iprcparm(f_type_))
case(f_ilu_n_,f_ilu_e_)
do i=1, prec%iprcparm(jac_sweeps_)

@ -35,18 +35,12 @@
!!$
!!$
subroutine psb_dbldaggrmat(a,desc_a,ac,desc_ac,p,info)
use psb_serial_mod
use psb_penv_mod
use psb_base_mod
use psb_prec_type
use psb_descriptor_type
use psb_spmat_type
use psb_tools_mod
use psb_psblas_mod
use psb_error_mod
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_dspmat_type), intent(inout), target :: ac
type(psb_dspmat_type), intent(out), target :: ac
type(psb_desc_type), intent(in) :: desc_a
type(psb_desc_type), intent(inout) :: desc_ac
type(psb_dbaseprc_type), intent(inout), target :: p
@ -104,11 +98,8 @@ subroutine psb_dbldaggrmat(a,desc_a,ac,desc_ac,p,info)
contains
subroutine raw_aggregate(info)
use psb_base_mod
use psb_prec_type
use psb_const_mod
use psb_psblas_mod
use psb_error_mod
use psb_penv_mod
implicit none
include 'mpif.h'
@ -147,8 +138,8 @@ contains
do i=1, nrow
p%mlia(i) = p%mlia(i) + naggrm1
end do
call psb_halo(p%mlia,desc_a,info)
end if
call psb_halo(p%mlia,desc_a,info)
if(info /= 0) then
call psb_errpush(4010,name,a_err='psb_halo')
@ -178,16 +169,10 @@ contains
nzt = psb_sp_get_nnzeros(b)
j = 0
do i=1, nzt
if ((1<=b%ia2(i)).and.(b%ia2(i)<=nrow)) then
j = j + 1
b%aspk(j) = b%aspk(i)
b%ia1(j) = p%mlia(b%ia1(i))
b%ia2(j) = p%mlia(b%ia2(i))
end if
b%ia1(i) = p%mlia(b%ia1(i))
b%ia2(i) = p%mlia(b%ia2(i))
enddo
b%infoa(psb_nnz_)=j
call psb_fixcoo(b,info)
nzt = psb_sp_get_nnzeros(b)
@ -339,14 +324,10 @@ contains
subroutine smooth_aggregate(info)
use psb_serial_mod
use psb_const_mod
use psb_comm_mod
use psb_tools_mod
use psb_error_mod
use psb_penv_mod
use psb_base_mod
use psb_prec_type
use mpi
implicit none
include 'mpif.h'
integer, intent(out) :: info

@ -35,17 +35,8 @@
!!$
!!$
subroutine psb_ddiagsc_bld(a,desc_a,p,upd,info)
use psb_serial_mod
Use psb_spmat_type
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_tools_mod
use psb_comm_mod
use psb_const_mod
use psb_psblas_mod
use psb_error_mod
use psb_penv_mod
Implicit None
type(psb_dspmat_type), target :: a

@ -35,11 +35,8 @@
!!$
!!$
subroutine psb_dgenaggrmap(aggr_type,a,desc_a,nlaggr,ilaggr,info)
use psb_spmat_type
use psb_serial_mod
use psb_descriptor_type
use psb_error_mod
use psb_penv_mod
use psb_base_mod
use psb_prec_type
implicit none
integer, intent(in) :: aggr_type
type(psb_dspmat_type), intent(in) :: a

@ -50,17 +50,8 @@
!* *
!*****************************************************************************
subroutine psb_dilu_bld(a,desc_a,p,upd,info)
use psb_serial_mod
use psb_const_mod
use psb_base_mod
use psb_prec_type
use psb_descriptor_type
use psb_spmat_type
use psb_tools_mod
use psb_psblas_mod
use psb_error_mod
use psb_realloc_mod
use psb_penv_mod
use psb_prec_mod, only : psb_as_matbld, psb_ilu_fct, psb_sp_renum
implicit none
!
! .. Scalar Arguments ..
@ -84,6 +75,45 @@ subroutine psb_dilu_bld(a,desc_a,p,upd,info)
integer :: ictxt,np,me
character(len=20) :: name, ch_err
interface psb_ilu_fct
subroutine psb_dilu_fct(a,l,u,d,info,blck)
use psb_base_mod
integer, intent(out) :: info
type(psb_dspmat_type),intent(in) :: a
type(psb_dspmat_type),intent(inout) :: l,u
type(psb_dspmat_type),intent(in), optional, target :: blck
real(kind(1.d0)), intent(inout) :: d(:)
end subroutine psb_dilu_fct
end interface
interface psb_asmatbld
Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_dspmat_type), Intent(in) :: a
Type(psb_dspmat_type), Intent(inout) :: blk
Type(psb_desc_type), Intent(inout) :: desc_p
Type(psb_desc_type), Intent(in) :: desc_data
Character, Intent(in) :: upd
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_dasmatbld
end interface
interface psb_sp_renum
subroutine psb_dsp_renum(a,desc_a,blck,p,atmp,info)
use psb_base_mod
use psb_prec_type
implicit none
type(psb_dspmat_type), intent(in) :: a,blck
type(psb_dspmat_type), intent(inout) :: atmp
type(psb_dbaseprc_type), intent(inout) :: p
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
end subroutine psb_dsp_renum
end interface
if(psb_get_errstatus().ne.0) return
info=0
name='psb_ilu_bld'
@ -115,18 +145,18 @@ subroutine psb_dilu_bld(a,desc_a,p,upd,info)
t1= mpi_wtime()
if(debug) write(0,*)me,': calling psb_as_matbld',p%iprcparm(p_type_),p%iprcparm(n_ovr_)
if(debug) write(0,*)me,': calling psb_asmatbld',p%iprcparm(p_type_),p%iprcparm(n_ovr_)
if (debug) call psb_barrier(ictxt)
call psb_as_matbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info)
if(info/=0) then
info=4010
ch_err='psb_as_matbld'
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
t2= mpi_wtime()
if (debug) write(0,*)me,': out of psb_as_matbld'
if (debug) write(0,*)me,': out of psb_asmatbld'
if (debug) call psb_barrier(ictxt)
if (allocated(p%av)) then
@ -178,7 +208,7 @@ subroutine psb_dilu_bld(a,desc_a,p,upd,info)
if (debug) then
write(0,*) me,'Done psb_as_matbld'
write(0,*) me,'Done psb_asmatbld'
call psb_barrier(ictxt)
endif

@ -41,11 +41,7 @@ subroutine psb_dilu_fct(a,l,u,d,info,blck)
! into L/D/U.
!
!
use psb_spmat_type
use psb_serial_mod
use psb_tools_mod
use psb_error_mod
use psb_const_mod
use psb_base_mod
implicit none
! .. Scalar Arguments ..
integer, intent(out) :: info

@ -83,15 +83,8 @@ subroutine psb_dmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
! 6. baseprecv(ilev)%nlaggr Number of aggregates on the various procs.
!
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_psblas_mod
use psb_penv_mod
use psb_const_mod
use psb_error_mod
use psb_penv_mod
use psb_prec_mod, only : psb_baseprc_aply
implicit none
type(psb_desc_type),intent(in) :: desc_data
@ -120,6 +113,20 @@ subroutine psb_dmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
end type psb_mlprec_wrk_type
type(psb_mlprec_wrk_type), allocatable :: mlprec_wrk(:)
interface psb_baseprc_aply
subroutine psb_dbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: prec
real(kind(0.d0)),intent(inout) :: x(:), y(:)
real(kind(0.d0)),intent(in) :: alpha,beta
character(len=1) :: trans
real(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_dbaseprc_aply
end interface
name='psb_mlprc_aply'
info = 0
call psb_erractionsave(err_act)

@ -36,14 +36,8 @@
!!$
subroutine psb_dmlprc_bld(a,desc_a,p,info)
use psb_serial_mod
use psb_tools_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_const_mod
use psb_error_mod
use psb_penv_mod
use psb_prec_mod, only : psb_genaggrmap, psb_bldaggrmat, psb_baseprc_bld
implicit none
type(psb_dspmat_type), intent(in), target :: a
@ -58,6 +52,43 @@ subroutine psb_dmlprc_bld(a,desc_a,p,info)
logical, parameter :: debug=.false.
type(psb_dspmat_type) :: ac
interface psb_baseprc_bld
subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
use psb_base_mod
use psb_prec_type
type(psb_dspmat_type), target :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dbaseprc_type),intent(inout) :: p
integer, intent(out) :: info
character, intent(in), optional :: upd
end subroutine psb_dbaseprc_bld
end interface
interface psb_genaggrmap
subroutine psb_dgenaggrmap(aggr_type,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use psb_prec_type
implicit none
integer, intent(in) :: aggr_type
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer, allocatable :: ilaggr(:),nlaggr(:)
integer, intent(out) :: info
end subroutine psb_dgenaggrmap
end interface
interface psb_bldaggrmat
subroutine psb_dbldaggrmat(a,desc_a,ac,desc_ac,p,info)
use psb_base_mod
use psb_prec_type
type(psb_dspmat_type), intent(in), target :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dspmat_type), intent(out),target :: ac
type(psb_desc_type), intent(inout) :: desc_ac
type(psb_dbaseprc_type), intent(inout), target :: p
integer, intent(out) :: info
end subroutine psb_dbldaggrmat
end interface
integer :: ictxt, np, me

@ -36,14 +36,8 @@
!!$
subroutine psb_dprc_aply(prec,x,y,desc_data,info,trans, work)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_psblas_mod
use psb_const_mod
use psb_error_mod
use psb_penv_mod
use psb_prec_mod, only: psb_mlprc_aply, psb_baseprc_aply
implicit none
type(psb_desc_type),intent(in) :: desc_data
@ -61,6 +55,34 @@ subroutine psb_dprc_aply(prec,x,y,desc_data,info,trans, work)
external mpi_wtime
character(len=20) :: name
interface psb_baseprc_aply
subroutine psb_dbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: prec
real(kind(0.d0)),intent(inout) :: x(:), y(:)
real(kind(0.d0)),intent(in) :: alpha,beta
character(len=1) :: trans
real(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_dbaseprc_aply
end interface
interface psb_mlprc_aply
subroutine psb_dmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: baseprecv(:)
real(kind(0.d0)),intent(in) :: alpha,beta
real(kind(0.d0)),intent(inout) :: x(:), y(:)
character :: trans
real(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_dmlprc_aply
end interface
name='psb_dprc_aply'
info = 0
call psb_erractionsave(err_act)
@ -159,18 +181,20 @@ end subroutine psb_dprc_aply
!!$
subroutine psb_dprc_aply1(prec,x,desc_data,info,trans)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_psblas_mod
use psb_const_mod
use psb_error_mod
use psb_penv_mod
implicit none
interface psb_prc_aply
type(psb_desc_type),intent(in) :: desc_data
type(psb_dprec_type), intent(in) :: prec
real(kind(0.d0)),intent(inout) :: x(:)
integer, intent(out) :: info
character(len=1), optional :: trans
logical,parameter :: debug=.false., debugprt=.false.
interface
subroutine psb_dprc_aply(prec,x,y,desc_data,info,trans, work)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
implicit none
@ -182,13 +206,6 @@ subroutine psb_dprc_aply1(prec,x,desc_data,info,trans)
real(kind(0.d0)), optional, target :: work(:)
end subroutine psb_dprc_aply
end interface
type(psb_desc_type),intent(in) :: desc_data
type(psb_dprec_type), intent(in) :: prec
real(kind(0.d0)),intent(inout) :: x(:)
integer, intent(out) :: info
character(len=1), optional :: trans
logical,parameter :: debug=.false., debugprt=.false.
! Local variables
character :: trans_
@ -214,7 +231,7 @@ subroutine psb_dprc_aply1(prec,x,desc_data,info,trans)
goto 9999
end if
if (debug) write(0,*) 'Prc_aply1 Size(x) ',size(x), size(ww),size(w1)
call psb_prc_aply(prec,x,ww,desc_data,info,trans_,work=w1)
call psb_dprc_aply(prec,x,ww,desc_data,info,trans_,work=w1)
if(info /=0) goto 9999
x(:) = ww(:)
deallocate(ww,W1)

@ -36,17 +36,9 @@
!!$
subroutine psb_dprecbld(a,desc_a,p,info,upd)
use psb_serial_mod
Use psb_spmat_type
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_tools_mod
use psb_comm_mod
use psb_const_mod
use psb_psblas_mod
use psb_error_mod
use psb_penv_mod
use psb_prec_mod, only: psb_mlprc_bld, psb_baseprc_bld
use psb_prec_mod
Implicit None
type(psb_dspmat_type), target :: a
@ -60,9 +52,9 @@ subroutine psb_dprecbld(a,desc_a,p,info,upd)
integer :: int_err(5)
character :: iupd
logical, parameter :: debug=.false., filedump=.false.
logical, parameter :: debug=.false.
integer,parameter :: iroot=0,iout=60,ilout=40
character(len=20) :: name, ch_err,dumpname
character(len=20) :: name, ch_err
if(psb_get_errstatus().ne.0) return
info=0
@ -142,18 +134,7 @@ subroutine psb_dprecbld(a,desc_a,p,info,upd)
if (debug) then
write(0,*) 'Return from ',i-1,' call to mlprcbld ',info
endif
if (filedump) then
write(dumpname,'(a,i2.2,a,i2.2,a)'),'ac_lev_',i,'.',me,'.out'
open(20,file=dumpname)
call psb_csprt(20,p%baseprecv(i)%av(ac_))
call flush(20)
close(20)
write(dumpname,'(a,i2.2,a,i2.2,a)'),'nd_lev_',i,'.',me,'.out'
open(20,file=dumpname)
call psb_csprt(20,p%baseprecv(i)%av(ap_nd_))
call flush(20)
close(20)
end if
end do
endif

@ -35,16 +35,9 @@
!!$
!!$
subroutine psb_dprecfree(p,info)
!...free sparse matrix structure...
use psb_descriptor_type
use psb_serial_mod
use psb_const_mod
use psb_base_mod
use psb_prec_type
use psb_tools_mod
use psb_error_mod
implicit none
!....parameters...
type(psb_dprec_type), intent(inout) :: p
integer, intent(out) :: info

@ -36,10 +36,8 @@
!!$
subroutine psb_dprecset(p,ptype,info,iv,rs,rv,ilev,nlev)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_string_mod
implicit none
type(psb_dprec_type), intent(inout) :: p
character(len=*), intent(in) :: ptype

@ -35,13 +35,9 @@
!!$
!!$
subroutine psb_dslu_bld(a,desc_a,p,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_tools_mod
use psb_const_mod
use psb_penv_mod
use psb_prec_mod, only: psb_as_matbld
implicit none
type(psb_dspmat_type), intent(inout) :: a
@ -57,12 +53,27 @@ subroutine psb_dslu_bld(a,desc_a,p,info)
logical, parameter :: debug=.false.
character(len=20) :: name, ch_err
interface psb_asmatbld
Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_dspmat_type), Intent(in) :: a
Type(psb_dspmat_type), Intent(inout) :: blk
Type(psb_desc_type), Intent(inout) :: desc_p
Type(psb_desc_type), Intent(in) :: desc_data
Character, Intent(in) :: upd
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_dasmatbld
end interface
if(psb_get_errstatus().ne.0) return
info=0
name='psb_slu_bld'
call psb_erractionsave(err_act)
ictxt = desc_A%matrix_data(psb_ctxt_)
ictxt = desc_a%matrix_data(psb_ctxt_)
call psb_info(ictxt, me, np)
@ -88,18 +99,18 @@ subroutine psb_dslu_bld(a,desc_a,p,info)
write(0,*) me, 'SPLUBLD: Done csdp',info,nza,atmp%m,atmp%k
call psb_barrier(ictxt)
endif
call psb_as_matbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info,outfmt=fmt)
if(info /= 0) then
info=4010
ch_err='psb_as_matbld'
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
nzb = blck%infoa(psb_nnz_)
if (Debug) then
write(0,*) me, 'SPLUBLD: Done as_matbld',info,nzb,blck%fida
write(0,*) me, 'SPLUBLD: Done asmatbld',info,nzb,blck%fida
call psb_barrier(ictxt)
endif
if (nzb > 0 ) then

@ -35,15 +35,8 @@
!!$
!!$
subroutine psb_dsp_renum(a,desc_a,blck,p,atmp,info)
use psb_serial_mod
use psb_const_mod
use psb_base_mod
use psb_prec_type
use psb_descriptor_type
use psb_spmat_type
use psb_tools_mod
use psb_psblas_mod
use psb_error_mod
use psb_penv_mod
implicit none
! .. array Arguments ..
@ -377,7 +370,6 @@ contains
integer,dimension(:,:),allocatable::NDstk
integer,dimension(:),allocatable::iOld,renum,ndeg,lvl,lvls1,lvls2,ccstor
!--- Per la common area.
character(len=20) :: name, ch_err

@ -35,13 +35,8 @@
!!$
!!$
subroutine psb_dumf_bld(a,desc_a,p,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_tools_mod
use psb_const_mod
use psb_penv_mod
use psb_prec_mod, only: psb_as_matbld
implicit none
type(psb_dspmat_type), intent(inout) :: a
@ -58,6 +53,21 @@ subroutine psb_dumf_bld(a,desc_a,p,info)
logical, parameter :: debug=.false.
character(len=20) :: name, ch_err
interface psb_asmatbld
Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_dspmat_type), Intent(in) :: a
Type(psb_dspmat_type), Intent(inout) :: blk
Type(psb_desc_type), Intent(inout) :: desc_p
Type(psb_desc_type), Intent(in) :: desc_data
Character, Intent(in) :: upd
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_dasmatbld
end interface
info=0
name='psb_umf_bld'
call psb_erractionsave(err_act)
@ -89,18 +99,18 @@ subroutine psb_dumf_bld(a,desc_a,p,info)
write(0,*) me, 'UMFBLD: Done csdp',info,nza,atmp%m,atmp%k,nzb
call psb_barrier(ictxt)
endif
call psb_as_matbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info,outfmt=fmt)
if(info /= 0) then
info=4010
ch_err='psb_as_matbld'
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
nzb = psb_sp_get_nnzeros(blck)
if (Debug) then
write(0,*) me, 'UMFBLD: Done as_matbld',info,nzb,blck%fida
write(0,*) me, 'UMFBLD: Done asmatbld',info,nzb,blck%fida
call psb_barrier(ictxt)
endif
if (nzb > 0 ) then

@ -1,7 +1,13 @@
!!$
!!$
!!$ MD2P4
!!$ Multilevel Domain Decomposition Parallel Preconditioner Package for PSBLAS
!!$ for
!!$ Parallel Sparse BLAS v2.0
!!$ (C) Copyright 2006 Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$ Daniela di Serafino Second University of Naples
!!$ Pasqua D'Ambra ICAR-CNR
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
@ -11,14 +17,14 @@
!!$ 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
!!$ 3. The name of the MD2P4 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
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE MD2P4 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
@ -28,14 +34,13 @@
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
module psb_prec_mod
use psb_prec_type
interface psb_precbld
subroutine psb_dprecbld(a,desc_a,prec,info,upd)
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
implicit none
type(psb_dspmat_type), intent(in), target :: a
@ -45,7 +50,7 @@ module psb_prec_mod
character, intent(in),optional :: upd
end subroutine psb_dprecbld
subroutine psb_zprecbld(a,desc_a,prec,info,upd)
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
implicit none
type(psb_zspmat_type), intent(in), target :: a
@ -58,8 +63,7 @@ module psb_prec_mod
interface psb_precset
subroutine psb_dprecset(prec,ptype,info,iv,rs,rv,ilev,nlev)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
implicit none
type(psb_dprec_type), intent(inout) :: prec
@ -71,8 +75,7 @@ module psb_prec_mod
real(kind(1.d0)), optional, intent(in) :: rv(:)
end subroutine psb_dprecset
subroutine psb_zprecset(prec,ptype,info,iv,rs,rv,ilev,nlev)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
implicit none
type(psb_zprec_type), intent(inout) :: prec
@ -88,17 +91,13 @@ module psb_prec_mod
interface psb_precfree
subroutine psb_dprecfree(p,info)
use psb_descriptor_type
use psb_serial_mod
use psb_const_mod
use psb_base_mod
use psb_prec_type
type(psb_dprec_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_dprecfree
subroutine psb_zprecfree(p,info)
use psb_descriptor_type
use psb_serial_mod
use psb_const_mod
use psb_base_mod
use psb_prec_type
type(psb_zprec_type), intent(inout) :: p
integer, intent(out) :: info
@ -107,8 +106,7 @@ module psb_prec_mod
interface psb_prc_aply
subroutine psb_dprc_aply(prec,x,y,desc_data,info,trans,work)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_dprec_type), intent(in) :: prec
@ -118,8 +116,7 @@ module psb_prec_mod
real(kind(0.d0)),intent(inout), optional, target :: work(:)
end subroutine psb_dprc_aply
subroutine psb_dprc_aply1(prec,x,desc_data,info,trans)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_dprec_type), intent(in) :: prec
@ -128,8 +125,7 @@ module psb_prec_mod
character(len=1), optional :: trans
end subroutine psb_dprc_aply1
subroutine psb_zprc_aply(prec,x,y,desc_data,info,trans,work)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_zprec_type), intent(in) :: prec
@ -139,8 +135,7 @@ module psb_prec_mod
complex(kind(0.d0)),intent(inout), optional, target :: work(:)
end subroutine psb_zprc_aply
subroutine psb_zprc_aply1(prec,x,desc_data,info,trans)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_zprec_type), intent(in) :: prec
@ -152,8 +147,7 @@ module psb_prec_mod
interface psb_baseprc_bld
subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
Use psb_spmat_type
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_dspmat_type), target :: a
type(psb_desc_type), intent(in), target :: desc_a
@ -162,8 +156,7 @@ module psb_prec_mod
character, intent(in), optional :: upd
end subroutine psb_dbaseprc_bld
subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
Use psb_spmat_type
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_zspmat_type), target :: a
type(psb_desc_type), intent(in), target :: desc_a
@ -175,24 +168,16 @@ module psb_prec_mod
interface psb_mlprc_bld
subroutine psb_dmlprc_bld(a,desc_a,p,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_const_mod
implicit none
type(psb_dspmat_type), intent(in), target :: a
type(psb_desc_type), intent(in), target :: desc_a
type(psb_dbaseprc_type), intent(inout), target :: p
integer, intent(out) :: info
end subroutine psb_dmlprc_bld
subroutine psb_zmlprc_bld(a,desc_a,p,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_const_mod
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_desc_type), intent(in), target :: desc_a
type(psb_zbaseprc_type), intent(inout),target :: p
@ -203,7 +188,7 @@ module psb_prec_mod
interface psb_baseprc_aply
subroutine psb_dbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: prec
@ -215,7 +200,7 @@ module psb_prec_mod
end subroutine psb_dbaseprc_aply
subroutine psb_zbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(in) :: prec
@ -229,7 +214,7 @@ module psb_prec_mod
interface psb_mlprc_aply
subroutine psb_dmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: baseprecv(:)
@ -240,7 +225,7 @@ module psb_prec_mod
integer, intent(out) :: info
end subroutine psb_dmlprc_aply
subroutine psb_zmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(in) :: baseprecv(:)
@ -254,7 +239,7 @@ module psb_prec_mod
interface psb_bjac_aply
subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_desc_type), intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: prec
@ -266,7 +251,7 @@ module psb_prec_mod
end subroutine psb_dbjac_aply
subroutine psb_zbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
type(psb_desc_type), intent(in) :: desc_data
type(psb_zbaseprc_type), intent(in) :: prec
@ -281,8 +266,7 @@ module psb_prec_mod
interface psb_diagsc_bld
subroutine psb_ddiagsc_bld(a,desc_data,p,upd,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_dspmat_type), intent(in), target :: a
@ -291,8 +275,7 @@ module psb_prec_mod
character, intent(in) :: upd
end subroutine psb_ddiagsc_bld
subroutine psb_zdiagsc_bld(a,desc_data,p,upd,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_zspmat_type), intent(in), target :: a
@ -304,8 +287,7 @@ module psb_prec_mod
interface psb_ilu_bld
subroutine psb_dilu_bld(a,desc_data,p,upd,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_dspmat_type), intent(in), target :: a
@ -314,8 +296,7 @@ module psb_prec_mod
character, intent(in) :: upd
end subroutine psb_dilu_bld
subroutine psb_zilu_bld(a,desc_data,p,upd,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_zspmat_type), intent(in), target :: a
@ -327,24 +308,16 @@ module psb_prec_mod
interface psb_slu_bld
subroutine psb_dslu_bld(a,desc_a,p,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_const_mod
implicit none
type(psb_dspmat_type), intent(inout) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_dslu_bld
subroutine psb_zslu_bld(a,desc_a,p,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_const_mod
implicit none
type(psb_zspmat_type), intent(inout) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_zbaseprc_type), intent(inout) :: p
@ -354,24 +327,16 @@ module psb_prec_mod
interface psb_umf_bld
subroutine psb_dumf_bld(a,desc_a,p,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_const_mod
implicit none
type(psb_dspmat_type), intent(inout) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_dumf_bld
subroutine psb_zumf_bld(a,desc_a,p,info)
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_const_mod
implicit none
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_zbaseprc_type), intent(inout) :: p
@ -382,7 +347,7 @@ module psb_prec_mod
interface psb_ilu_fct
subroutine psb_dilu_fct(a,l,u,d,info,blck)
use psb_spmat_type
use psb_base_mod
integer, intent(out) :: info
type(psb_dspmat_type),intent(in) :: a
type(psb_dspmat_type),intent(inout) :: l,u
@ -390,7 +355,7 @@ module psb_prec_mod
real(kind(1.d0)), intent(inout) :: d(:)
end subroutine psb_dilu_fct
subroutine psb_zilu_fct(a,l,u,d,info,blck)
use psb_spmat_type
use psb_base_mod
integer, intent(out) :: info
type(psb_zspmat_type),intent(in) :: a
type(psb_zspmat_type),intent(inout) :: l,u
@ -401,9 +366,8 @@ module psb_prec_mod
interface psb_as_matbld
Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_serial_mod
Use psb_descriptor_type
Use psb_prec_type
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_dspmat_type), Intent(in) :: a
Type(psb_dspmat_type), Intent(inout) :: blk
@ -414,9 +378,8 @@ module psb_prec_mod
character(len=5), optional :: outfmt
end Subroutine psb_dasmatbld
Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_serial_mod
Use psb_descriptor_type
Use psb_prec_type
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_zspmat_type), Intent(in) :: a
Type(psb_zspmat_type), Intent(inout) :: blk
@ -430,12 +393,8 @@ module psb_prec_mod
interface psb_sp_renum
subroutine psb_dsp_renum(a,desc_a,blck,p,atmp,info)
use psb_base_mod
use psb_prec_type
use psb_descriptor_type
use psb_spmat_type
implicit none
! .. array Arguments ..
type(psb_dspmat_type), intent(in) :: a,blck
type(psb_dspmat_type), intent(inout) :: atmp
type(psb_dbaseprc_type), intent(inout) :: p
@ -443,12 +402,8 @@ module psb_prec_mod
integer, intent(out) :: info
end subroutine psb_dsp_renum
subroutine psb_zsp_renum(a,desc_a,blck,p,atmp,info)
use psb_base_mod
use psb_prec_type
use psb_descriptor_type
use psb_spmat_type
implicit none
! .. array Arguments ..
type(psb_zspmat_type), intent(in) :: a,blck
type(psb_zspmat_type), intent(inout) :: atmp
type(psb_zbaseprc_type), intent(inout) :: p
@ -460,9 +415,8 @@ module psb_prec_mod
interface psb_genaggrmap
subroutine psb_dgenaggrmap(aggr_type,a,desc_a,nlaggr,ilaggr,info)
use psb_spmat_type
use psb_descriptor_type
implicit none
use psb_base_mod
use psb_prec_type
integer, intent(in) :: aggr_type
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
@ -470,9 +424,8 @@ module psb_prec_mod
integer, intent(out) :: info
end subroutine psb_dgenaggrmap
subroutine psb_zgenaggrmap(aggr_type,a,desc_a,nlaggr,ilaggr,info)
use psb_spmat_type
use psb_descriptor_type
implicit none
use psb_base_mod
use psb_prec_type
integer, intent(in) :: aggr_type
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
@ -483,9 +436,8 @@ module psb_prec_mod
interface psb_bldaggrmat
subroutine psb_dbldaggrmat(a,desc_a,ac,desc_ac,p,info)
use psb_base_mod
use psb_prec_type
use psb_descriptor_type
use psb_spmat_type
type(psb_dspmat_type), intent(in), target :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dspmat_type), intent(inout),target :: ac
@ -494,9 +446,8 @@ module psb_prec_mod
integer, intent(out) :: info
end subroutine psb_dbldaggrmat
subroutine psb_zbldaggrmat(a,desc_a,ac,desc_ac,p,info)
use psb_base_mod
use psb_prec_type
use psb_descriptor_type
use psb_spmat_type
type(psb_zspmat_type), intent(in), target :: a
type(psb_zbaseprc_type), intent(inout),target :: p
type(psb_zspmat_type), intent(inout),target :: ac

@ -1,7 +1,13 @@
!!$
!!$
!!$ MD2P4
!!$ Multilevel Domain Decomposition Parallel Preconditioner Package for PSBLAS
!!$ for
!!$ Parallel Sparse BLAS v2.0
!!$ (C) Copyright 2006 Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$ Daniela di Serafino Second University of Naples
!!$ Pasqua D'Ambra ICAR-CNR
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
@ -11,14 +17,14 @@
!!$ 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
!!$ 3. The name of the MD2P4 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
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE MD2P4 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
@ -27,17 +33,13 @@
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
module psb_prec_type
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!! Module to define PREC_DATA, !!
!! structure for preconditioning. !!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
module psb_prec_type
use psb_const_mod
use psb_spmat_type
use psb_descriptor_type
use psb_base_mod
integer, parameter :: min_prec_=0, noprec_=0, diagsc_=1, bja_=2,&
& asm_=3, ras_=5, ash_=4, rash_=6, ras2lv_=7, ras2lvm_=8,&
@ -653,9 +655,8 @@ contains
end subroutine psb_dcheck_def
subroutine psb_dbase_precfree(p,info)
use psb_serial_mod
use psb_descriptor_type
use psb_tools_mod
use psb_base_mod
type(psb_dbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
integer :: i
@ -728,7 +729,8 @@ contains
end subroutine psb_dbase_precfree
subroutine psb_nullify_dbaseprec(p)
use psb_descriptor_type
use psb_base_mod
type(psb_dbaseprc_type), intent(inout) :: p
nullify(p%base_a)
@ -739,9 +741,7 @@ contains
end subroutine psb_nullify_dbaseprec
subroutine psb_zbase_precfree(p,info)
use psb_serial_mod
use psb_descriptor_type
use psb_tools_mod
use psb_base_mod
type(psb_zbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
integer :: i
@ -809,7 +809,8 @@ contains
end subroutine psb_zbase_precfree
subroutine psb_nullify_zbaseprec(p)
use psb_descriptor_type
use psb_base_mod
type(psb_zbaseprc_type), intent(inout) :: p

@ -35,7 +35,8 @@
*
*/
/* This file is an interface to the UMFPACK routines for sparse
factorization.
factorization. It was obtained by adapting umfpack_di_demo
under the original copyright terms reproduced below.
PSBLAS v 2.0 */

@ -53,13 +53,8 @@
!*****************************************************************************
Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_serial_mod
use psb_descriptor_type
Use psb_prec_type
use psb_tools_mod
use psb_const_mod
use psb_error_mod
use psb_penv_mod
use psb_base_mod
use psb_prec_type
Implicit None
! .. Array Arguments ..

@ -39,15 +39,8 @@ subroutine psb_zbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
! Compute Y <- beta*Y + alpha*K^-1 X
! where K is a a basic preconditioner stored in prec
!
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_psblas_mod
use psb_const_mod
use psb_error_mod
use psb_penv_mod
use psb_prec_mod, only : psb_bjac_aply
implicit none
type(psb_desc_type),intent(in) :: desc_data
@ -68,6 +61,20 @@ subroutine psb_zbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
external mpi_wtime
character(len=20) :: name, ch_err
interface psb_bjac_aply
subroutine psb_zbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type), intent(in) :: desc_data
type(psb_zbaseprc_type), intent(in) :: prec
complex(kind(0.d0)),intent(inout) :: x(:), y(:)
complex(kind(0.d0)),intent(in) :: alpha,beta
character(len=1) :: trans
complex(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_zbjac_aply
end interface
name='psb_zbaseprc_aply'
info = 0
call psb_erractionsave(err_act)

@ -36,18 +36,8 @@
!!$
subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
use psb_serial_mod
Use psb_spmat_type
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_tools_mod
use psb_comm_mod
use psb_const_mod
use psb_psblas_mod
use psb_error_mod
use psb_penv_mod
use psb_prec_mod, only: psb_diagsc_bld, psb_ilu_bld, &
& psb_slu_bld,psb_umf_bld
Implicit None
type(psb_zspmat_type), target :: a
@ -56,6 +46,51 @@ subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
integer, intent(out) :: info
character, intent(in), optional :: upd
interface psb_diagsc_bld
subroutine psb_zdiagsc_bld(a,desc_data,p,upd,info)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_zspmat_type), intent(in), target :: a
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
end subroutine psb_zdiagsc_bld
end interface
interface psb_ilu_bld
subroutine psb_zilu_bld(a,desc_data,p,upd,info)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_zspmat_type), intent(in), target :: a
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
end subroutine psb_zilu_bld
end interface
interface psb_slu_bld
subroutine psb_zslu_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
type(psb_zspmat_type), intent(inout) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_zbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_zslu_bld
end interface
interface psb_umf_bld
subroutine psb_zumf_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_zbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_zumf_bld
end interface
! Local scalars
Integer :: err, nnzero, n_row, n_col,I,j,k,ictxt,&

@ -42,13 +42,8 @@ subroutine psb_zbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
! but since both are INTENT(IN) this should be legal.
!
use psb_serial_mod
use psb_descriptor_type
use psb_base_mod
use psb_prec_type
use psb_psblas_mod
use psb_const_mod
use psb_error_mod
use psb_penv_mod
implicit none
type(psb_desc_type), intent(in) :: desc_data

@ -35,19 +35,13 @@
!!$
!!$
subroutine psb_zbldaggrmat(a,desc_a,ac,desc_ac,p,info)
use psb_serial_mod
use psb_penv_mod
use psb_base_mod
use psb_prec_type
use psb_descriptor_type
use psb_spmat_type
use psb_tools_mod
use psb_psblas_mod
use psb_error_mod
implicit none
type(psb_zspmat_type), intent(in), target :: a
type(psb_zbaseprc_type), intent(inout),target :: p
type(psb_zspmat_type), intent(inout), target :: ac
type(psb_zspmat_type), intent(out), target :: ac
type(psb_desc_type), intent(in) :: desc_a
type(psb_desc_type), intent(inout) :: desc_ac
integer, intent(out) :: info
@ -103,16 +97,13 @@ subroutine psb_zbldaggrmat(a,desc_a,ac,desc_ac,p,info)
contains
subroutine raw_aggregate(info)
use psb_base_mod
use psb_prec_type
use psb_const_mod
use psb_psblas_mod
use psb_error_mod
use psb_penv_mod
implicit none
include 'mpif.h'
integer, intent(out) :: info
type(psb_zspmat_type) :: b
type(psb_zspmat_type) :: b, tmp
integer, pointer :: nzbr(:), idisp(:)
integer :: ictxt, nrow, nglob, ncol, ntaggr, nzac, ip, ndx,&
& naggr, np, me, nzt,jl,nzl,nlr,&
@ -146,8 +137,8 @@ contains
do i=1, nrow
p%mlia(i) = p%mlia(i) + naggrm1
end do
call psb_halo(p%mlia,desc_a,info)
end if
call psb_halo(p%mlia,desc_a,info)
if(info /= 0) then
call psb_errpush(4010,name,a_err='psb_halo')
@ -177,16 +168,10 @@ contains
nzt = psb_sp_get_nnzeros(b)
j = 0
do i=1, nzt
if ((1<=b%ia2(i)).and.(b%ia2(i)<=nrow)) then
j = j + 1
b%aspk(j) = b%aspk(i)
b%ia1(j) = p%mlia(b%ia1(i))
b%ia2(j) = p%mlia(b%ia2(i))
end if
b%ia1(i) = p%mlia(b%ia1(i))
b%ia2(i) = p%mlia(b%ia2(i))
enddo
b%infoa(psb_nnz_)=j
call psb_fixcoo(b,info)
nzt = psb_sp_get_nnzeros(b)
@ -338,14 +323,10 @@ contains
subroutine smooth_aggregate(info)
use psb_serial_mod
use psb_const_mod
use psb_comm_mod
use psb_tools_mod
use psb_error_mod
use psb_penv_mod
use psb_base_mod
use psb_prec_type
use mpi
implicit none
include 'mpif.h'
integer, intent(out) :: info

Some files were not shown because too many files have changed in this diff Show More

Loading…
Cancel
Save