base/newserial/psbn_d_base_mat_mod.f03
 base/newserial/psbn_d_csr_mat_mod.f03
 prec/psb_ddiagsc_bld.f90
 prec/psb_dprecbld.f90
 prec/psb_dprecinit.f90
 prec/psb_prec_mod.f90
 test/pargen/runs/ppde.inp
 test/serial/psbn_d_cxx_impl.f03
 test/serial/psbn_d_cxx_mat_mod.f03
 util/psb_mat_dist_mod.f90

Moved csgetblk and csclip to the base level, the only specific
implementation is that of csgetrow.
psblas3-type-indexed
Salvatore Filippone 16 years ago
parent 9bfb2980e3
commit ffe5ab739d

@ -71,8 +71,6 @@ module psbn_d_base_mat_mod
procedure, pass(a) :: free => d_coo_free procedure, pass(a) :: free => d_coo_free
procedure, pass(a) :: trim => d_coo_trim procedure, pass(a) :: trim => d_coo_trim
procedure, pass(a) :: d_csgetrow => d_coo_csgetrow procedure, pass(a) :: d_csgetrow => d_coo_csgetrow
procedure, pass(a) :: d_csgetblk => d_coo_csgetblk
procedure, pass(a) :: csclip => d_coo_csclip
procedure, pass(a) :: print => d_coo_print procedure, pass(a) :: print => d_coo_print
procedure, pass(a) :: get_fmt => d_coo_get_fmt procedure, pass(a) :: get_fmt => d_coo_get_fmt
@ -84,8 +82,7 @@ module psbn_d_base_mat_mod
& d_fix_coo, d_coo_free, d_coo_print, d_coo_get_fmt, & & d_fix_coo, d_coo_free, d_coo_print, d_coo_get_fmt, &
& d_cp_coo_to_coo, d_cp_coo_from_coo, & & d_cp_coo_to_coo, d_cp_coo_from_coo, &
& d_cp_coo_to_fmt, d_cp_coo_from_fmt, & & d_cp_coo_to_fmt, d_cp_coo_from_fmt, &
& d_coo_scals, d_coo_scal, d_coo_csgetrow, d_coo_csgetblk, & & d_coo_scals, d_coo_scal, d_coo_csgetrow
& d_coo_csclip
interface interface
@ -568,6 +565,7 @@ contains
end subroutine d_csgetrow end subroutine d_csgetrow
subroutine d_csgetblk(imin,imax,a,b,info,& subroutine d_csgetblk(imin,imax,a,b,info,&
& jmin,jmax,iren,append,rscale,cscale) & jmin,jmax,iren,append,rscale,cscale)
! Output is always in COO format ! Output is always in COO format
@ -583,19 +581,44 @@ contains
integer, intent(in), optional :: iren(:) integer, intent(in), optional :: iren(:)
integer, intent(in), optional :: jmin,jmax integer, intent(in), optional :: jmin,jmax
logical, intent(in), optional :: rscale,cscale logical, intent(in), optional :: rscale,cscale
Integer :: err_act Integer :: err_act, nzin, nzout
character(len=20) :: name='csget' character(len=20) :: name='csget'
logical :: append_
logical, parameter :: debug=.false. logical, parameter :: debug=.false.
call psb_get_erraction(err_act) call psb_erractionsave(err_act)
! This is the base version. If we get here info = 0
! it means the derived class is incomplete,
! so we throw an error.
info = 700
call psb_errpush(info,name,a_err=a%get_fmt())
if (err_act /= psb_act_ret_) then if (present(append)) then
append_ = append
else
append_ = .false.
endif
if (append_) then
nzin = a%get_nzeros()
else
nzin = 0
endif
call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
& jmin=jmin, jmax=jmax, iren=iren, append=append_, &
& nzin=nzin, rscale=rscale, cscale=cscale)
if (info /= 0) goto 9999
call b%set_nzeros(nzin+nzout)
call b%fix(info)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error() call psb_error()
return
end if end if
return return
@ -614,19 +637,78 @@ contains
integer,intent(out) :: info integer,intent(out) :: info
integer, intent(in), optional :: imin,imax,jmin,jmax integer, intent(in), optional :: imin,imax,jmin,jmax
logical, intent(in), optional :: rscale,cscale logical, intent(in), optional :: rscale,cscale
Integer :: err_act
character(len=20) :: name='csclip' Integer :: err_act, nzin, nzout, imin_, imax_, jmin_, jmax_, mb,nb
character(len=20) :: name='csget'
logical :: rscale_, cscale_
logical, parameter :: debug=.false. logical, parameter :: debug=.false.
call psb_get_erraction(err_act) call psb_erractionsave(err_act)
! This is the base version. If we get here info = 0
! it means the derived class is incomplete,
! so we throw an error.
info = 700
call psb_errpush(info,name,a_err=a%get_fmt())
if (err_act /= psb_act_ret_) then nzin = 0
if (present(imin)) then
imin_ = imin
else
imin_ = 1
end if
if (present(imax)) then
imax_ = imax
else
imax_ = a%get_nrows()
end if
if (present(jmin)) then
jmin_ = jmin
else
jmin_ = 1
end if
if (present(jmax)) then
jmax_ = jmax
else
jmax_ = a%get_ncols()
end if
if (present(rscale)) then
rscale_ = rscale
else
rscale_ = .true.
end if
if (present(cscale)) then
cscale_ = cscale
else
cscale_ = .true.
end if
if (rscale_) then
mb = imax_ - imin_ +1
else
mb = a%get_nrows() ! Should this be imax_ ??
endif
if (cscale_) then
nb = jmax_ - jmin_ +1
else
nb = a%get_ncols() ! Should this be jmax_ ??
endif
call b%allocate(mb,nb)
call a%csget(imin_,imax_,nzout,b%ia,b%ja,b%val,info,&
& jmin=jmin_, jmax=jmax_, append=.false., &
& nzin=nzin, rscale=rscale_, cscale=cscale_)
if (info /= 0) goto 9999
call b%set_nzeros(nzin+nzout)
call b%fix(info)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error() call psb_error()
return
end if end if
return return
@ -1392,156 +1474,6 @@ contains
end subroutine d_coo_csgetrow end subroutine d_coo_csgetrow
subroutine d_coo_csgetblk(imin,imax,a,b,info,&
& jmin,jmax,iren,append,rscale,cscale)
! Output is always in COO format
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_coo_sparse_mat), intent(in) :: a
class(psbn_d_coo_sparse_mat), intent(inout) :: b
integer, intent(in) :: imin,imax
integer,intent(out) :: info
logical, intent(in), optional :: append
integer, intent(in), optional :: iren(:)
integer, intent(in), optional :: jmin,jmax
logical, intent(in), optional :: rscale,cscale
Integer :: err_act, nzin, nzout
character(len=20) :: name='csget'
logical :: append_
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
info = 0
if (present(append)) then
append_ = append
else
append_ = .false.
endif
if (append_) then
nzin = a%get_nzeros()
else
nzin = 0
endif
call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
& jmin=jmin, jmax=jmax, iren=iren, append=append_, &
& nzin=nzin, rscale=rscale, cscale=cscale)
if (info /= 0) goto 9999
call b%set_nzeros(nzin+nzout)
call b%fix(info)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_coo_csgetblk
subroutine d_coo_csclip(a,b,info,&
& imin,imax,jmin,jmax,rscale,cscale)
! Output is always in COO format
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_coo_sparse_mat), intent(in) :: a
class(psbn_d_coo_sparse_mat), intent(out) :: b
integer,intent(out) :: info
integer, intent(in), optional :: imin,imax,jmin,jmax
logical, intent(in), optional :: rscale,cscale
Integer :: err_act, nzin, nzout, imin_, imax_, jmin_, jmax_, mb,nb
character(len=20) :: name='csget'
logical :: rscale_, cscale_
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
info = 0
nzin = 0
if (present(imin)) then
imin_ = imin
else
imin_ = 1
end if
if (present(imax)) then
imax_ = imax
else
imax_ = a%get_nrows()
end if
if (present(jmin)) then
jmin_ = jmin
else
jmin_ = 1
end if
if (present(jmax)) then
jmax_ = jmax
else
jmax_ = a%get_ncols()
end if
if (present(rscale)) then
rscale_ = rscale
else
rscale_ = .true.
end if
if (present(cscale)) then
cscale_ = cscale
else
cscale_ = .true.
end if
if (rscale_) then
mb = imax_ - imin_ +1
else
mb = a%get_nrows() ! Should this be imax_ ??
endif
if (cscale_) then
nb = jmax_ - jmin_ +1
else
nb = a%get_ncols() ! Should this be jmax_ ??
endif
call b%allocate(mb,nb)
call a%csget(imin_,imax_,nzout,b%ia,b%ja,b%val,info,&
& jmin=jmin_, jmax=jmax_, append=.false., &
& nzin=nzin, rscale=rscale_, cscale=cscale_)
if (info /= 0) goto 9999
call b%set_nzeros(nzin+nzout)
call b%fix(info)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_coo_csclip
!!$ !!$
!!$ subroutine d_coo_csget(irw,a,nz,ia,ja,val,info,iren,lrw,append,nzin) !!$ subroutine d_coo_csget(irw,a,nz,ia,ja,val,info,iren,lrw,append,nzin)
!!$ ! Output is always in COO format !!$ ! Output is always in COO format

