Merged changes for serial version.

psblas3-type-indexed
Salvatore Filippone 18 years ago
parent 11cf3b128d
commit 719d69246e

@ -2,6 +2,10 @@ Changelog. A lot less detailed than usual, at least for past
history. history.
2007/02/01: Merged serial version: we provide a minimal fake mpi to
allow compiling and running without mpi and blacs. Only
tested with gnu42 so far.
2007/01/23: Defined new field ext_index in desc_type, and 2007/01/23: Defined new field ext_index in desc_type, and
fixed long standing inconsistency in usage of overlap for fixed long standing inconsistency in usage of overlap for
AS preconditioners. Modified halo to accept selector for AS preconditioners. Modified halo to accept selector for

@ -1,6 +1,6 @@
.mod=.mod .mod=.mod
.fh=.fh .fh=.fh
.SUFFIXES: .f90 $(.mod) .F90 .SUFFIXES: .f90 $(.mod) .F90 .F
####################### Section 1 ####################### ####################### Section 1 #######################
@ -78,8 +78,12 @@ $(.mod).o:
$(F90) $(F90COPT) $(INCDIRS) -c $< $(F90) $(F90COPT) $(INCDIRS) -c $<
.f90.o: .f90.o:
$(F90) $(F90COPT) $(INCDIRS) -c $< $(F90) $(F90COPT) $(INCDIRS) -c $<
.F.o:
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<
.F90.o: .F90.o:
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $< $(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<
.F90$(.mod):
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<

@ -0,0 +1,92 @@
.mod=.mod
.fh=.fh
.SUFFIXES: .f90 $(.mod) .F90 .F
####################### Section 1 #######################
# Define your compilers and compiler flags here #
##########################################################
F90=/usr/local/gcc42/bin/gfortran
FC=$(F90)
F77=$(FC)
CC=/usr/local/gcc42/bin/gcc
F90COPT=-O3 -ggdb
FCOPT=-O3 -ggdb
CCOPT=-O3 -ggdb
####################### Section 2 #######################
# Define your linker and linker flags here #
##########################################################
F90LINK=$(FC)
FLINK=$(FC)
MPF90=$(FC)
MPF77=$(FC)
MPCC=$(CC)
####################### Section 3 #######################
# Specify paths to libraries #
##########################################################
BLAS=-lblas-gcc42 -L$(HOME)/LIB
# No BLACS in serialMPI. But we need the fakempi.o
#BLACS=-lmpiblacs-gcc42 -L$(HOME)/LIB
EXTRA_COBJS=fakempi.o
####################### Section 4 #######################
# Other useful tools&defines #
##########################################################
SLUDIR=/usr/local/SuperLU_3.0
SLU=-lslu_lx_gcc42 -L$(SLUDIR)
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)
FDEFINES=-DHAVE_MOVE_ALLOC -DSERIAL_MPI
AR=ar -cur
RANLIB=ranlib
####################### Section 5 #######################
# Do not edit this #
##########################################################
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:
.f.o:
$(FC) $(FCOPT) $(INCDIRS) -c $<
.c.o:
$(CC) $(CCOPT) $(INCDIRS) $(CDEFINES) -c $<
.f$(.mod):
$(F90) $(FCOPT) $(INCDIRS) -c $<
.f90$(.mod):
$(F90) $(F90COPT) $(INCDIRS) -c $<
.f90.o:
$(F90) $(F90COPT) $(INCDIRS) -c $<
.F.o:
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<
.F90.o:
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<
.F90$(.mod):
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<