@ -30,8 +30,6 @@ module psbn_d_csr_mat_mod
procedure, pass(a) :: mv_to_fmt => d_mv_csr_to_fmt procedure, pass(a) :: mv_to_fmt => d_mv_csr_to_fmt
procedure, pass(a) :: mv_from_fmt => d_mv_csr_from_fmt procedure, pass(a) :: mv_from_fmt => d_mv_csr_from_fmt
procedure, pass(a) :: d_csgetrow => d_csr_csgetrow procedure, pass(a) :: d_csgetrow => d_csr_csgetrow
procedure, pass(a) :: d_csgetblk => d_csr_csgetblk
procedure, pass(a) :: csclip => d_csr_csclip
procedure, pass(a) :: get_size => d_csr_get_size procedure, pass(a) :: get_size => d_csr_get_size
procedure, pass(a) :: free => d_csr_free procedure, pass(a) :: free => d_csr_free
procedure, pass(a) :: trim => d_csr_trim procedure, pass(a) :: trim => d_csr_trim
@ -44,8 +42,7 @@ module psbn_d_csr_mat_mod
& d_mv_csr_to_coo, d_mv_csr_from_coo, & & d_mv_csr_to_coo, d_mv_csr_from_coo, &
& d_cp_csr_to_fmt, d_cp_csr_from_fmt, & & d_cp_csr_to_fmt, d_cp_csr_from_fmt, &
& d_mv_csr_to_fmt, d_mv_csr_from_fmt, & & d_mv_csr_to_fmt, d_mv_csr_from_fmt, &
& d_csr_scals, d_csr_scal, d_csr_trim, d_csr_csgetrow, d_csr_csgetblk, & & d_csr_scals, d_csr_scal, d_csr_trim, d_csr_csgetrow, d_csr_get_size
& d_csr_csclip, d_csr_get_size
interface interface

@ -32,10 +32,11 @@
subroutine psb_ddiagsc_bld(a,desc_a,p,upd,info) subroutine psb_ddiagsc_bld(a,desc_a,p,upd,info)
use psb_base_mod use psb_base_mod
use psbn_d_mat_mod
use psb_prec_mod, psb_protect_name => psb_ddiagsc_bld use psb_prec_mod, psb_protect_name => psb_ddiagsc_bld
Implicit None Implicit None
type(psb_dspmat_type), intent(in), target :: a type(psbn_d_sparse_mat), intent(in), target :: a
type(psb_desc_type), intent(in) :: desc_a type(psb_desc_type), intent(in) :: desc_a
type(psb_dprec_type),intent(inout) :: p type(psb_dprec_type),intent(inout) :: p
character, intent(in) :: upd character, intent(in) :: upd
@ -76,7 +77,7 @@ subroutine psb_ddiagsc_bld(a,desc_a,p,upd,info)
! !
! Retrieve the diagonal entries of the matrix A ! Retrieve the diagonal entries of the matrix A
! !
call psb_sp_getdiag(a,p%d,info) call a%get_diag(p%d,info)
if(info /= 0) then if(info /= 0) then
info=4010 info=4010
ch_err='psb_sp_getdiag' ch_err='psb_sp_getdiag'
@ -105,18 +106,18 @@ subroutine psb_ddiagsc_bld(a,desc_a,p,upd,info)
endif endif
end do end do
if (a%pl(1) /= 0) then !!$ if (a%pl(1) /= 0) then
! !!$ !
! Apply the same row permutation as in the sparse matrix A !!$ ! Apply the same row permutation as in the sparse matrix A
! !!$ !
call psb_gelp('n',a%pl,p%d,info) !!$ call psb_gelp('n',a%pl,p%d,info)
if(info /= 0) then !!$ if(info /= 0) then
info=4010 !!$ info=4010
ch_err='psb_dgelp' !!$ ch_err='psb_dgelp'
call psb_errpush(info,name,a_err=ch_err) !!$ call psb_errpush(info,name,a_err=ch_err)
goto 9999 !!$ goto 9999
end if !!$ end if
endif !!$ endif
call psb_erractionrestore(err_act) call psb_erractionrestore(err_act)
return return

@ -102,7 +102,7 @@ subroutine psb_dprecbld(aa,desc_a,p,info,upd)
case (psb_diag_) case (psb_diag_)
call psb_diagsc_bld(a,desc_a,p,upd_,info) call psb_diagsc_bld(aa,desc_a,p,upd_,info)
if(info /= 0) then if(info /= 0) then
info=4010 info=4010
ch_err='psb_diagsc_bld' ch_err='psb_diagsc_bld'