@ -1,6 +1,6 @@
.mod=.mod .mod=.mod
.fh=.fh .fh=.fh
.SUFFIXES: .f90 $(.mod) .F90 .SUFFIXES: .f90 $(.mod) .F90 .F
####################### Section 1 ####################### ####################### Section 1 #######################
@ -80,8 +80,12 @@ $(.mod).o:
$(F90) $(F90COPT) $(INCDIRS) -c $< $(F90) $(F90COPT) $(INCDIRS) -c $<
.f90.o: .f90.o:
$(F90) $(F90COPT) $(INCDIRS) -c $< $(F90) $(F90COPT) $(INCDIRS) -c $<
.F.o:
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<
.F90.o: .F90.o:
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $< $(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<
.F90$(.mod):
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<

@ -1,13 +1,13 @@
.mod=.mod .mod=.mod
.fh=.fh .fh=.fh
.SUFFIXES: .f90 $(.mod) .F90 .SUFFIXES: .f90 $(.mod) .F90 .F
####################### Section 1 ####################### ####################### Section 1 #######################
# Define your compilers and compiler flags here # # Define your compilers and compiler flags here #
########################################################## ##########################################################
F90=xlf95 -qsuffix=f=f90:cpp=F90 F90=xlf95 -qsuffix=f=f90:cpp=F90
FC=xlf FC=xlf -qsuffix=cpp=F
F77=$(FC) F77=$(FC)
CC=xlc CC=xlc
F90COPT= -O3 F90COPT= -O3
@ -78,8 +78,12 @@ $(.mod).o:
$(F90) $(F90COPT) $(INCDIRS) -c $< $(F90) $(F90COPT) $(INCDIRS) -c $<
.f90.o: .f90.o:
$(F90) $(F90COPT) $(INCDIRS) -c $< $(F90) $(F90COPT) $(INCDIRS) -c $<
.F.o:
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<
.F90.o: .F90.o:
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $< $(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<
.F90$(.mod):
$(F90) $(F90COPT) $(INCDIRS) $(FDEFINES) -c $<

@ -1,6 +1,6 @@
include Make.inc include Make.inc
#PREC=../mld2p4 PREC=../mld2p4
PREC=baseprec #PREC=baseprec
library: library:
( [ -d lib ] || mkdir lib) ( [ -d lib ] || mkdir lib)

@ -121,7 +121,6 @@ subroutine psb_iscatterm(globx, locx, desc_a, info, iroot)
n = psb_cd_get_global_cols(desc_a) n = psb_cd_get_global_cols(desc_a)
call psb_bcast(ictxt,k,root=iiroot) call psb_bcast(ictxt,k,root=iiroot)
! there should be a global check on k here!!! ! there should be a global check on k here!!!
call psb_chkglobvect(m,n,size(globx),iglobx,jglobx,desc_a,info) call psb_chkglobvect(m,n,size(globx),iglobx,jglobx,desc_a,info)

@ -11,8 +11,7 @@ COBJS = avltree.o srcht.o
MPFOBJS = psi_dswapdata.o psi_dswaptran.o psi_iswapdata.o \ MPFOBJS = psi_dswapdata.o psi_dswaptran.o psi_iswapdata.o \
psi_iswaptran.o psi_desc_index.o \ psi_iswaptran.o psi_desc_index.o \
psi_zswapdata.o psi_zswaptran.o psi_zswapdata.o psi_zswaptran.o psi_extrct_dl.o
MPFOBJS2 = psi_extrct_dl.o
LIBDIR = .. LIBDIR = ..
MODDIR = ../modules MODDIR = ../modules
INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I . INCDIRS = -I $(LIBDIR) -I $(MODDIR) -I .
@ -25,7 +24,6 @@ lib: mpfobjs $(FOBJS) $(FOBJS2) $(COBJS) $(MPFOBJS2)
mpfobjs: psi_gthsct.o mpfobjs: psi_gthsct.o
(make $(MPFOBJS) F90="$(MPF90)" FC="$(MPF90)" FCOPT="$(F90COPT)") (make $(MPFOBJS) F90="$(MPF90)" FC="$(MPF90)" FCOPT="$(F90COPT)")
(make $(MPFOBJS2) F90="$(MPF77)" FC="$(MPF77)" FCOPT="$(FCOPT)")
(make $(FOBJS2) F90="$(MPF77)" FC="$(MPF77)" FCOPT="$(FCOPT)") (make $(FOBJS2) F90="$(MPF77)" FC="$(MPF77)" FCOPT="$(FCOPT)")
clean: clean:
/bin/rm -f $(MPFOBJS) $(FOBJS) $(COBJS) $(FOBJS2) $(MPFOBJS2) /bin/rm -f $(MPFOBJS) $(FOBJS) $(COBJS) $(FOBJS2) $(MPFOBJS2)

@ -1,282 +0,0 @@
C
C Parallel Sparse BLAS v2.0
C (C) Copyright 2006 Salvatore Filippone University of Rome Tor Vergata
C Alfredo Buttari University of Rome Tor Vergata
C
C Redistribution and use in source and binary forms, with or without
C modification, are permitted provided that the following conditions
C are met:
C 1. Redistributions of source code must retain the above copyright
C notice, this list of conditions and the following disclaimer.
C 2. Redistributions in binary form must reproduce the above copyright
C notice, this list of conditions, and the following disclaimer in the
C documentation and/or other materials provided with the distribution.
C 3. The name of the PSBLAS group or the names of its contributors may
C not be used to endorse or promote products derived from this
C software without specific written permission.
C
C THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
C ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
C TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
C PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
C BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
C CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
C SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
C INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
C CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
C ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
C POSSIBILITY OF SUCH DAMAGE.
C
C
subroutine psi_extract_dep_list(desc_data,
+ desc_str,dep_list,
+ length_dl,np,dl_lda,mode,info)
c internal routine
c ================
c
c _____called by psi_crea_halo and psi_crea_ovrlap ______
c
c purpose
c =======
c process root (pid=0) extracts for each process "k" the ordered list of process
c to which "k" must communicate. this list with its order is extracted from
c desc_str list
c
c
c input
c =======
c desc_data :integer array
c explanation:
c name explanation
c ------------------ -------------------------------------------------------
c desc_data array of integer that contains some local and global
c information of matrix.
c
c
c now we explain each of the above vectors.
c
c let a be a generic sparse matrix. we denote with matdata_a the matrix_data
c array for matrix a.
c data stored in matrix_data array are:
c
c notation stored in explanation
c --------------- ---------------------- -------------------------------------
c dec_type matdata_a[psb_dec_type_] decomposition type
c m matdata_a[m_] total number of equations
c n matdata_a[n_] total number of variables
c n_row matdata_a[psb_n_row_] number of local equations
c n_col matdata_a[psb_n_col_] number of local variables
c psb_ctxt_a matdata_a[ctxt_] the blacs context handle, indicating
c the global context of the operation
c on the matrix.
c the context itself is global.
c desc_str integer array
c explanation:
c let desc_str_p be the array desc_str for local process.
c this is composed of variable dimension blocks for each process to
c communicate to.
c each block contain indexes of local halo elements to exchange with other
c process.
c let p be the pointer to the first element of a block in desc_str_p.
c this block is stored in desc_str_p as :
c
c notation stored in explanation
c --------------- --------------------------- -----------------------------------
c process_id desc_str_p[p+psb_proc_id_] identifier of process which exchange
c data with.
c n_elements_recv desc_str_p[p+n_elem_recv_] number of elements to receive.
c elements_recv desc_str_p[p+elem_recv_+i] indexes of local elements to
c receive. these are stored in the
c array from location p+elem_recv_ to
c location p+elem_recv_+
c desc_str_p[p+n_elem_recv_]-1.
c if desc_data(psb_dec_type_) == 0
c then also will be:
c n_elements_send desc_str_p[p+n_elem_send_] number of elements to send.
c elements_send desc_str_p[p+elem_send_+i] indexes of local elements to
c send. these are stored in the
c array from location p+elem_send_ to
c location p+elem_send_+
c desc_str_p[p+n_elem_send_]-1.
c list is ended by -1 value
c
c np integer (global input)
c number of grid process.
c
c mode integer (global input)
c if mode =0 then will be inserted also duplicate element in
c a same dependence list
c if mode =1 then not will be inserted duplicate element in
c a same dependence list
c output
c =====
c only for root (pid=0) process:
c dep_list integer array(dl_lda,0:np)
c dependence list dep_list(*,i) is the list of process identifiers to which process i
c must communicate with. this list with its order is extracted from
c desc_str list.
c length_dl integer array(0:np)
c length_dl(i) is the length of dep_list(*,i) list
use psb_penv_mod
use psb_const_mod
use psb_error_mod
use psb_descriptor_type
implicit none
include 'mpif.h'
c ....scalar parameters...
integer np,dl_lda,mode, info
c ....array parameters....
integer desc_str(*),desc_data(*),
+ dep_list(dl_lda,0:np),length_dl(0:np)
integer, pointer :: itmp(:)
c .....local arrays....
integer int_err(5)
double precision real_err(5)
c .....local scalars...
integer i,nprow,npcol,me,mycol,pointer_dep_list,proc,j,err_act
integer ictxt, err, icomm
logical debug
parameter (debug=.false.)
character name*20
name='psi_extrct_dl'
call psb_erractionsave(err_act)
info = 0
ictxt = desc_data(psb_ctxt_)
call psb_info(ictxt,me,nprow)
do i=0,np
length_dl(i) = 0
enddo
i=1
if (debug) write(0,*) 'extract: info ',info,
+ desc_data(psb_dec_type_)
pointer_dep_list=1
c$$$ if (desc_data(psb_dec_type_).eq.psb_desc_bld_) then
if (psb_is_bld_dec(desc_data(psb_dec_type_))) then
do while (desc_str(i).ne.-1)
if (debug) write(0,*) me,' extract: looping ',i,
+ desc_str(i),desc_str(i+1),desc_str(i+2)
c ...with different decomposition type we have different
c structure of indices lists............................
if ((desc_str(i+1).ne.0).or.(desc_str(i+2).ne.0)) then
c ..if number of element to be exchanged !=0
proc=desc_str(i)
if ((proc.lt.0).or.(proc.ge.nprow)) then
if (debug) write(0,*) 'extract error ',i,desc_str(i)
info = 9999
int_err(1) = i
int_err(2) = desc_str(i)
goto 998
endif
! if((me.eq.1).and.(proc.eq.3))write(0,*)'found 3'
if (mode.eq.1) then
c ...search if already exist proc
c in dep_list(*,me)...
j=1
do while ((j.lt.pointer_dep_list).and.
+ (dep_list(j,me).ne.proc))
j=j+1
enddo
if (j.eq.pointer_dep_list) then
c ...if not found.....
dep_list(pointer_dep_list,me)=proc
pointer_dep_list=pointer_dep_list+1
endif
else if (mode.eq.0) then
if (pointer_dep_list.gt.dl_lda) then
info = 4000
goto 998
endif
dep_list(pointer_dep_list,me)=proc
pointer_dep_list=pointer_dep_list+1
endif
endif
i=i+desc_str(i+1)+2
enddo
c$$$ else if (desc_data(psb_dec_type_).eq.psb_desc_upd_) then
else if (psb_is_upd_dec(desc_data(psb_dec_type_))) then
do while (desc_str(i).ne.-1)
if (debug) write(0,*) 'extract: looping ',i,desc_str(i)
c ...with different decomposition type we have different
c structure of indices lists............................
if (desc_str(i+1).ne.0) then
proc=desc_str(i)
c ..if number of element to be exchanged !=0
if (mode.eq.1) then
c ...search if already exist proc....
j=1
do while ((j.lt.pointer_dep_list).and.
+ (dep_list(j,me).ne.proc))
j=j+1
enddo
if (j.eq.pointer_dep_list) then
c ...if not found.....
if (pointer_dep_list.gt.dl_lda) then
info = 4000
goto 998
endif
dep_list(pointer_dep_list,me)=proc
pointer_dep_list=pointer_dep_list+1
endif
else if (mode.eq.0) then
if (pointer_dep_list.gt.dl_lda) then
info = 4000
goto 998
endif
dep_list(pointer_dep_list,me)=proc
pointer_dep_list=pointer_dep_list+1
endif
endif
i=i+desc_str(i+1)+2
enddo
else
write(0,*) 'invalid dec_type',desc_data(psb_dec_type_)
info = 2020
goto 9999
endif
length_dl(me)=pointer_dep_list-1
c ... check for errors...
998 continue
if (debug) write(0,*) 'extract: info ',info
err = info
if (err.ne.0) goto 9999
call psb_sum(ictxt,length_dl(0:np))
call psb_get_mpicomm(ictxt,icomm )
allocate(itmp(dl_lda),stat=info)
if (info /= 0) goto 9999
itmp(1:dl_lda) = dep_list(1:dl_lda,me)
call mpi_allgather(itmp,dl_lda,mpi_integer,
+ dep_list,dl_lda,mpi_integer,icomm,info)
deallocate(itmp)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_errpush(info,name,i_err=int_err)
call psb_erractionrestore(err_act)
if (err_act.eq.psb_act_ret_) then
return
else
call psb_error()
endif
return
end

@ -0,0 +1,276 @@
!!$
!!$ 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.
!!$
!!$
subroutine psi_extract_dep_list(desc_data,desc_str,dep_list,&
& length_dl,np,dl_lda,mode,info)
! internal routine
! ================
!
! _____called by psi_crea_halo and psi_crea_ovrlap ______
!
! purpose
! =======
! process root (pid=0) extracts for each process "k" the ordered list of process
! to which "k" must communicate. this list with its order is extracted from
! desc_str list
!
!
! input
! =======
! desc_data :integer array
! explanation:
! name explanation
! ------------------ -------------------------------------------------------
! desc_data array of integer that contains some local and global
! information of matrix.
!
!
! now we explain each of the above vectors.
!
! let a be a generic sparse matrix. we denote with matdata_a the matrix_data
! array for matrix a.
! data stored in matrix_data array are:
!
! notation stored in explanation
! --------------- ---------------------- -------------------------------------
! dec_type matdata_a[psb_dec_type_] decomposition type
! m matdata_a[m_] total number of equations
! n matdata_a[n_] total number of variables
! n_row matdata_a[psb_n_row_] number of local equations
! n_col matdata_a[psb_n_col_] number of local variables
! psb_ctxt_a matdata_a[ctxt_] the blacs context handle, indicating
! the global context of the operation
! on the matrix.
! the context itself is global.
! desc_str integer array
! explanation:
! let desc_str_p be the array desc_str for local process.
!! this is composed of variable dimension blocks for each process to
! communicate to.
! each block contain indexes of local halo elements to exchange with other
! process.
! let p be the pointer to the first element of a block in desc_str_p.
! this block is stored in desc_str_p as :
!
! notation stored in explanation
! --------------- --------------------------- -----------------------------------
! process_id desc_str_p[p+psb_proc_id_] identifier of process which exchange
! data with.
! n_elements_recv desc_str_p[p+n_elem_recv_] number of elements to receive.
! elements_recv desc_str_p[p+elem_recv_+i] indexes of local elements to
! receive. these are stored in the
! array from location p+elem_recv_ to
! location p+elem_recv_+
! desc_str_p[p+n_elem_recv_]-1.
! if desc_data(psb_dec_type_) == 0
! then also will be:
! n_elements_send desc_str_p[p+n_elem_send_] number of elements to send.
! elements_send desc_str_p[p+elem_send_+i] indexes of local elements to
! send. these are stored in the
! array from location p+elem_send_ to
! location p+elem_send_+
! desc_str_p[p+n_elem_send_]-1.
! list is ended by -1 value
!
! np integer (global input)
! number of grid process.
!
! mode integer (global input)
! if mode =0 then will be inserted also duplicate element in
! a same dependence list
! if mode =1 then not will be inserted duplicate element in
! a same dependence list
! output
! =====
! only for root (pid=0) process:
! dep_list integer array(dl_lda,0:np)
! dependence list dep_list(*,i) is the list of process identifiers to which process i
! must communicate with. this list with its order is extracted from
! desc_str list.
! length_dl integer array(0:np)
! length_dl(i) is the length of dep_list(*,i) list
use mpi
use psb_penv_mod
use psb_const_mod
use psb_error_mod
use psb_descriptor_type
implicit none
! ....scalar parameters...
integer np,dl_lda,mode, info
! ....array parameters....
integer :: desc_str(*),desc_data(*),dep_list(dl_lda,0:np),length_dl(0:np)
integer, pointer :: itmp(:)
! .....local arrays....
integer int_err(5)
double precision real_err(5)
! .....local scalars...
integer i,nprow,npcol,me,mycol,pointer_dep_list,proc,j,err_act
integer ictxt, err, icomm
logical, parameter :: debug=.false.
character name*20
name='psi_extrct_dl'
call psb_erractionsave(err_act)
info = 0
ictxt = desc_data(psb_ctxt_)
call psb_info(ictxt,me,nprow)
do i=0,np
length_dl(i) = 0
enddo
i=1
if (debug) write(0,*) 'extract: info ',info,desc_data(psb_dec_type_)
pointer_dep_list=1
if (psb_is_bld_dec(desc_data(psb_dec_type_))) then
do while (desc_str(i) /= -1)
if (debug) write(0,*) me,' extract: looping ',i,&
& desc_str(i),desc_str(i+1),desc_str(i+2)
! ...with different decomposition type we have different
! structure of indices lists............................
if ((desc_str(i+1) /= 0).or.(desc_str(i+2) /= 0)) then
! ..if number of element to be exchanged !=0
proc=desc_str(i)
if ((proc < 0).or.(proc.ge.nprow)) then
if (debug) write(0,*) 'extract error ',i,desc_str(i)
info = 9999
int_err(1) = i
int_err(2) = desc_str(i)
goto 998
endif
! if((me == 1).and.(proc == 3))write(0,*)'found 3'
if (mode == 1) then
! ...search if already exist proc
! in dep_list(*,me)...
j=1
do while ((j < pointer_dep_list).and.&
& (dep_list(j,me) /= proc))
j=j+1
enddo
if (j == pointer_dep_list) then
! ...if not found.....
dep_list(pointer_dep_list,me)=proc
pointer_dep_list=pointer_dep_list+1
endif
else if (mode == 0) then
if (pointer_dep_list.gt.dl_lda) then
info = 4000
goto 998
endif
dep_list(pointer_dep_list,me)=proc
pointer_dep_list=pointer_dep_list+1
endif
endif
i=i+desc_str(i+1)+2
enddo
else if (psb_is_upd_dec(desc_data(psb_dec_type_))) then
do while (desc_str(i) /= -1)
if (debug) write(0,*) 'extract: looping ',i,desc_str(i)
! ...with different decomposition type we have different
! structure of indices lists............................
if (desc_str(i+1) /= 0) then
proc=desc_str(i)
! ..if number of element to be exchanged !=0
if (mode == 1) then
! ...search if already exist proc....
j=1
do while ((j < pointer_dep_list).and.&
& (dep_list(j,me) /= proc))
j=j+1
enddo
if (j == pointer_dep_list) then
! ...if not found.....
if (pointer_dep_list.gt.dl_lda) then
info = 4000
goto 998
endif
dep_list(pointer_dep_list,me)=proc
pointer_dep_list=pointer_dep_list+1
endif
else if (mode == 0) then
if (pointer_dep_list.gt.dl_lda) then
info = 4000
goto 998
endif
dep_list(pointer_dep_list,me)=proc
pointer_dep_list=pointer_dep_list+1
endif
endif
i=i+desc_str(i+1)+2
enddo
else
write(0,*) 'invalid dec_type',desc_data(psb_dec_type_)
info = 2020
goto 9999
endif
length_dl(me)=pointer_dep_list-1
! ... check for errors...
998 continue
if (debug) write(0,*) 'extract: info ',info
err = info
if (err /= 0) goto 9999
call psb_sum(ictxt,length_dl(0:np))
call psb_get_mpicomm(ictxt,icomm )
allocate(itmp(dl_lda),stat=info)
if (info /= 0) goto 9999
itmp(1:dl_lda) = dep_list(1:dl_lda,me)
call mpi_allgather(itmp,dl_lda,mpi_integer,&
& dep_list,dl_lda,mpi_integer,icomm,info)
deallocate(itmp)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_errpush(info,name,i_err=int_err)
call psb_erractionrestore(err_act)
if (err_act == psb_act_ret_) then
return
else
call psb_error()
endif
return
end subroutine psi_extract_dep_list

@ -12,7 +12,7 @@ MODULES = psb_realloc_mod.o psb_string_mod.o psb_spmat_type.o \
LIBMOD=psb_base_mod$(.mod) LIBMOD=psb_base_mod$(.mod)
MPFOBJS=psb_penv_mod.o MPFOBJS=psb_penv_mod.o
OBJS = error.o psb_base_mod.o OBJS = error.o psb_base_mod.o $(EXTRA_COBJS)
LIBDIR = .. LIBDIR = ..
INCDIRS = -I . INCDIRS = -I .

@ -23,15 +23,23 @@ subroutine psb_restore_coher(ictxt,isvch)
end subroutine psb_restore_coher end subroutine psb_restore_coher
subroutine psb_get_mpicomm(ictxt,comm) subroutine psb_get_mpicomm(ictxt,comm)
integer :: ictxt, comm integer :: ictxt, comm
#if !defined(SERIAL_MPI)
call blacs_get(ictxt,10,comm) call blacs_get(ictxt,10,comm)
#else
comm = ictxt
#endif
end subroutine psb_get_mpicomm end subroutine psb_get_mpicomm
subroutine psb_get_rank(rank,ictxt,id) subroutine psb_get_rank(rank,ictxt,id)
integer :: rank,ictxt, id integer :: rank,ictxt, id
integer :: blacs_pnum integer :: blacs_pnum
#if defined(SERIAL_MPI)
rank = 0
#else
rank = blacs_pnum(ictxt,id,0) rank = blacs_pnum(ictxt,id,0)
#endif
end subroutine psb_get_rank end subroutine psb_get_rank
#ifdef ESSL_BLACS #if defined(ESSL_BLACS) || defined(SERIAL_MPI)
! !
! Need these, as they are not in the ESSL implementation ! Need these, as they are not in the ESSL implementation
! of the BLACS. ! of the BLACS.

@ -101,7 +101,7 @@ subroutine FCpsb_get_errverbosity(v)
integer, intent(out) :: v integer, intent(out) :: v
call psb_get_errverbosity(v) v = psb_get_errverbosity()
end subroutine FCpsb_get_errverbosity end subroutine FCpsb_get_errverbosity

@ -0,0 +1,134 @@
#include <sys/time.h>
#include <stdio.h>
#include <string.h>
#ifdef Add_
#define mpi_wtime mpi_wtime_
#define mpi_send mpi_send_
#define mpi_irecv mpi_irecv_
#define mpi_wait mpi_wait_
#define mpi_alltoall mpi_alltoall_
#define mpi_alltoallv mpi_alltoallv_
#define mpi_allgather mpi_allgather_
#define mpi_allgatherv mpi_allgatherv_
#endif
#define mpi_integer 1
#define mpi_double 3
#define mpi_double_complex 5
double mpi_wtime()
{
struct timeval tt;
struct timezone tz;
double temp;
if (gettimeofday(&tt,&tz) != 0) {
fprintf(stderr,"Fatal error for gettimeofday ??? \n");
temp=0.0;
} else {
temp = ((double)tt.tv_sec) + ((double)tt.tv_usec)*1.0e-6;
}
return(temp);
}
void mpi_wait()
{
return;
}
void mpi_send()
{
return;
}
void mpi_irecv()
{
return;
}
void mpi_alltoall(void* sdb, int* sdc, int* sdt,
void* rvb, int* rvc, int* rvt, int* comm, int* ierr)
{
int i,j,k;
if (*sdt == mpi_integer) {
memcpy(rvb,sdb, (*sdc)*sizeof(int));
}
if (*sdt == mpi_double) {
memcpy(rvb,sdb, (*sdc)*sizeof(double));
}
if (*sdt == mpi_double_complex) {
memcpy(rvb,sdb, (*sdc)*2*sizeof(double));
}
*ierr = 0;
}
void mpi_alltoallv(void* sdb, int* sdc, int* sdspl, int* sdt,
void* rvb, int* rvc, int* rdspl, int* rvt, int* comm, int* ierr)
{
int i,j,k;
if (*sdt == mpi_integer) {
memcpy((rvb+rdspl[0]*sizeof(int)),
(sdb+sdspl[0]*sizeof(int)),(*sdc)*sizeof(int));
}
if (*sdt == mpi_double) {
memcpy((rvb+rdspl[0]*sizeof(double)),
(sdb+sdspl[0]*sizeof(double)),(*sdc)*sizeof(double));
}
if (*sdt == mpi_double_complex) {
memcpy((rvb+rdspl[0]*2*sizeof(double)),
(sdb+sdspl[0]*2*sizeof(double)),(*sdc)*2*sizeof(double));
}
*ierr = 0;
}
void mpi_allgather(void* sdb, int* sdc, int* sdt,
void* rvb, int* rvc, int* rvt, int* comm, int* ierr)
{
int i,j,k;
if (*sdt == mpi_integer) {
memcpy(rvb,sdb, (*sdc)*sizeof(int));
}
if (*sdt == mpi_double) {
memcpy(rvb,sdb, (*sdc)*sizeof(double));
}
if (*sdt == mpi_double_complex) {
memcpy(rvb,sdb, (*sdc)*2*sizeof(double));
}
*ierr = 0;
}
void mpi_allgatherv(void* sdb, int* sdc, int* sdt,
void* rvb, int* rvc, int* rdspl,
int* rvt, int* comm, int* ierr)
{
int i,j,k;
if (*sdt == mpi_integer) {
memcpy(rvb,sdb, (*sdc)*sizeof(int));
}
if (*sdt == mpi_double) {
memcpy(rvb,sdb, (*sdc)*sizeof(double));
}
if (*sdt == mpi_double_complex) {
memcpy(rvb,sdb, (*sdc)*2*sizeof(double));
}
if (*sdt == mpi_integer) {
memcpy((rvb+rdspl[0]*sizeof(int)),
(sdb),(*sdc)*sizeof(int));
}
if (*sdt == mpi_double) {
memcpy((rvb+rdspl[0]*sizeof(double)),
(sdb),(*sdc)*sizeof(double));
}
if (*sdt == mpi_double_complex) {
memcpy((rvb+rdspl[0]*2*sizeof(double)),
(sdb),(*sdc)*2*sizeof(double));
}
*ierr = 0;
}

@ -112,8 +112,10 @@ contains
integer :: temp(2) integer :: temp(2)
integer, parameter :: ione=1 integer, parameter :: ione=1
! Cannot use psb_amx or otherwise we have a recursion in module usage ! Cannot use psb_amx or otherwise we have a recursion in module usage
#if !defined(SERIAL_MPI)
call igamx2d(ictxt, 'A', ' ', ione, ione, err, ione,& call igamx2d(ictxt, 'A', ' ', ione, ione, err, ione,&
&temp ,temp,-ione ,-ione,-ione) &temp ,temp,-ione ,-ione,-ione)
#endif
end subroutine psb_errcomm end subroutine psb_errcomm
@ -127,10 +129,11 @@ contains
! returns verbosity of the error message ! returns verbosity of the error message
subroutine psb_get_errverbosity(v) function psb_get_errverbosity()
integer, intent(out) :: v integer :: psb_get_errverbosity
v=verbosity_level
end subroutine psb_get_errverbosity psb_get_errverbosity=verbosity_level
end function psb_get_errverbosity
@ -203,7 +206,6 @@ contains
integer :: nprow, npcol, me, mypcol integer :: nprow, npcol, me, mypcol
integer, parameter :: ione=1, izero=0 integer, parameter :: ione=1, izero=0
call blacs_gridinfo(ictxt, nprow, npcol, me, mypcol)
if(error_status.gt.0) then if(error_status.gt.0) then
if(verbosity_level.gt.1) then if(verbosity_level.gt.1) then
@ -214,7 +216,11 @@ contains
call psb_errmsg(err_c, r_name, i_e_d, a_e_d,me) call psb_errmsg(err_c, r_name, i_e_d, a_e_d,me)
! write(0,'(50("="))') ! write(0,'(50("="))')
end do end do
#if defined(SERIAL_MPI)
stop
#else
call blacs_abort(ictxt,-1) call blacs_abort(ictxt,-1)
#endif
else else
call psb_errpop(err_c, r_name, i_e_d, a_e_d) call psb_errpop(err_c, r_name, i_e_d, a_e_d)
@ -222,12 +228,20 @@ contains
do while (error_stack%n_elems.gt.0) do while (error_stack%n_elems.gt.0)
call psb_errpop(err_c, r_name, i_e_d, a_e_d) call psb_errpop(err_c, r_name, i_e_d, a_e_d)
end do end do
#if defined(SERIAL_MPI)
stop
#else
call blacs_abort(ictxt,-1) call blacs_abort(ictxt,-1)
#endif
end if end if
end if end if
if(error_status.gt.izero) then if(error_status.gt.izero) then
#if defined(SERIAL_MPI)
stop
#else
call blacs_abort(ictxt,err_c) call blacs_abort(ictxt,err_c)
#endif
end if end if

@ -100,8 +100,10 @@ C
DOUBLE PRECISION HIS( 2 ) DOUBLE PRECISION HIS( 2 )
* .. * ..
* .. External Subroutines .. * .. External Subroutines ..
#if !defined(SERIAL_MPI)
EXTERNAL BLACS_GRIDINFO, DGEBR2D, DGEBS2D, EXTERNAL BLACS_GRIDINFO, DGEBR2D, DGEBS2D,
$ DGERV2D, DGESD2D $ DGERV2D, DGESD2D
#endif
* .. * ..
* .. External Functions .. * .. External Functions ..
LOGICAL LSAME LOGICAL LSAME
@ -122,6 +124,8 @@ C
TRDEST = RDEST0 TRDEST = RDEST0
TCDEST = CDEST0 TCDEST = CDEST0
END IF END IF
#if !defined(SERIAL_MPI)
* *
* Get grid parameters. * Get grid parameters.
* *
@ -234,6 +238,7 @@ C
$ TRDEST, TCDEST ) $ TRDEST, TCDEST )
END IF END IF
END IF END IF
#endif
* *
RETURN RETURN
* *

@ -70,11 +70,17 @@ subroutine psb_cd_inloc(v, ictxt, desc_a, info)
if (debug) write(*,*) 'psb_cdall: ',np,me if (debug) write(*,*) 'psb_cdall: ',np,me
if (.false.) then
loc_row = size(v) loc_row = size(v)
m = loc_row m = loc_row
call psb_sum(ictxt,m) call psb_sum(ictxt,m)
else
m = maxval(v)
call psb_max(ictxt,m)
end if
n = m n = m
!... check m and n parameters.... !... check m and n parameters....
if (m < 1) then if (m < 1) then
info = 10 info = 10

@ -37,6 +37,7 @@
! desc_a - type(<psb_desc_type>). The communication descriptor. ! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code. ! info - integer. Eventually returns an error code.
subroutine psb_icdasb(desc_a,info,ext_hv) subroutine psb_icdasb(desc_a,info,ext_hv)
use mpi
use psb_descriptor_type use psb_descriptor_type
use psb_serial_mod use psb_serial_mod
use psb_const_mod use psb_const_mod
@ -44,8 +45,6 @@ subroutine psb_icdasb(desc_a,info,ext_hv)
use psb_error_mod use psb_error_mod
use psb_penv_mod use psb_penv_mod
implicit none implicit none
include 'mpif.h'
!...Parameters.... !...Parameters....
type(psb_desc_type), intent(inout) :: desc_a type(psb_desc_type), intent(inout) :: desc_a
integer, intent(out) :: info integer, intent(out) :: info