@ -43,12 +43,7 @@ subroutine psb_dprecinit(p,ptype,info)
call psb_realloc(psb_ifpsz,p%iprcparm,info) call psb_realloc(psb_ifpsz,p%iprcparm,info)
if (info == 0) call psb_realloc(psb_rfpsz,p%rprcparm,info) if (info == 0) call psb_realloc(psb_rfpsz,p%rprcparm,info)
if (info /= 0) return if (info /= 0) return
p%iprcparm(:) = 0
p%iprcparm(:) = 0
p%iprcparm(psb_p_type_) = psb_noprec_
p%iprcparm(psb_f_type_) = psb_f_none_
return
select case(psb_toupper(ptype(1:len_trim(ptype)))) select case(psb_toupper(ptype(1:len_trim(ptype))))
case ('NONE','NOPREC') case ('NONE','NOPREC')
p%iprcparm(:) = 0 p%iprcparm(:) = 0
@ -60,12 +55,12 @@ subroutine psb_dprecinit(p,ptype,info)
p%iprcparm(psb_p_type_) = psb_diag_ p%iprcparm(psb_p_type_) = psb_diag_
p%iprcparm(psb_f_type_) = psb_f_none_ p%iprcparm(psb_f_type_) = psb_f_none_
case ('BJAC') !!$ case ('BJAC')
p%iprcparm(:) = 0 !!$ p%iprcparm(:) = 0
p%iprcparm(psb_p_type_) = psb_bjac_ !!$ p%iprcparm(psb_p_type_) = psb_bjac_
p%iprcparm(psb_f_type_) = psb_f_ilu_n_ !!$ p%iprcparm(psb_f_type_) = psb_f_ilu_n_
p%iprcparm(psb_ilu_fill_in_) = 0 !!$ p%iprcparm(psb_ilu_fill_in_) = 0
!!$
case default case default
write(0,*) 'Unknown preconditioner type request "',ptype,'"' write(0,*) 'Unknown preconditioner type request "',ptype,'"'
info = 2 info = 2

@ -406,8 +406,9 @@ module psb_prec_mod
subroutine psb_ddiagsc_bld(a,desc_a,p,upd,info) subroutine psb_ddiagsc_bld(a,desc_a,p,upd,info)
use psb_base_mod, only : psb_desc_type, psb_dspmat_type, psb_dpk_ use psb_base_mod, only : psb_desc_type, psb_dspmat_type, psb_dpk_
use psb_prec_type, only : psb_dprec_type use psb_prec_type, only : psb_dprec_type
use psbn_d_mat_mod
integer, intent(out) :: info integer, intent(out) :: info
type(psb_dspmat_type), intent(in), target :: a type(psbn_d_sparse_mat), intent(in), target :: a
type(psb_dprec_type), intent(inout) :: p type(psb_dprec_type), intent(inout) :: p
type(psb_desc_type), intent(in) :: desc_a type(psb_desc_type), intent(in) :: desc_a
character, intent(in) :: upd character, intent(in) :: upd

@ -1,11 +1,11 @@
7 Number of entries below this 7 Number of entries below this
BICGSTAB Iterative method BICGSTAB CGS BICG BICGSTABL RGMRES BICGSTAB Iterative method BICGSTAB CGS BICG BICGSTABL RGMRES
NONE Preconditioner NONE DIAG BJAC DIAG Preconditioner NONE DIAG BJAC
CSR Storage format for matrix A: CSR COO JAD CSR Storage format for matrix A: CSR COO JAD
050 Domain size (acutal system is this**3) 060 Domain size (acutal system is this**3)
1 Stopping criterion 1 Stopping criterion
1000 MAXIT 0400 MAXIT
040 ITRACE -40 ITRACE
20 IRST restart for RGMRES and BiCGSTABL 20 IRST restart for RGMRES and BiCGSTABL

@ -1,16 +1,16 @@
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !=====================================
! !
! !
! !
! Computational routines ! Computational routines
! !
! !
! !
! !
! !
! !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !=====================================
subroutine d_cxx_csmv_impl(alpha,a,x,beta,y,info,trans) subroutine d_cxx_csmv_impl(alpha,a,x,beta,y,info,trans)
use psb_error_mod use psb_error_mod
@ -980,21 +980,199 @@ function d_cxx_csnmi_impl(a) result(res)
end function d_cxx_csnmi_impl end function d_cxx_csnmi_impl
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !=====================================
! !
! !
! !
! Data management ! Data management
! !
! !
! !
! !
! !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !=====================================
subroutine d_cxx_csgetrow_impl(imin,imax,a,nz,ia,ja,val,info,&
& jmin,jmax,iren,append,nzin,rscale,cscale)
! Output is always in COO format
use psb_error_mod
use psb_const_mod
use psb_error_mod
use psbn_d_base_mat_mod
use psbn_d_cxx_mat_mod, psb_protect_name => d_cxx_csgetrow_impl
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
integer, intent(in) :: imin,imax
integer, intent(out) :: nz
integer, allocatable, intent(inout) :: ia(:), ja(:)
real(psb_dpk_), allocatable, intent(inout) :: val(:)
integer,intent(out) :: info
logical, intent(in), optional :: append
integer, intent(in), optional :: iren(:)
integer, intent(in), optional :: jmin,jmax, nzin
logical, intent(in), optional :: rscale,cscale
logical :: append_, rscale_, cscale_
integer :: nzin_, jmin_, jmax_, err_act, i
character(len=20) :: name='csget'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
info = 0
if (present(jmin)) then
jmin_ = jmin
else
jmin_ = 1
endif
if (present(jmax)) then
jmax_ = jmax
else
jmax_ = a%get_ncols()
endif
if ((imax<imin).or.(jmax_<jmin_)) return
if (present(append)) then
append_=append
else
append_=.false.
endif
if ((append_).and.(present(nzin))) then
nzin_ = nzin
else
nzin_ = 0
endif
if (present(rscale)) then
rscale_ = rscale
else
rscale_ = .false.
endif
if (present(cscale)) then
cscale_ = cscale
else
cscale_ = .false.
endif
if ((rscale_.or.cscale_).and.(present(iren))) then
info = 583
call psb_errpush(info,name,a_err='iren (rscale.or.cscale)')
goto 9999
end if
call cxx_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,info,&
& iren)
if (rscale_) then
do i=nzin_+1, nzin_+nz
ia(i) = ia(i) - imin + 1
end do
end if
if (cscale_) then
do i=nzin_+1, nzin_+nz
ja(i) = ja(i) - jmin_ + 1
end do
end if
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
subroutine d_cxx_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
contains
subroutine cxx_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,info,&
& iren)
use psb_const_mod
use psb_error_mod
use psb_realloc_mod
use psb_sort_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
integer :: imin,imax,jmin,jmax
integer, intent(out) :: nz
integer, allocatable, intent(inout) :: ia(:), ja(:)
real(psb_dpk_), allocatable, intent(inout) :: val(:)
integer, intent(in) :: nzin
logical, intent(in) :: append
integer :: info
integer, optional :: iren(:)
integer :: nzin_, nza, idx,i,j,k, nzt, irw, lrw
integer :: debug_level, debug_unit
character(len=20) :: name='coo_getrow'
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
nza = a%get_nzeros()
irw = imin
lrw = min(imax,a%get_nrows())
if (irw<0) then
write(debug_unit,*) ' spgtrow Error : idx no good ',irw
info = 2
return
end if
if (append) then
nzin_ = nzin
else
nzin_ = 0
endif
nzt = a%irp(lrw)-a%irp(irw)
nz = 0
call psb_ensure_size(nzin_+nzt,ia,info)
if (info==0) call psb_ensure_size(nzin_+nzt,ja,info)
if (info==0) call psb_ensure_size(nzin_+nzt,val,info)
if (info /= 0) return
if (present(iren)) then
do i=irw, lrw
do j=a%irp(i), a%irp(i+1) - 1
if ((jmin <= a%ja(j)).and.(a%ja(j)<=jmax)) then
nzin_ = nzin_ + 1
nz = nz + 1
val(nzin_) = a%val(j)
ia(nzin_) = iren(i)
ja(nzin_) = iren(a%ja(j))
end if
enddo
end do
else
do i=irw, lrw
do j=a%irp(i), a%irp(i+1) - 1
if ((jmin <= a%ja(j)).and.(a%ja(j)<=jmax)) then
nzin_ = nzin_ + 1
nz = nz + 1
val(nzin_) = a%val(j)
ia(nzin_) = (i)
ja(nzin_) = (a%ja(j))
end if
enddo
end do
end if
end subroutine cxx_getrow
end subroutine d_cxx_csgetrow_impl
subroutine d_cxx_csput_impl(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl)
use psb_error_mod use psb_error_mod
use psb_realloc_mod use psb_realloc_mod
use psbn_d_cxx_mat_mod, psb_protect_name => d_cxx_csput_impl use psbn_d_cxx_mat_mod, psb_protect_name => d_cxx_csput_impl
@ -1058,7 +1236,7 @@ contains
use psb_const_mod use psb_const_mod
use psb_realloc_mod use psb_realloc_mod
use psb_string_mod use psb_string_mod
use psb_serial_mod use psb_sort_mod
implicit none implicit none
class(psbn_d_cxx_sparse_mat), intent(inout) :: a class(psbn_d_cxx_sparse_mat), intent(inout) :: a
@ -1596,3 +1774,4 @@ subroutine d_cp_cxx_from_fmt_impl(a,b,info)
if (info == 0) call a%mv_from_coo(tmp,info) if (info == 0) call a%mv_from_coo(tmp,info)
end select end select
end subroutine d_cp_cxx_from_fmt_impl end subroutine d_cp_cxx_from_fmt_impl