@ -121,14 +121,14 @@ Subroutine psb_zsphalo(a,desc_a,blk,info,rwcnv,clcnv,outfmt)
sdsz(:)=0 sdsz(:)=0
rvsz(:)=0 rvsz(:)=0
ipx = 1 ipx = 1
blk%k = a%k
blk%m = 0
brvindx(ipx) = 0 brvindx(ipx) = 0
bsdindx(ipx) = 0 bsdindx(ipx) = 0
counter=1 counter=1
idx = 0 idx = 0
idxs = 0 idxs = 0
idxr = 0 idxr = 0
blk%k = a%k
blk%m = 0
! For all rows in the halo descriptor, extract and send/receive. ! For all rows in the halo descriptor, extract and send/receive.
Do Do
proc=desc_a%halo_index(counter) proc=desc_a%halo_index(counter)

@ -52,7 +52,7 @@ subroutine psb_zbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
character ::diagl, diagu character ::diagl, diagu
integer :: ictxt,np,me,i, isz, nrg, err_act integer :: ictxt,np,me,i, isz, nrg, err_act
real(kind(1.d0)) :: t1, t2, t3, t4, t5, t6, t7 real(kind(1.d0)) :: t1, t2, t3, t4, t5, t6, t7
logical,parameter :: debug=.false., debugprt=. logical,parameter :: debug=.false., debugprt=.false.
character(len=20) :: name, ch_err character(len=20) :: name, ch_err
interface psb_bjac_aply interface psb_bjac_aply

@ -1,5 +1,5 @@
10 Number of inputs 10 Number of inputs
young1r.mtx This (and others) from: http://math.nist.gov/MatrixMarket/ or lapl600x600.mtx young1r.mtx This (and others) from: http://math.nist.gov/MatrixMarket/ or
NONE http://www.cise.ufl.edu/research/sparse/matrices/index.html NONE http://www.cise.ufl.edu/research/sparse/matrices/index.html
BICGSTAB BICGSTAB
CSR CSR
@ -7,6 +7,6 @@ CSR
1 ISTOPC 1 ISTOPC
00800 ITMAX 00800 ITMAX
-1 ITRACE -1 ITRACE
4 IPREC 0:NONE 1:DIAGSC 2:ILU 3: AS 4: RAS 5: ASH 6: RASH 2 IPREC 0:NONE 1:DIAGSC 2:ILU
1 ML 1 ML
1.d-6 EPS 1.d-6 EPS

Loading…
Cancel
Save