@ -10,10 +10,13 @@ module psbn_d_cxx_mat_mod
contains contains
procedure, pass(a) :: get_nzeros => d_cxx_get_nzeros procedure, pass(a) :: get_nzeros => d_cxx_get_nzeros
procedure, pass(a) :: get_fmt => d_cxx_get_fmt procedure, pass(a) :: get_fmt => d_cxx_get_fmt
procedure, pass(a) :: get_diag => d_cxx_get_diag
procedure, pass(a) :: d_base_csmm => d_cxx_csmm procedure, pass(a) :: d_base_csmm => d_cxx_csmm
procedure, pass(a) :: d_base_csmv => d_cxx_csmv procedure, pass(a) :: d_base_csmv => d_cxx_csmv
procedure, pass(a) :: d_base_cssm => d_cxx_cssm procedure, pass(a) :: d_base_cssm => d_cxx_cssm
procedure, pass(a) :: d_base_cssv => d_cxx_cssv procedure, pass(a) :: d_base_cssv => d_cxx_cssv
procedure, pass(a) :: d_scals => d_cxx_scals
procedure, pass(a) :: d_scal => d_cxx_scal
procedure, pass(a) :: csnmi => d_cxx_csnmi procedure, pass(a) :: csnmi => d_cxx_csnmi
procedure, pass(a) :: reallocate_nz => d_cxx_reallocate_nz procedure, pass(a) :: reallocate_nz => d_cxx_reallocate_nz
procedure, pass(a) :: csput => d_cxx_csput procedure, pass(a) :: csput => d_cxx_csput
@ -26,16 +29,20 @@ module psbn_d_cxx_mat_mod
procedure, pass(a) :: mv_from_coo => d_mv_cxx_from_coo procedure, pass(a) :: mv_from_coo => d_mv_cxx_from_coo
procedure, pass(a) :: mv_to_fmt => d_mv_cxx_to_fmt procedure, pass(a) :: mv_to_fmt => d_mv_cxx_to_fmt
procedure, pass(a) :: mv_from_fmt => d_mv_cxx_from_fmt procedure, pass(a) :: mv_from_fmt => d_mv_cxx_from_fmt
procedure, pass(a) :: d_csgetrow => d_cxx_csgetrow
procedure, pass(a) :: get_size => d_cxx_get_size
procedure, pass(a) :: free => d_cxx_free procedure, pass(a) :: free => d_cxx_free
procedure, pass(a) :: trim => d_cxx_trim
procedure, pass(a) :: print => d_cxx_print procedure, pass(a) :: print => d_cxx_print
end type psbn_d_cxx_sparse_mat end type psbn_d_cxx_sparse_mat
private :: d_cxx_get_nzeros, d_cxx_csmm, d_cxx_csmv, d_cxx_cssm, d_cxx_cssv, & private :: d_cxx_get_nzeros, d_cxx_csmm, d_cxx_csmv, d_cxx_cssm, d_cxx_cssv, &
& d_cxx_csput, d_cxx_reallocate_nz, d_cxx_allocate_mnnz, & & d_cxx_csput, d_cxx_reallocate_nz, d_cxx_allocate_mnnz, &
& d_cxx_free, d_cxx_print, d_cxx_get_fmt, d_cxx_csnmi, & & d_cxx_free, d_cxx_print, d_cxx_get_fmt, d_cxx_csnmi, get_diag, &
& d_cp_cxx_to_coo, d_cp_cxx_from_coo, & & d_cp_cxx_to_coo, d_cp_cxx_from_coo, &
& d_mv_cxx_to_coo, d_mv_cxx_from_coo, & & d_mv_cxx_to_coo, d_mv_cxx_from_coo, &
& d_cp_cxx_to_fmt, d_cp_cxx_from_fmt, & & d_cp_cxx_to_fmt, d_cp_cxx_from_fmt, &
& d_mv_cxx_to_fmt, d_mv_cxx_from_fmt & d_mv_cxx_to_fmt, d_mv_cxx_from_fmt, &
& d_cxx_scals, d_cxx_scal, d_cxx_trim, d_cxx_csgetrow, d_cxx_get_size
interface interface
@ -140,6 +147,26 @@ module psbn_d_cxx_mat_mod
end subroutine d_cxx_csput_impl end subroutine d_cxx_csput_impl
end interface end interface
interface
subroutine d_cxx_csgetrow_impl(imin,imax,a,nz,ia,ja,val,info,&
& jmin,jmax,iren,append,nzin,rscale,cscale)
use psb_const_mod
import psbn_d_cxx_sparse_mat
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
integer, intent(in) :: imin,imax
integer, intent(out) :: nz
integer, allocatable, intent(inout) :: ia(:), ja(:)
real(psb_dpk_), allocatable, intent(inout) :: val(:)
integer,intent(out) :: info
logical, intent(in), optional :: append
integer, intent(in), optional :: iren(:)
integer, intent(in), optional :: jmin,jmax, nzin
logical, intent(in), optional :: rscale,cscale
end subroutine d_cxx_csgetrow_impl
end interface
interface d_cxx_cssm_impl interface d_cxx_cssm_impl
subroutine d_cxx_cssv_impl(alpha,a,x,beta,y,info,trans) subroutine d_cxx_cssv_impl(alpha,a,x,beta,y,info,trans)
use psb_const_mod use psb_const_mod
@ -195,7 +222,7 @@ module psbn_d_cxx_mat_mod
contains contains
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !=====================================
! !
! !
! !
@ -205,7 +232,7 @@ contains
! !
! !
! !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !=====================================
function d_cxx_get_fmt(a) result(res) function d_cxx_get_fmt(a) result(res)
implicit none implicit none
@ -221,8 +248,31 @@ contains
res = a%irp(a%m+1)-1 res = a%irp(a%m+1)-1
end function d_cxx_get_nzeros end function d_cxx_get_nzeros
function d_cxx_get_size(a) result(res)
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
integer :: res
res = -1
if (allocated(a%ja)) then
if (res >= 0) then
res = min(res,size(a%ja))
else
res = size(a%ja)
end if
end if
if (allocated(a%val)) then
if (res >= 0) then
res = min(res,size(a%val))
else
res = size(a%val)
end if
end if
end function d_cxx_get_size
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !=====================================
! !
! !
! !
@ -232,7 +282,7 @@ contains
! !
! !
! !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !=====================================
subroutine d_cxx_reallocate_nz(nz,a) subroutine d_cxx_reallocate_nz(nz,a)
@ -269,7 +319,7 @@ contains
end subroutine d_cxx_reallocate_nz end subroutine d_cxx_reallocate_nz
subroutine d_cxx_csput(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) subroutine d_cxx_csput(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl)
use psb_const_mod use psb_const_mod
use psb_error_mod use psb_error_mod
implicit none implicit none
@ -316,7 +366,7 @@ contains
if (nz == 0) return if (nz == 0) return
call d_cxx_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) call d_cxx_csput_impl(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl)
if (info /= 0) goto 9999 if (info /= 0) goto 9999
call psb_erractionrestore(err_act) call psb_erractionrestore(err_act)
@ -332,6 +382,198 @@ contains
return return
end subroutine d_cxx_csput end subroutine d_cxx_csput
subroutine d_cxx_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
& jmin,jmax,iren,append,nzin,rscale,cscale)
! Output is always in COO format
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
integer, intent(in) :: imin,imax
integer, intent(out) :: nz
integer, allocatable, intent(inout) :: ia(:), ja(:)
real(psb_dpk_), allocatable, intent(inout) :: val(:)
integer,intent(out) :: info
logical, intent(in), optional :: append
integer, intent(in), optional :: iren(:)
integer, intent(in), optional :: jmin,jmax, nzin
logical, intent(in), optional :: rscale,cscale
Integer :: err_act
character(len=20) :: name='csget'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
info = 0
call d_cxx_csgetrow_impl(imin,imax,a,nz,ia,ja,val,info,&
& jmin,jmax,iren,append,nzin,rscale,cscale)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_csgetrow
subroutine d_cxx_csgetblk(imin,imax,a,b,info,&
& jmin,jmax,iren,append,rscale,cscale)
! Output is always in COO format
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
class(psbn_d_coo_sparse_mat), intent(inout) :: b
integer, intent(in) :: imin,imax
integer,intent(out) :: info
logical, intent(in), optional :: append
integer, intent(in), optional :: iren(:)
integer, intent(in), optional :: jmin,jmax
logical, intent(in), optional :: rscale,cscale
Integer :: err_act, nzin, nzout
character(len=20) :: name='csget'
logical :: append_
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
info = 0
if (present(append)) then
append_ = append
else
append_ = .false.
endif
if (append_) then
nzin = a%get_nzeros()
else
nzin = 0
endif
call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
& jmin=jmin, jmax=jmax, iren=iren, append=append_, &
& nzin=nzin, rscale=rscale, cscale=cscale)
if (info /= 0) goto 9999
call b%set_nzeros(nzin+nzout)
call b%fix(info)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_csgetblk
subroutine d_cxx_csclip(a,b,info,&
& imin,imax,jmin,jmax,rscale,cscale)
! Output is always in COO format
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
class(psbn_d_coo_sparse_mat), intent(out) :: b
integer,intent(out) :: info
integer, intent(in), optional :: imin,imax,jmin,jmax
logical, intent(in), optional :: rscale,cscale
Integer :: err_act, nzin, nzout, imin_, imax_, jmin_, jmax_, mb,nb
character(len=20) :: name='csget'
logical :: rscale_, cscale_
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
info = 0
nzin = 0
if (present(imin)) then
imin_ = imin
else
imin_ = 1
end if
if (present(imax)) then
imax_ = imax
else
imax_ = a%get_nrows()
end if
if (present(jmin)) then
jmin_ = jmin
else
jmin_ = 1
end if
if (present(jmax)) then
jmax_ = jmax
else
jmax_ = a%get_ncols()
end if
if (present(rscale)) then
rscale_ = rscale
else
rscale_ = .true.
end if
if (present(cscale)) then
cscale_ = cscale
else
cscale_ = .true.
end if
if (rscale_) then
mb = imax_ - imin_ +1
else
mb = a%get_nrows() ! Should this be imax_ ??
endif
if (cscale_) then
nb = jmax_ - jmin_ +1
else
nb = a%get_ncols() ! Should this be jmax_ ??
endif
call b%allocate(mb,nb)
call a%csget(imin_,imax_,nzout,b%ia,b%ja,b%val,info,&
& jmin=jmin_, jmax=jmax_, append=.false., &
& nzin=nzin, rscale=rscale_, cscale=cscale_)
if (info /= 0) goto 9999
call b%set_nzeros(nzin+nzout)
call b%fix(info)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_csclip
subroutine d_cxx_free(a) subroutine d_cxx_free(a)
implicit none implicit none
@ -349,6 +591,38 @@ contains
end subroutine d_cxx_free end subroutine d_cxx_free
subroutine d_cxx_trim(a)
use psb_realloc_mod
use psb_error_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(inout) :: a
Integer :: err_act, info, nz, m
character(len=20) :: name='trim'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
info = 0
m = a%get_nrows()
nz = a%get_nzeros()
if (info == 0) call psb_realloc(m+1,a%irp,info)
if (info == 0) call psb_realloc(nz,a%ja,info)
if (info == 0) call psb_realloc(nz,a%val,info)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_trim
subroutine d_cp_cxx_to_coo(a,b,info) subroutine d_cp_cxx_to_coo(a,b,info)
use psb_error_mod use psb_error_mod
@ -671,7 +945,6 @@ contains
subroutine d_cxx_print(iout,a,iv,eirs,eics,head,ivr,ivc) subroutine d_cxx_print(iout,a,iv,eirs,eics,head,ivr,ivc)
use psb_spmat_type
use psb_string_mod use psb_string_mod
implicit none implicit none
@ -752,7 +1025,7 @@ contains
end subroutine d_cxx_print end subroutine d_cxx_print
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !=====================================
! !
! !
! !
@ -763,7 +1036,7 @@ contains
! !
! !
! !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !=====================================
subroutine d_cxx_csmv(alpha,a,x,beta,y,info,trans) subroutine d_cxx_csmv(alpha,a,x,beta,y,info,trans)
@ -972,6 +1245,129 @@ contains
end function d_cxx_csnmi end function d_cxx_csnmi
subroutine d_cxx_get_diag(a,d,info)
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(out) :: d(:)
integer, intent(out) :: info
Integer :: err_act, mnm, i, j, k
character(len=20) :: name='get_diag'
logical, parameter :: debug=.false.
info = 0
call psb_erractionsave(err_act)
mnm = min(a%get_nrows(),a%get_ncols())
if (size(d) < mnm) then
info=35
call psb_errpush(info,name,i_err=(/2,size(d),0,0,0/))
goto 9999
end if
do i=1, mnm
do k=a%irp(i),a%irp(i+1)-1
j=a%ja(k)
if ((j==i) .and.(j <= mnm )) then
d(i) = a%val(k)
endif
enddo
end do
do i=mnm+1,size(d)
d(i) = dzero
end do
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_get_diag
subroutine d_cxx_scal(d,a,info)
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(inout) :: a
real(psb_dpk_), intent(in) :: d(:)
integer, intent(out) :: info
Integer :: err_act,mnm, i, j, m
character(len=20) :: name='scal'
logical, parameter :: debug=.false.
info = 0
call psb_erractionsave(err_act)
m = a%get_nrows()
if (size(d) < m) then
info=35
call psb_errpush(info,name,i_err=(/2,size(d),0,0,0/))
goto 9999
end if
do i=1, m
do j = a%irp(i), a%irp(i+1) -1
a%val(j) = a%val(j) * d(i)
end do
enddo
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_scal
subroutine d_cxx_scals(d,a,info)
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(inout) :: a
real(psb_dpk_), intent(in) :: d
integer, intent(out) :: info
Integer :: err_act,mnm, i, j, m
character(len=20) :: name='scal'
logical, parameter :: debug=.false.
info = 0
call psb_erractionsave(err_act)
do i=1,a%get_nzeros()
a%val(i) = a%val(i) * d
enddo
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_scals
end module psbn_d_cxx_mat_mod end module psbn_d_cxx_mat_mod

@ -547,13 +547,14 @@ contains
! on exit : unchanged. ! on exit : unchanged.
! !
use psb_base_mod use psb_base_mod
use psbn_d_mat_mod
implicit none implicit none
! parameters ! parameters
type(psb_dspmat_type) :: a_glob type(psb_dspmat_type) :: a_glob
real(psb_dpk_) :: b_glob(:) real(psb_dpk_) :: b_glob(:)
integer :: ictxt integer :: ictxt
type(psb_dspmat_type) :: a type(psbn_d_sparse_mat) :: a
real(psb_dpk_), allocatable :: b(:) real(psb_dpk_), allocatable :: b(:)
type(psb_desc_type) :: desc_a type(psb_desc_type) :: desc_a
integer, intent(out) :: info integer, intent(out) :: info

Loading…
Cancel
Save