From 607a0aa94946257fbde85ac561da8b3a6eb4e4f8 Mon Sep 17 00:00:00 2001 From: Salvatore Filippone Date: Thu, 10 Sep 2009 16:43:24 +0000 Subject: [PATCH] psblas3: base/Makefile base/newserial/psbn_base_mat_mod.f03 base/newserial/psbn_d_base_mat_mod.f03 base/newserial/psbn_d_coo_impl.f03 base/newserial/psbn_d_csr_impl.f03 base/newserial/psbn_d_csr_mat_mod.f03 base/newserial/psbn_mat_impl.f03 base/newserial/psbn_mat_mod.f03 test/serial/Makefile test/serial/d_matgen.f03 test/serial/psbn_d_cxx_impl.f03 test/serial/psbn_d_cxx_mat_mod.f03 Improved error messages. Various cosmetic fixes. --- base/Makefile | 1 + base/newserial/psbn_base_mat_mod.f03 | 50 +- base/newserial/psbn_d_base_mat_mod.f03 | 1129 +++++++++++++----------- base/newserial/psbn_d_coo_impl.f03 | 632 +++++-------- base/newserial/psbn_d_csr_impl.f03 | 52 ++ base/newserial/psbn_d_csr_mat_mod.f03 | 501 ++++++----- base/newserial/psbn_mat_impl.f03 | 11 +- base/newserial/psbn_mat_mod.f03 | 80 +- test/serial/Makefile | 3 +- test/serial/d_matgen.f03 | 10 +- test/serial/psbn_d_cxx_impl.f03 | 90 +- test/serial/psbn_d_cxx_mat_mod.f03 | 533 ++++++----- 12 files changed, 1659 insertions(+), 1433 deletions(-) diff --git a/base/Makefile b/base/Makefile index ac3edb53..ccaa913c 100644 --- a/base/Makefile +++ b/base/Makefile @@ -22,6 +22,7 @@ clean: (cd tools; make clean) (cd serial; make clean) (cd psblas; make clean) + (cd newserial; make clean) veryclean: clean /bin/rm -f $(HERE)/$(LIBNAME) $(LIBMOD) *$(.mod) diff --git a/base/newserial/psbn_base_mat_mod.f03 b/base/newserial/psbn_base_mat_mod.f03 index 396e9f10..04c5a964 100644 --- a/base/newserial/psbn_base_mat_mod.f03 +++ b/base/newserial/psbn_base_mat_mod.f03 @@ -30,21 +30,13 @@ module psbn_base_mat_mod integer, private :: state, duplicate logical, private :: triangle, unitd, upper, sorted contains - procedure, pass(a) :: set_nrows - procedure, pass(a) :: set_ncols - procedure, pass(a) :: set_dupl - procedure, pass(a) :: set_state - procedure, pass(a) :: set_null - procedure, pass(a) :: set_bld - procedure, pass(a) :: set_upd - procedure, pass(a) :: set_asb - procedure, pass(a) :: set_sorted - procedure, pass(a) :: set_upper - procedure, pass(a) :: set_lower - procedure, pass(a) :: set_triangle - procedure, pass(a) :: set_unit - + ! ==================================== + ! + ! Getters + ! + ! + ! ==================================== procedure, pass(a) :: get_nrows procedure, pass(a) :: get_ncols procedure, pass(a) :: get_nzeros @@ -52,8 +44,6 @@ module psbn_base_mat_mod procedure, pass(a) :: get_state procedure, pass(a) :: get_dupl procedure, pass(a) :: get_fmt - - procedure, pass(a) :: is_null procedure, pass(a) :: is_bld procedure, pass(a) :: is_upd @@ -63,13 +53,39 @@ module psbn_base_mat_mod procedure, pass(a) :: is_lower procedure, pass(a) :: is_triangle procedure, pass(a) :: is_unit + + ! ==================================== + ! + ! Setters + ! + ! ==================================== + procedure, pass(a) :: set_nrows + procedure, pass(a) :: set_ncols + procedure, pass(a) :: set_dupl + procedure, pass(a) :: set_state + procedure, pass(a) :: set_null + procedure, pass(a) :: set_bld + procedure, pass(a) :: set_upd + procedure, pass(a) :: set_asb + procedure, pass(a) :: set_sorted + procedure, pass(a) :: set_upper + procedure, pass(a) :: set_lower + procedure, pass(a) :: set_triangle + procedure, pass(a) :: set_unit + + + + ! ==================================== + ! + ! Data management + ! + ! ==================================== procedure, pass(a) :: get_neigh procedure, pass(a) :: allocate_mnnz procedure, pass(a) :: reallocate_nz procedure, pass(a) :: free generic, public :: allocate => allocate_mnnz generic, public :: reallocate => reallocate_nz - procedure, pass(a) :: print => sparse_print end type psbn_base_sparse_mat diff --git a/base/newserial/psbn_d_base_mat_mod.f03 b/base/newserial/psbn_d_base_mat_mod.f03 index bf3a832d..8f51d99f 100644 --- a/base/newserial/psbn_d_base_mat_mod.f03 +++ b/base/newserial/psbn_d_base_mat_mod.f03 @@ -1,7 +1,7 @@ module psbn_d_base_mat_mod - + use psbn_base_mat_mod - + type, extends(psbn_base_sparse_mat) :: psbn_d_base_sparse_mat contains procedure, pass(a) :: d_base_csmv @@ -10,6 +10,7 @@ module psbn_d_base_mat_mod procedure, pass(a) :: d_base_cssv procedure, pass(a) :: d_base_cssm generic, public :: cssm => d_base_cssm, d_base_cssv + procedure, pass(a) :: csnmi procedure, pass(a) :: csput procedure, pass(a) :: cp_to_coo procedure, pass(a) :: cp_from_coo @@ -21,7 +22,7 @@ module psbn_d_base_mat_mod procedure, pass(a) :: mv_from_fmt end type psbn_d_base_sparse_mat private :: d_base_csmv, d_base_csmm, d_base_cssv, d_base_cssm,& - & csput, cp_to_coo, cp_from_coo, cp_to_fmt, cp_from_fmt, & + & csnmi, csput, cp_to_coo, cp_from_coo, cp_to_fmt, cp_from_fmt, & & mv_to_coo, mv_from_coo, mv_to_fmt, mv_from_fmt @@ -33,32 +34,33 @@ module psbn_d_base_mat_mod contains - procedure, pass(a) :: get_size => d_coo_get_size - procedure, pass(a) :: get_nzeros => d_coo_get_nzeros - procedure, pass(a) :: set_nzeros => d_coo_set_nzeros - procedure, pass(a) :: d_base_csmm => d_coo_csmm - procedure, pass(a) :: d_base_csmv => d_coo_csmv - procedure, pass(a) :: d_base_cssm => d_coo_cssm - procedure, pass(a) :: d_base_cssv => d_coo_cssv - procedure, pass(a) :: csput => d_coo_csput - procedure, pass(a) :: reallocate_nz => d_coo_reallocate_nz - procedure, pass(a) :: allocate_mnnz => d_coo_allocate_mnnz - procedure, pass(a) :: cp_to_coo => d_cp_coo_to_coo - procedure, pass(a) :: cp_from_coo => d_cp_coo_from_coo - procedure, pass(a) :: cp_to_fmt => d_cp_coo_to_fmt - procedure, pass(a) :: cp_from_fmt => d_cp_coo_from_fmt - procedure, pass(a) :: mv_to_coo => d_mv_coo_to_coo - procedure, pass(a) :: mv_from_coo => d_mv_coo_from_coo - procedure, pass(a) :: mv_to_fmt => d_mv_coo_to_fmt - procedure, pass(a) :: mv_from_fmt => d_mv_coo_from_fmt - procedure, pass(a) :: fix => d_fix_coo - procedure, pass(a) :: free => d_coo_free - procedure, pass(a) :: print => d_coo_print - procedure, pass(a) :: get_fmt => d_coo_get_fmt + procedure, pass(a) :: get_size => d_coo_get_size + procedure, pass(a) :: get_nzeros => d_coo_get_nzeros + procedure, pass(a) :: set_nzeros => d_coo_set_nzeros + procedure, pass(a) :: d_base_csmm => d_coo_csmm + procedure, pass(a) :: d_base_csmv => d_coo_csmv + procedure, pass(a) :: d_base_cssm => d_coo_cssm + procedure, pass(a) :: d_base_cssv => d_coo_cssv + procedure, pass(a) :: csnmi => d_coo_csnmi + procedure, pass(a) :: csput => d_coo_csput + procedure, pass(a) :: reallocate_nz => d_coo_reallocate_nz + procedure, pass(a) :: allocate_mnnz => d_coo_allocate_mnnz + procedure, pass(a) :: cp_to_coo => d_cp_coo_to_coo + procedure, pass(a) :: cp_from_coo => d_cp_coo_from_coo + procedure, pass(a) :: cp_to_fmt => d_cp_coo_to_fmt + procedure, pass(a) :: cp_from_fmt => d_cp_coo_from_fmt + procedure, pass(a) :: mv_to_coo => d_mv_coo_to_coo + procedure, pass(a) :: mv_from_coo => d_mv_coo_from_coo + procedure, pass(a) :: mv_to_fmt => d_mv_coo_to_fmt + procedure, pass(a) :: mv_from_fmt => d_mv_coo_from_fmt + procedure, pass(a) :: fix => d_fix_coo + procedure, pass(a) :: free => d_coo_free + procedure, pass(a) :: print => d_coo_print + procedure, pass(a) :: get_fmt => d_coo_get_fmt end type psbn_d_coo_sparse_mat private :: d_coo_get_nzeros, d_coo_set_nzeros, & - & d_coo_csmm, d_coo_csmv, d_coo_cssm, d_coo_cssv, & + & d_coo_csmm, d_coo_csmv, d_coo_cssm, d_coo_cssv, d_coo_csnmi, & & d_coo_csput, d_coo_reallocate_nz, d_coo_allocate_mnnz, & & d_fix_coo, d_coo_free, d_coo_print, d_coo_get_fmt, & & d_cp_coo_to_coo, d_cp_coo_from_coo, & @@ -115,7 +117,7 @@ module psbn_d_base_mat_mod integer, intent(out) :: info end subroutine d_cp_coo_to_fmt_impl end interface - + interface subroutine d_cp_coo_from_fmt_impl(a,b,info) use psb_const_mod @@ -135,7 +137,7 @@ module psbn_d_base_mat_mod integer, intent(out) :: info end subroutine d_mv_coo_to_coo_impl end interface - + interface subroutine d_mv_coo_from_coo_impl(a,b,info) use psb_const_mod @@ -155,7 +157,7 @@ module psbn_d_base_mat_mod integer, intent(out) :: info end subroutine d_mv_coo_to_fmt_impl end interface - + interface subroutine d_mv_coo_from_fmt_impl(a,b,info) use psb_const_mod @@ -166,7 +168,7 @@ module psbn_d_base_mat_mod end subroutine d_mv_coo_from_fmt_impl end interface - + interface subroutine d_coo_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) use psb_const_mod @@ -178,7 +180,7 @@ module psbn_d_base_mat_mod integer, intent(in), optional :: gtl(:) end subroutine d_coo_csput_impl end interface - + interface d_coo_cssm_impl subroutine d_coo_cssv_impl(alpha,a,x,beta,y,info,trans) use psb_const_mod @@ -199,7 +201,7 @@ module psbn_d_base_mat_mod character, optional, intent(in) :: trans end subroutine d_coo_cssm_impl end interface - + interface d_coo_csmm_impl subroutine d_coo_csmv_impl(alpha,a,x,beta,y,info,trans) use psb_const_mod @@ -220,10 +222,34 @@ module psbn_d_base_mat_mod character, optional, intent(in) :: trans end subroutine d_coo_csmm_impl end interface - + + + interface d_coo_csnmi_impl + function d_coo_csnmi_impl(a) result(res) + use psb_const_mod + import psbn_d_coo_sparse_mat + class(psbn_d_coo_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + end function d_coo_csnmi_impl + end interface + + contains - - + + + !==================================== + ! + ! + ! + ! Data management + ! + ! + ! + ! + ! + !==================================== + + subroutine cp_to_coo(a,b,info) use psb_error_mod use psb_realloc_mod @@ -231,25 +257,25 @@ contains class(psbn_d_base_sparse_mat), intent(in) :: a class(psbn_d_coo_sparse_mat), intent(out) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='to_coo' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine cp_to_coo - + subroutine cp_from_coo(a,b,info) use psb_error_mod use psb_realloc_mod @@ -257,26 +283,26 @@ contains class(psbn_d_base_sparse_mat), intent(inout) :: a class(psbn_d_coo_sparse_mat), intent(in) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='from_coo' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine cp_from_coo - - + + subroutine cp_to_fmt(a,b,info) use psb_error_mod use psb_realloc_mod @@ -284,25 +310,25 @@ contains class(psbn_d_base_sparse_mat), intent(in) :: a class(psbn_d_base_sparse_mat), intent(out) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='to_fmt' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine cp_to_fmt - + subroutine cp_from_fmt(a,b,info) use psb_error_mod use psb_realloc_mod @@ -310,26 +336,26 @@ contains class(psbn_d_base_sparse_mat), intent(inout) :: a class(psbn_d_base_sparse_mat), intent(in) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='from_fmt' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine cp_from_fmt - - + + subroutine mv_to_coo(a,b,info) use psb_error_mod use psb_realloc_mod @@ -337,25 +363,25 @@ contains class(psbn_d_base_sparse_mat), intent(inout) :: a class(psbn_d_coo_sparse_mat), intent(out) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='to_coo' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine mv_to_coo - + subroutine mv_from_coo(a,b,info) use psb_error_mod use psb_realloc_mod @@ -363,26 +389,26 @@ contains class(psbn_d_base_sparse_mat), intent(inout) :: a class(psbn_d_coo_sparse_mat), intent(inout) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='from_coo' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine mv_from_coo - - + + subroutine mv_to_fmt(a,b,info) use psb_error_mod use psb_realloc_mod @@ -390,25 +416,25 @@ contains class(psbn_d_base_sparse_mat), intent(inout) :: a class(psbn_d_base_sparse_mat), intent(out) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='to_fmt' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine mv_to_fmt - + subroutine mv_from_fmt(a,b,info) use psb_error_mod use psb_realloc_mod @@ -416,67 +442,26 @@ contains class(psbn_d_base_sparse_mat), intent(inout) :: a class(psbn_d_base_sparse_mat), intent(inout) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='from_fmt' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine mv_from_fmt - - function d_coo_get_fmt(a) result(res) - implicit none - class(psbn_d_coo_sparse_mat), intent(in) :: a - character(len=5) :: res - res = 'COO' - end function d_coo_get_fmt - - - subroutine d_fix_coo(a,info,idir) - use psb_error_mod - use psb_const_mod - implicit none - class(psbn_d_coo_sparse_mat), intent(inout) :: a - integer, intent(out) :: info - integer, intent(in), optional :: idir - Integer :: err_act - character(len=20) :: name='fix_coo' - logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - info = 0 - call d_fix_coo_impl(a,info,idir) - - if (info /= 0) goto 9999 - - call psb_erractionrestore(err_act) - return - -9999 continue - call psb_erractionrestore(err_act) - - call psb_errpush(info,name) - - if (err_act /= psb_act_ret_) then - call psb_error() - end if - return - - - end subroutine d_fix_coo - - + + subroutine csput(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) use psb_error_mod use psb_realloc_mod @@ -486,25 +471,38 @@ contains integer, intent(in) :: nz, ia(:), ja(:), imin,imax,jmin,jmax integer, intent(out) :: info integer, intent(in), optional :: gtl(:) - + Integer :: err_act character(len=20) :: name='csput' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine csput - + + !==================================== + ! + ! + ! + ! Computational routines + ! + ! + ! + ! + ! + ! + !==================================== + subroutine d_base_csmm(alpha,a,x,beta,y,info,trans) use psb_error_mod implicit none @@ -513,25 +511,25 @@ contains real(psb_dpk_), intent(inout) :: y(:,:) integer, intent(out) :: info character, optional, intent(in) :: trans - + Integer :: err_act character(len=20) :: name='d_base_csmm' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine d_base_csmm - + subroutine d_base_csmv(alpha,a,x,beta,y,info,trans) use psb_error_mod implicit none @@ -540,26 +538,26 @@ contains real(kind(1.d0)), intent(inout) :: y(:) integer, intent(out) :: info character, optional, intent(in) :: trans - + Integer :: err_act character(len=20) :: name='d_base_csmv' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - - + + end subroutine d_base_csmv - + subroutine d_base_cssm(alpha,a,x,beta,y,info,trans) use psb_error_mod implicit none @@ -568,25 +566,25 @@ contains real(psb_dpk_), intent(inout) :: y(:,:) integer, intent(out) :: info character, optional, intent(in) :: trans - + Integer :: err_act character(len=20) :: name='d_base_cssm' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - + end subroutine d_base_cssm - + subroutine d_base_cssv(alpha,a,x,beta,y,info,trans) use psb_error_mod implicit none @@ -595,59 +593,211 @@ contains real(psb_dpk_), intent(inout) :: y(:) integer, intent(out) :: info character, optional, intent(in) :: trans - + Integer :: err_act character(len=20) :: name='d_base_cssv' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) + + call psb_get_erraction(err_act) ! This is the base version. If we get here ! 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 call psb_error() end if return - - + + end subroutine d_base_cssv - - - subroutine d_cp_coo_to_coo(a,b,info) + + function csnmi(a) result(res) use psb_error_mod - use psb_realloc_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 :: err_act - character(len=20) :: name='to_coo' + class(psbn_d_base_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + + Integer :: err_act, info + character(len=20) :: name='csnmi' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - info = 0 - call d_cp_coo_to_coo_impl(a,b,info) - if (info /= 0) goto 9999 - - call psb_erractionrestore(err_act) - return - -9999 continue - call psb_erractionrestore(err_act) - - call psb_errpush(info,name) - - if (err_act /= psb_act_ret_) then + + call psb_get_erraction(err_act) + ! This is the base version. If we get here + ! it means the derived class is incomplete, + ! so we throw an error. + info = 700 + call psb_errpush(info,name,a_err=a%get_fmt()) + write(0,*) 'Got into error path',err_act,psb_act_ret_ + if (err_act /= psb_act_ret_) then call psb_error() end if + res = -done + return - + + end function csnmi + + + + + !==================================== + ! + ! + ! + ! Getters + ! + ! + ! + ! + ! + !==================================== + + + + function d_coo_get_fmt(a) result(res) + implicit none + class(psbn_d_coo_sparse_mat), intent(in) :: a + character(len=5) :: res + res = 'COO' + end function d_coo_get_fmt + + + function d_coo_get_size(a) result(res) + implicit none + class(psbn_d_coo_sparse_mat), intent(in) :: a + integer :: res + res = -1 + + if (allocated(a%ia)) res = size(a%ia) + 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_coo_get_size + + + function d_coo_get_nzeros(a) result(res) + implicit none + class(psbn_d_coo_sparse_mat), intent(in) :: a + integer :: res + res = a%nnz + end function d_coo_get_nzeros + + + !==================================== + ! + ! + ! + ! Setters + ! + ! + ! + ! + ! + ! + !==================================== + + subroutine d_coo_set_nzeros(nz,a) + implicit none + integer, intent(in) :: nz + class(psbn_d_coo_sparse_mat), intent(inout) :: a + + a%nnz = nz + + end subroutine d_coo_set_nzeros + + !==================================== + ! + ! + ! + ! Data management + ! + ! + ! + ! + ! + !==================================== + + + subroutine d_fix_coo(a,info,idir) + use psb_error_mod + use psb_const_mod + implicit none + class(psbn_d_coo_sparse_mat), intent(inout) :: a + integer, intent(out) :: info + integer, intent(in), optional :: idir + Integer :: err_act + character(len=20) :: name='fix_coo' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + info = 0 + call d_fix_coo_impl(a,info,idir) + + if (info /= 0) goto 9999 + + call psb_erractionrestore(err_act) + return + +9999 continue + call psb_erractionrestore(err_act) + + call psb_errpush(info,name) + + if (err_act /= psb_act_ret_) then + call psb_error() + end if + return + + + end subroutine d_fix_coo + + + subroutine d_cp_coo_to_coo(a,b,info) + use psb_error_mod + use psb_realloc_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 :: err_act + character(len=20) :: name='to_coo' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + info = 0 + call d_cp_coo_to_coo_impl(a,b,info) + if (info /= 0) goto 9999 + + call psb_erractionrestore(err_act) + return + +9999 continue + call psb_erractionrestore(err_act) + + call psb_errpush(info,name) + + if (err_act /= psb_act_ret_) then + call psb_error() + end if + return + end subroutine d_cp_coo_to_coo - + subroutine d_cp_coo_from_coo(a,b,info) use psb_error_mod use psb_realloc_mod @@ -655,31 +805,31 @@ contains class(psbn_d_coo_sparse_mat), intent(inout) :: a class(psbn_d_coo_sparse_mat), intent(in) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='from_coo' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) info = 0 call d_cp_coo_from_coo_impl(a,b,info) if (info /= 0) goto 9999 - + call psb_erractionrestore(err_act) return - + 9999 continue call psb_erractionrestore(err_act) - + call psb_errpush(info,name) - + if (err_act /= psb_act_ret_) then call psb_error() end if return - + end subroutine d_cp_coo_from_coo - + subroutine d_cp_coo_to_fmt(a,b,info) use psb_error_mod use psb_realloc_mod @@ -687,31 +837,31 @@ contains class(psbn_d_coo_sparse_mat), intent(in) :: a class(psbn_d_base_sparse_mat), intent(out) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='to_coo' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) info = 0 call d_cp_coo_to_fmt_impl(a,b,info) if (info /= 0) goto 9999 - + call psb_erractionrestore(err_act) return - + 9999 continue call psb_erractionrestore(err_act) - + call psb_errpush(info,name) - + if (err_act /= psb_act_ret_) then call psb_error() end if return - + end subroutine d_cp_coo_to_fmt - + subroutine d_cp_coo_from_fmt(a,b,info) use psb_error_mod use psb_realloc_mod @@ -719,33 +869,33 @@ contains class(psbn_d_coo_sparse_mat), intent(inout) :: a class(psbn_d_base_sparse_mat), intent(in) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='from_coo' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) info = 0 call d_cp_coo_from_fmt_impl(a,b,info) if (info /= 0) goto 9999 - + call psb_erractionrestore(err_act) return - + 9999 continue call psb_erractionrestore(err_act) - + call psb_errpush(info,name) - + if (err_act /= psb_act_ret_) then call psb_error() end if return - + end subroutine d_cp_coo_from_fmt - + subroutine d_mv_coo_to_coo(a,b,info) use psb_error_mod use psb_realloc_mod @@ -753,31 +903,31 @@ contains class(psbn_d_coo_sparse_mat), intent(inout) :: a class(psbn_d_coo_sparse_mat), intent(out) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='to_coo' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) info = 0 call d_mv_coo_to_coo_impl(a,b,info) if (info /= 0) goto 9999 - + call psb_erractionrestore(err_act) return - + 9999 continue call psb_erractionrestore(err_act) - + call psb_errpush(info,name) - + if (err_act /= psb_act_ret_) then call psb_error() end if return - + end subroutine d_mv_coo_to_coo - + subroutine d_mv_coo_from_coo(a,b,info) use psb_error_mod use psb_realloc_mod @@ -785,31 +935,31 @@ contains class(psbn_d_coo_sparse_mat), intent(inout) :: a class(psbn_d_coo_sparse_mat), intent(inout) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='from_coo' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) info = 0 call d_mv_coo_from_coo_impl(a,b,info) if (info /= 0) goto 9999 - + call psb_erractionrestore(err_act) return - + 9999 continue call psb_erractionrestore(err_act) - + call psb_errpush(info,name) - + if (err_act /= psb_act_ret_) then call psb_error() end if return - + end subroutine d_mv_coo_from_coo - + subroutine d_mv_coo_to_fmt(a,b,info) use psb_error_mod use psb_realloc_mod @@ -817,31 +967,31 @@ contains class(psbn_d_coo_sparse_mat), intent(inout) :: a class(psbn_d_base_sparse_mat), intent(out) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='to_coo' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) info = 0 call d_mv_coo_to_fmt_impl(a,b,info) if (info /= 0) goto 9999 - + call psb_erractionrestore(err_act) return - + 9999 continue call psb_erractionrestore(err_act) - + call psb_errpush(info,name) - + if (err_act /= psb_act_ret_) then call psb_error() end if return - + end subroutine d_mv_coo_to_fmt - + subroutine d_mv_coo_from_fmt(a,b,info) use psb_error_mod use psb_realloc_mod @@ -849,33 +999,33 @@ contains class(psbn_d_coo_sparse_mat), intent(inout) :: a class(psbn_d_base_sparse_mat), intent(inout) :: b integer, intent(out) :: info - + Integer :: err_act character(len=20) :: name='from_coo' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) info = 0 call d_mv_coo_from_fmt_impl(a,b,info) if (info /= 0) goto 9999 - + call psb_erractionrestore(err_act) return - + 9999 continue call psb_erractionrestore(err_act) - + call psb_errpush(info,name) - + if (err_act /= psb_act_ret_) then call psb_error() end if return - + end subroutine d_mv_coo_from_fmt - + subroutine d_coo_reallocate_nz(nz,a) use psb_error_mod use psb_realloc_mod @@ -885,73 +1035,31 @@ contains Integer :: err_act, info character(len=20) :: name='d_coo_reallocate_nz' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) - + call psb_realloc(nz,a%ia,a%ja,a%val,info) - + if (info /= 0) then call psb_errpush(4000,name) goto 9999 end if - + 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_reallocate_nz - - - function d_coo_get_size(a) result(res) - implicit none - class(psbn_d_coo_sparse_mat), intent(in) :: a - integer :: res - res = -1 - if (allocated(a%ia)) res = size(a%ia) - 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_coo_get_size - - - function d_coo_get_nzeros(a) result(res) - implicit none - class(psbn_d_coo_sparse_mat), intent(in) :: a - integer :: res - res = a%nnz - end function d_coo_get_nzeros - - - subroutine d_coo_set_nzeros(nz,a) - implicit none - integer, intent(in) :: nz - class(psbn_d_coo_sparse_mat), intent(inout) :: a - - a%nnz = nz - - end subroutine d_coo_set_nzeros - - + subroutine d_coo_csput(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) use psb_error_mod use psb_realloc_mod @@ -961,16 +1069,16 @@ contains integer, intent(in) :: nz, ia(:), ja(:), imin,imax,jmin,jmax integer, intent(out) :: info integer, intent(in), optional :: gtl(:) - - + + Integer :: err_act character(len=20) :: name='d_coo_csput' logical, parameter :: debug=.false. integer :: nza, i,j,k, nzl, isza, int_err(5) - + call psb_erractionsave(err_act) info = 0 - + if (nz <= 0) then info = 10 int_err(1)=1 @@ -983,7 +1091,7 @@ contains call psb_errpush(info,name,i_err=int_err) goto 9999 end if - + if (size(ja) < nz) then info = 35 int_err(1)=3 @@ -996,27 +1104,190 @@ contains call psb_errpush(info,name,i_err=int_err) goto 9999 end if - + if (nz == 0) return nza = a%get_nzeros() call d_coo_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) 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_csput - - + + + subroutine d_coo_free(a) + implicit none + + class(psbn_d_coo_sparse_mat), intent(inout) :: a + + if (allocated(a%ia)) deallocate(a%ia) + if (allocated(a%ja)) deallocate(a%ja) + if (allocated(a%val)) deallocate(a%val) + call a%set_null() + call a%set_nrows(0) + call a%set_ncols(0) + + return + + end subroutine d_coo_free + + subroutine d_coo_allocate_mnnz(m,n,a,nz) + use psb_error_mod + use psb_realloc_mod + implicit none + integer, intent(in) :: m,n + class(psbn_d_coo_sparse_mat), intent(inout) :: a + integer, intent(in), optional :: nz + Integer :: err_act, info, nz_ + character(len=20) :: name='allocate_mnz' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + info = 0 + if (m < 0) then + info = 10 + call psb_errpush(info,name,i_err=(/1,0,0,0,0/)) + goto 9999 + endif + if (n < 0) then + info = 10 + call psb_errpush(info,name,i_err=(/2,0,0,0,0/)) + goto 9999 + endif + if (present(nz)) then + nz_ = nz + else + nz_ = max(7*m,7*n,1) + end if + if (nz_ < 0) then + info = 10 + call psb_errpush(info,name,i_err=(/3,0,0,0,0/)) + goto 9999 + endif + + if (info == 0) call psb_realloc(nz_,a%ia,info) + if (info == 0) call psb_realloc(nz_,a%ja,info) + if (info == 0) call psb_realloc(nz_,a%val,info) + if (info == 0) then + call a%set_nrows(m) + call a%set_ncols(n) + call a%set_nzeros(0) + call a%set_bld() + call a%set_triangle(.false.) + call a%set_dupl(psbn_dupl_def_) + end if + 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_allocate_mnnz + + + subroutine d_coo_print(iout,a,iv,eirs,eics,head,ivr,ivc) + use psb_spmat_type + use psb_string_mod + implicit none + + integer, intent(in) :: iout + class(psbn_d_coo_sparse_mat), intent(in) :: a + integer, intent(in), optional :: iv(:) + integer, intent(in), optional :: eirs,eics + character(len=*), optional :: head + integer, intent(in), optional :: ivr(:), ivc(:) + + Integer :: err_act + character(len=20) :: name='d_coo_print' + logical, parameter :: debug=.false. + + character(len=80) :: frmtv + integer :: irs,ics,i,j, nmx, ni, nr, nc, nz + + if (present(eirs)) then + irs = eirs + else + irs = 0 + endif + if (present(eics)) then + ics = eics + else + ics = 0 + endif + + if (present(head)) then + write(iout,'(a)') '%%MatrixMarket matrix coordinate real general' + write(iout,'(a,a)') '% ',head + write(iout,'(a)') '%' + write(iout,'(a,a)') '% COO' + endif + + nr = a%get_nrows() + nc = a%get_ncols() + nz = a%get_nzeros() + nmx = max(nr,nc,1) + ni = floor(log10(1.0*nmx)) + 1 + + write(frmtv,'(a,i3.3,a,i3.3,a)') '(2(i',ni,',1x),es26.18,1x,2(i',ni,',1x))' + write(iout,*) nr, nc, nz + if(present(iv)) then + do j=1,a%get_nzeros() + write(iout,frmtv) iv(a%ia(j)),iv(a%ja(j)),a%val(j) + enddo + else + if (present(ivr).and..not.present(ivc)) then + do j=1,a%get_nzeros() + write(iout,frmtv) ivr(a%ia(j)),a%ja(j),a%val(j) + enddo + else if (present(ivr).and.present(ivc)) then + do j=1,a%get_nzeros() + write(iout,frmtv) ivr(a%ia(j)),ivc(a%ja(j)),a%val(j) + enddo + else if (.not.present(ivr).and.present(ivc)) then + do j=1,a%get_nzeros() + write(iout,frmtv) a%ia(j),ivc(a%ja(j)),a%val(j) + enddo + else if (.not.present(ivr).and..not.present(ivc)) then + do j=1,a%get_nzeros() + write(iout,frmtv) a%ia(j),a%ja(j),a%val(j) + enddo + endif + endif + + end subroutine d_coo_print + + + + !==================================== + ! + ! + ! + ! Computational routines + ! + ! + ! + ! + ! + ! + !==================================== + subroutine d_coo_csmv(alpha,a,x,beta,y,info,trans) use psb_error_mod implicit none @@ -1025,7 +1296,7 @@ contains real(psb_dpk_), intent(inout) :: y(:) integer, intent(out) :: info character, optional, intent(in) :: trans - + character :: trans_ integer :: i,j,k,m,n, nnz, ir, jc real(psb_dpk_) :: acc @@ -1033,34 +1304,34 @@ contains Integer :: err_act character(len=20) :: name='d_coo_csmv' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) - + if (.not.a%is_asb()) then info = 1121 call psb_errpush(info,name) goto 9999 endif - - + + call d_coo_csmm_impl(alpha,a,x,beta,y,info,trans) - + 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_csmv - + subroutine d_coo_csmm(alpha,a,x,beta,y,info,trans) use psb_error_mod implicit none @@ -1069,7 +1340,7 @@ contains real(psb_dpk_), intent(inout) :: y(:,:) integer, intent(out) :: info character, optional, intent(in) :: trans - + character :: trans_ integer :: i,j,k,m,n, nnz, ir, jc, nc real(psb_dpk_), allocatable :: acc(:) @@ -1077,30 +1348,30 @@ contains Integer :: err_act character(len=20) :: name='d_coo_csmm' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) - - - + + + call d_coo_csmm_impl(alpha,a,x,beta,y,info,trans) - + 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_csmm - - + + subroutine d_coo_cssv(alpha,a,x,beta,y,info,trans) use psb_error_mod implicit none @@ -1109,7 +1380,7 @@ contains real(psb_dpk_), intent(inout) :: y(:) integer, intent(out) :: info character, optional, intent(in) :: trans - + character :: trans_ integer :: i,j,k,m,n, nnz, ir, jc real(psb_dpk_) :: acc @@ -1118,43 +1389,43 @@ contains Integer :: err_act character(len=20) :: name='d_coo_cssv' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) - + if (.not.a%is_asb()) then info = 1121 call psb_errpush(info,name) goto 9999 endif - - + + if (.not. (a%is_triangle())) then write(0,*) 'Called SM on a non-triangular mat!' info = 1121 call psb_errpush(info,name) goto 9999 end if - + call d_coo_cssm_impl(alpha,a,x,beta,y,info,trans) - + 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_cssv - - - + + + subroutine d_coo_cssm(alpha,a,x,beta,y,info,trans) use psb_error_mod implicit none @@ -1163,7 +1434,7 @@ contains real(psb_dpk_), intent(inout) :: y(:,:) integer, intent(out) :: info character, optional, intent(in) :: trans - + character :: trans_ integer :: i,j,k,m,n, nnz, ir, jc, nc real(psb_dpk_) :: acc @@ -1172,190 +1443,58 @@ contains Integer :: err_act character(len=20) :: name='d_coo_csmm' logical, parameter :: debug=.false. - + call psb_erractionsave(err_act) - + if (.not.a%is_asb()) then info = 1121 call psb_errpush(info,name) goto 9999 endif - - + + if (.not. (a%is_triangle())) then write(0,*) 'Called SM on a non-triangular mat!' info = 1121 call psb_errpush(info,name) goto 9999 end if - + call d_coo_cssm_impl(alpha,a,x,beta,y,info,trans) 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_cssm - - - subroutine d_coo_free(a) - implicit none - - class(psbn_d_coo_sparse_mat), intent(inout) :: a - - if (allocated(a%ia)) deallocate(a%ia) - if (allocated(a%ja)) deallocate(a%ja) - if (allocated(a%val)) deallocate(a%val) - call a%set_null() - call a%set_nrows(0) - call a%set_ncols(0) - - return - - end subroutine d_coo_free - - subroutine d_coo_allocate_mnnz(m,n,a,nz) + + function d_coo_csnmi(a) result(res) use psb_error_mod - use psb_realloc_mod + use psb_const_mod implicit none - integer, intent(in) :: m,n - class(psbn_d_coo_sparse_mat), intent(inout) :: a - integer, intent(in), optional :: nz - Integer :: err_act, info, nz_ - character(len=20) :: name='allocate_mnz' + class(psbn_d_coo_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + + Integer :: err_act + character(len=20) :: name='csnmi' logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - info = 0 - if (m < 0) then - info = 10 - call psb_errpush(info,name,i_err=(/1,0,0,0,0/)) - goto 9999 - endif - if (n < 0) then - info = 10 - call psb_errpush(info,name,i_err=(/2,0,0,0,0/)) - goto 9999 - endif - if (present(nz)) then - nz_ = nz - else - nz_ = max(7*m,7*n,1) - end if - if (nz_ < 0) then - info = 10 - call psb_errpush(info,name,i_err=(/3,0,0,0,0/)) - goto 9999 - endif - - if (info == 0) call psb_realloc(nz_,a%ia,info) - if (info == 0) call psb_realloc(nz_,a%ja,info) - if (info == 0) call psb_realloc(nz_,a%val,info) - if (info == 0) then - call a%set_nrows(m) - call a%set_ncols(n) - call a%set_nzeros(0) - call a%set_bld() - call a%set_triangle(.false.) - call a%set_dupl(psbn_dupl_def_) - end if - 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_allocate_mnnz - - subroutine d_coo_print(iout,a,iv,eirs,eics,head,ivr,ivc) - use psb_spmat_type - use psb_string_mod - implicit none - integer, intent(in) :: iout - class(psbn_d_coo_sparse_mat), intent(in) :: a - integer, intent(in), optional :: iv(:) - integer, intent(in), optional :: eirs,eics - character(len=*), optional :: head - integer, intent(in), optional :: ivr(:), ivc(:) + res = d_coo_csnmi_impl(a) - Integer :: err_act - character(len=20) :: name='d_coo_print' - logical, parameter :: debug=.false. + return - character(len=80) :: frmtv - integer :: irs,ics,i,j, nmx, ni, nr, nc, nz - - if (present(eirs)) then - irs = eirs - else - irs = 0 - endif - if (present(eics)) then - ics = eics - else - ics = 0 - endif - - if (present(head)) then - write(iout,'(a)') '%%MatrixMarket matrix coordinate real general' - write(iout,'(a,a)') '% ',head - write(iout,'(a)') '%' - write(iout,'(a,a)') '% COO' - endif - - nr = a%get_nrows() - nc = a%get_ncols() - nz = a%get_nzeros() - nmx = max(nr,nc,1) - ni = floor(log10(1.0*nmx)) + 1 - - write(frmtv,'(a,i3.3,a,i3.3,a)') '(2(i',ni,',1x),es26.18,1x,2(i',ni,',1x))' - write(iout,*) nr, nc, nz - if(present(iv)) then - do j=1,a%get_nzeros() - write(iout,frmtv) iv(a%ia(j)),iv(a%ja(j)),a%val(j) - enddo - else - if (present(ivr).and..not.present(ivc)) then - do j=1,a%get_nzeros() - write(iout,frmtv) ivr(a%ia(j)),a%ja(j),a%val(j) - enddo - else if (present(ivr).and.present(ivc)) then - do j=1,a%get_nzeros() - write(iout,frmtv) ivr(a%ia(j)),ivc(a%ja(j)),a%val(j) - enddo - else if (.not.present(ivr).and.present(ivc)) then - do j=1,a%get_nzeros() - write(iout,frmtv) a%ia(j),ivc(a%ja(j)),a%val(j) - enddo - else if (.not.present(ivr).and..not.present(ivc)) then - do j=1,a%get_nzeros() - write(iout,frmtv) a%ia(j),a%ja(j),a%val(j) - enddo - endif - endif + end function d_coo_csnmi - end subroutine d_coo_print - - end module psbn_d_base_mat_mod diff --git a/base/newserial/psbn_d_coo_impl.f03 b/base/newserial/psbn_d_coo_impl.f03 index 75d6bfe9..de79522e 100644 --- a/base/newserial/psbn_d_coo_impl.f03 +++ b/base/newserial/psbn_d_coo_impl.f03 @@ -2,7 +2,7 @@ subroutine d_coo_cssm_impl(alpha,a,x,beta,y,info,trans) use psb_error_mod use psbn_d_base_mat_mod, psb_protect_name => d_coo_cssm_impl - implicit none + implicit none class(psbn_d_coo_sparse_mat), intent(in) :: a real(psb_dpk_), intent(in) :: alpha, beta, x(:,:) real(psb_dpk_), intent(inout) :: y(:,:) @@ -852,6 +852,54 @@ subroutine d_coo_csmm_impl(alpha,a,x,beta,y,info,trans) end subroutine d_coo_csmm_impl +function d_coo_csnmi_impl(a) result(res) + use psb_error_mod + use psbn_d_base_mat_mod, psb_protect_name => d_coo_csnmi_impl + implicit none + class(psbn_d_coo_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + + integer :: i,j,k,m,n, nnz, ir, jc, nc + real(psb_dpk_) :: acc + real(psb_dpk_), allocatable :: tmp(:,:) + logical :: tra + Integer :: err_act + character(len=20) :: name='d_base_csnmi' + logical, parameter :: debug=.false. + + + res = dzero + nnz = a%get_nzeros() + i = 1 + j = i + do while (i<=nnz) + do while ((a%ia(j) == a%ia(i)).and. (j <= nnz)) + j = j+1 + enddo + acc = dzero + do k=i, j-1 + acc = acc + abs(a%val(k)) + end do + res = max(res,acc) + end do + +end function d_coo_csnmi_impl + + + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + ! + ! + ! + ! Data management + ! + ! + ! + ! + ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + subroutine d_coo_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) use psb_error_mod @@ -950,80 +998,6 @@ subroutine d_coo_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) contains -!!$ -!!$ subroutine psb_inner_upd(nz,ia,ja,val,nza,aspk,maxsz,& -!!$ & imin,imax,jmin,jmax,info,gtl) -!!$ implicit none -!!$ -!!$ integer, intent(in) :: nz, imin,imax,jmin,jmax,maxsz -!!$ integer, intent(in) :: ia(:),ja(:) -!!$ integer, intent(inout) :: nza -!!$ real(psb_dpk_), intent(in) :: val(:) -!!$ real(psb_dpk_), intent(inout) :: aspk(:) -!!$ integer, intent(out) :: info -!!$ integer, intent(in), optional :: gtl(:) -!!$ integer :: i,ir,ic, ng,nzl -!!$ character(len=20) :: name, ch_err -!!$ -!!$ -!!$ name='psb_inner_upd' -!!$ nzl = 0 -!!$ if (present(gtl)) then -!!$ ng = size(gtl) -!!$ if ((nza > nzl)) then -!!$ do i=1, nz -!!$ nza = nza + 1 -!!$ if (nza>maxsz) then -!!$ call psb_errpush(50,name,i_err=(/7,maxsz,5,0,nza /)) -!!$ info = -71 -!!$ return -!!$ endif -!!$ aspk(nza) = val(i) -!!$ end do -!!$ else -!!$ do i=1, nz -!!$ ir = ia(i) -!!$ ic = ja(i) -!!$ if ((ir >=1).and.(ir<=ng).and.(ic>=1).and.(ic<=ng)) then -!!$ ir = gtl(ir) -!!$ ic = gtl(ic) -!!$ if ((ir >=imin).and.(ir<=imax).and.(ic>=jmin).and.(ic<=jmax)) then -!!$ nza = nza + 1 -!!$ if (nza>maxsz) then -!!$ info = -72 -!!$ return -!!$ endif -!!$ aspk(nza) = val(i) -!!$ end if -!!$ end if -!!$ end do -!!$ end if -!!$ else -!!$ if ((nza >= nzl)) then -!!$ do i=1, nz -!!$ nza = nza + 1 -!!$ if (nza>maxsz) then -!!$ info = -73 -!!$ return -!!$ endif -!!$ aspk(nza) = val(i) -!!$ end do -!!$ else -!!$ do i=1, nz -!!$ ir = ia(i) -!!$ ic = ja(i) -!!$ if ((ir >=imin).and.(ir<=imax).and.(ic>=jmin).and.(ic<=jmax)) then -!!$ nza = nza + 1 -!!$ if (nza>maxsz) then -!!$ info = -74 -!!$ return -!!$ endif -!!$ aspk(nza) = val(i) -!!$ end if -!!$ end do -!!$ end if -!!$ end if -!!$ end subroutine psb_inner_upd subroutine psb_inner_ins(nz,ia,ja,val,nza,ia1,ia2,aspk,maxsz,& & imin,imax,jmin,jmax,info,gtl) @@ -1082,151 +1056,56 @@ contains end subroutine psb_inner_ins - subroutine d_coo_srch_upd(nz,ia,ja,val,a,& - & imin,imax,jmin,jmax,info,gtl) - - use psb_const_mod - use psb_realloc_mod - use psb_string_mod - use psb_serial_mod - implicit none - - class(psbn_d_coo_sparse_mat), intent(inout) :: a - integer, intent(in) :: nz, imin,imax,jmin,jmax - integer, intent(in) :: ia(:),ja(:) - real(psb_dpk_), intent(in) :: val(:) - integer, intent(out) :: info - integer, intent(in), optional :: gtl(:) - integer :: i,ir,ic, ilr, ilc, ip, & - & i1,i2,nc,nnz,dupl,ng - integer :: debug_level, debug_unit - character(len=20) :: name='d_coo_srch_upd' - - info = 0 - debug_unit = psb_get_debug_unit() - debug_level = psb_get_debug_level() - - dupl = a%get_dupl() - - if (.not.a%is_sorted()) then - info = -4 - return - end if - - ilr = -1 - ilc = -1 - nnz = a%get_nzeros() - - - if (present(gtl)) then - ng = size(gtl) - - select case(dupl) - case(psbn_dupl_ovwrt_,psbn_dupl_err_) - ! Overwrite. - ! Cannot test for error, should have been caught earlier. - do i=1, nz - ir = ia(i) - ic = ja(i) - if ((ir >=1).and.(ir<=ng).and.(ic>=1).and.(ic<=ng)) then - ir = gtl(ir) - if ((ir > 0).and.(ir <= a%m)) then - ic = gtl(ic) - if (ir /= ilr) then - i1 = psb_ibsrch(ir,nnz,a%ia) - i2 = i1 - do - if (i2+1 > nnz) exit - if (a%ia(i2+1) /= a%ia(i2)) exit - i2 = i2 + 1 - end do - do - if (i1-1 < 1) exit - if (a%ia(i1-1) /= a%ia(i1)) exit - i1 = i1 - 1 - end do - ilr = ir - else - i1 = 1 - i2 = 1 - end if - nc = i2-i1+1 - ip = psb_issrch(ic,nc,a%ja(i1:i2)) - if (ip>0) then - a%val(i1+ip-1) = val(i) - else - info = i - return - end if - else - if (debug_level >= psb_debug_serial_) & - & write(debug_unit,*) trim(name),& - & ': Discarding row that does not belong to us.' - endif - end if - end do - case(psbn_dupl_add_) - ! Add - do i=1, nz - ir = ia(i) - ic = ja(i) - if ((ir >=1).and.(ir<=ng).and.(ic>=1).and.(ic<=ng)) then - ir = gtl(ir) - ic = gtl(ic) - if ((ir > 0).and.(ir <= a%m)) then - - if (ir /= ilr) then - i1 = psb_ibsrch(ir,nnz,a%ia) - i2 = i1 - do - if (i2+1 > nnz) exit - if (a%ia(i2+1) /= a%ia(i2)) exit - i2 = i2 + 1 - end do - do - if (i1-1 < 1) exit - if (a%ia(i1-1) /= a%ia(i1)) exit - i1 = i1 - 1 - end do - ilr = ir - else - i1 = 1 - i2 = 1 - end if - nc = i2-i1+1 - ip = psb_issrch(ic,nc,a%ja(i1:i2)) - if (ip>0) then - a%val(i1+ip-1) = a%val(i1+ip-1) + val(i) - else - info = i - return - end if - else - if (debug_level >= psb_debug_serial_) & - & write(debug_unit,*) trim(name),& - & ': Discarding row that does not belong to us.' - end if - end if - end do - - case default - info = -3 - if (debug_level >= psb_debug_serial_) & - & write(debug_unit,*) trim(name),& - & ': Duplicate handling: ',dupl - end select - - else - - select case(dupl) - case(psbn_dupl_ovwrt_,psbn_dupl_err_) - ! Overwrite. - ! Cannot test for error, should have been caught earlier. - do i=1, nz - ir = ia(i) - ic = ja(i) - if ((ir > 0).and.(ir <= a%m)) then + subroutine d_coo_srch_upd(nz,ia,ja,val,a,& + & imin,imax,jmin,jmax,info,gtl) + + use psb_const_mod + use psb_realloc_mod + use psb_string_mod + use psb_serial_mod + implicit none + + class(psbn_d_coo_sparse_mat), intent(inout) :: a + integer, intent(in) :: nz, imin,imax,jmin,jmax + integer, intent(in) :: ia(:),ja(:) + real(psb_dpk_), intent(in) :: val(:) + integer, intent(out) :: info + integer, intent(in), optional :: gtl(:) + integer :: i,ir,ic, ilr, ilc, ip, & + & i1,i2,nc,nnz,dupl,ng + integer :: debug_level, debug_unit + character(len=20) :: name='d_coo_srch_upd' + + info = 0 + debug_unit = psb_get_debug_unit() + debug_level = psb_get_debug_level() + + dupl = a%get_dupl() + + if (.not.a%is_sorted()) then + info = -4 + return + end if + + ilr = -1 + ilc = -1 + nnz = a%get_nzeros() + + if (present(gtl)) then + ng = size(gtl) + + select case(dupl) + case(psbn_dupl_ovwrt_,psbn_dupl_err_) + ! Overwrite. + ! Cannot test for error, should have been caught earlier. + do i=1, nz + ir = ia(i) + ic = ja(i) + if ((ir >=1).and.(ir<=ng).and.(ic>=1).and.(ic<=ng)) then + ir = gtl(ir) + if ((ir > 0).and.(ir <= a%m)) then + ic = gtl(ic) if (ir /= ilr) then i1 = psb_ibsrch(ir,nnz,a%ia) i2 = i1 @@ -1253,14 +1132,21 @@ contains info = i return end if - end if - end do - - case(psbn_dupl_add_) - ! Add - do i=1, nz - ir = ia(i) - ic = ja(i) + else + if (debug_level >= psb_debug_serial_) & + & write(debug_unit,*) trim(name),& + & ': Discarding row that does not belong to us.' + endif + end if + end do + case(psbn_dupl_add_) + ! Add + do i=1, nz + ir = ia(i) + ic = ja(i) + if ((ir >=1).and.(ir<=ng).and.(ic>=1).and.(ic<=ng)) then + ir = gtl(ir) + ic = gtl(ic) if ((ir > 0).and.(ir <= a%m)) then if (ir /= ilr) then @@ -1289,19 +1175,107 @@ contains info = i return end if + else + if (debug_level >= psb_debug_serial_) & + & write(debug_unit,*) trim(name),& + & ': Discarding row that does not belong to us.' end if - end do + end if + end do - case default - info = -3 - if (debug_level >= psb_debug_serial_) & - & write(debug_unit,*) trim(name),& - & ': Duplicate handling: ',dupl - end select + case default + info = -3 + if (debug_level >= psb_debug_serial_) & + & write(debug_unit,*) trim(name),& + & ': Duplicate handling: ',dupl + end select - end if + else - end subroutine d_coo_srch_upd + select case(dupl) + case(psbn_dupl_ovwrt_,psbn_dupl_err_) + ! Overwrite. + ! Cannot test for error, should have been caught earlier. + do i=1, nz + ir = ia(i) + ic = ja(i) + if ((ir > 0).and.(ir <= a%m)) then + + if (ir /= ilr) then + i1 = psb_ibsrch(ir,nnz,a%ia) + i2 = i1 + do + if (i2+1 > nnz) exit + if (a%ia(i2+1) /= a%ia(i2)) exit + i2 = i2 + 1 + end do + do + if (i1-1 < 1) exit + if (a%ia(i1-1) /= a%ia(i1)) exit + i1 = i1 - 1 + end do + ilr = ir + else + i1 = 1 + i2 = 1 + end if + nc = i2-i1+1 + ip = psb_issrch(ic,nc,a%ja(i1:i2)) + if (ip>0) then + a%val(i1+ip-1) = val(i) + else + info = i + return + end if + end if + end do + + case(psbn_dupl_add_) + ! Add + do i=1, nz + ir = ia(i) + ic = ja(i) + if ((ir > 0).and.(ir <= a%m)) then + + if (ir /= ilr) then + i1 = psb_ibsrch(ir,nnz,a%ia) + i2 = i1 + do + if (i2+1 > nnz) exit + if (a%ia(i2+1) /= a%ia(i2)) exit + i2 = i2 + 1 + end do + do + if (i1-1 < 1) exit + if (a%ia(i1-1) /= a%ia(i1)) exit + i1 = i1 - 1 + end do + ilr = ir + else + i1 = 1 + i2 = 1 + end if + nc = i2-i1+1 + ip = psb_issrch(ic,nc,a%ja(i1:i2)) + if (ip>0) then + a%val(i1+ip-1) = a%val(i1+ip-1) + val(i) + else + info = i + return + end if + end if + end do + + case default + info = -3 + if (debug_level >= psb_debug_serial_) & + & write(debug_unit,*) trim(name),& + & ': Duplicate handling: ',dupl + end select + + end if + + end subroutine d_coo_srch_upd end subroutine d_coo_csput_impl @@ -1518,7 +1492,7 @@ subroutine d_fix_coo_impl(a,info,idir) if (nza < 2) return dupl_ = a%get_dupl() -!!$ + call d_fix_coo_inner(nza,dupl_,a%ia,a%ja,a%val,i,info,idir_) call a%set_sorted() @@ -1526,188 +1500,6 @@ subroutine d_fix_coo_impl(a,info,idir) call a%set_asb() -!!$ allocate(iaux(nza+2),stat=info) -!!$ if (info /= 0) return -!!$ -!!$ -!!$ select case(idir_) -!!$ -!!$ case(0) ! Row major order -!!$ -!!$ call msort_up(nza,a%ia(1),iaux(1),iret) -!!$ if (iret == 0) & -!!$ & call psb_ip_reord(nza,a%val,a%ia,a%ja,iaux) -!!$ i = 1 -!!$ j = i -!!$ do while (i <= nza) -!!$ do while ((a%ia(j) == a%ia(i))) -!!$ j = j+1 -!!$ if (j > nza) exit -!!$ enddo -!!$ nzl = j - i -!!$ call msort_up(nzl,a%ja(i),iaux(1),iret) -!!$ if (iret == 0) & -!!$ & call psb_ip_reord(nzl,a%val(i:i+nzl-1),& -!!$ & a%ia(i:i+nzl-1),a%ja(i:i+nzl-1),iaux) -!!$ i = j -!!$ enddo -!!$ -!!$ i = 1 -!!$ irw = a%ia(i) -!!$ icl = a%ja(i) -!!$ j = 1 -!!$ -!!$ select case(dupl_) -!!$ case(psbn_dupl_ovwrt_) -!!$ -!!$ do -!!$ j = j + 1 -!!$ if (j > nza) exit -!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then -!!$ a%val(i) = a%val(j) -!!$ else -!!$ i = i+1 -!!$ a%val(i) = a%val(j) -!!$ a%ia(i) = a%ia(j) -!!$ a%ja(i) = a%ja(j) -!!$ irw = a%ia(i) -!!$ icl = a%ja(i) -!!$ endif -!!$ enddo -!!$ -!!$ case(psbn_dupl_add_) -!!$ -!!$ do -!!$ j = j + 1 -!!$ if (j > nza) exit -!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then -!!$ a%val(i) = a%val(i) + a%val(j) -!!$ else -!!$ i = i+1 -!!$ a%val(i) = a%val(j) -!!$ a%ia(i) = a%ia(j) -!!$ a%ja(i) = a%ja(j) -!!$ irw = a%ia(i) -!!$ icl = a%ja(i) -!!$ endif -!!$ enddo -!!$ -!!$ case(psbn_dupl_err_) -!!$ do -!!$ j = j + 1 -!!$ if (j > nza) exit -!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then -!!$ call psb_errpush(130,name) -!!$ goto 9999 -!!$ else -!!$ i = i+1 -!!$ a%val(i) = a%val(j) -!!$ a%ia(i) = a%ia(j) -!!$ a%ja(i) = a%ja(j) -!!$ irw = a%ia(i) -!!$ icl = a%ja(i) -!!$ endif -!!$ enddo -!!$ case default -!!$ write(0,*) 'Error in fix_coo: unsafe dupl',dupl_ -!!$ -!!$ end select -!!$ -!!$ -!!$ if(debug_level >= psb_debug_serial_)& -!!$ & write(debug_unit,*) trim(name),': end second loop' -!!$ -!!$ case(1) ! Col major order -!!$ -!!$ call msort_up(nza,a%ja(1),iaux(1),iret) -!!$ if (iret == 0) & -!!$ & call psb_ip_reord(nza,a%val,a%ia,a%ja,iaux) -!!$ i = 1 -!!$ j = i -!!$ do while (i <= nza) -!!$ do while ((a%ja(j) == a%ja(i))) -!!$ j = j+1 -!!$ if (j > nza) exit -!!$ enddo -!!$ nzl = j - i -!!$ call msort_up(nzl,a%ia(i),iaux(1),iret) -!!$ if (iret == 0) & -!!$ & call psb_ip_reord(nzl,a%val(i:i+nzl-1),& -!!$ & a%ia(i:i+nzl-1),a%ja(i:i+nzl-1),iaux) -!!$ i = j -!!$ enddo -!!$ -!!$ i = 1 -!!$ irw = a%ia(i) -!!$ icl = a%ja(i) -!!$ j = 1 -!!$ -!!$ -!!$ select case(dupl_) -!!$ case(psbn_dupl_ovwrt_) -!!$ do -!!$ j = j + 1 -!!$ if (j > nza) exit -!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then -!!$ a%val(i) = a%val(j) -!!$ else -!!$ i = i+1 -!!$ a%val(i) = a%val(j) -!!$ a%ia(i) = a%ia(j) -!!$ a%ja(i) = a%ja(j) -!!$ irw = a%ia(i) -!!$ icl = a%ja(i) -!!$ endif -!!$ enddo -!!$ -!!$ case(psbn_dupl_add_) -!!$ do -!!$ j = j + 1 -!!$ if (j > nza) exit -!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then -!!$ a%val(i) = a%val(i) + a%val(j) -!!$ else -!!$ i = i+1 -!!$ a%val(i) = a%val(j) -!!$ a%ia(i) = a%ia(j) -!!$ a%ja(i) = a%ja(j) -!!$ irw = a%ia(i) -!!$ icl = a%ja(i) -!!$ endif -!!$ enddo -!!$ -!!$ case(psbn_dupl_err_) -!!$ do -!!$ j = j + 1 -!!$ if (j > nza) exit -!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then -!!$ call psb_errpush(130,name) -!!$ goto 9999 -!!$ else -!!$ i = i+1 -!!$ a%val(i) = a%val(j) -!!$ a%ia(i) = a%ia(j) -!!$ a%ja(i) = a%ja(j) -!!$ irw = a%ia(i) -!!$ icl = a%ja(i) -!!$ endif -!!$ enddo -!!$ case default -!!$ write(0,*) 'Error in fix_coo: unsafe dupl',dupl_ -!!$ end select -!!$ if (debug_level >= psb_debug_serial_)& -!!$ & write(debug_unit,*) trim(name),': end second loop' -!!$ case default -!!$ write(debug_unit,*) trim(name),': unknown direction ',idir_ -!!$ end select -!!$ -!!$ -!!$ call a%set_sorted() -!!$ call a%set_nzeros(i) -!!$ call a%set_asb() -!!$ -!!$ deallocate(iaux) -!!$ call psb_erractionrestore(err_act) return @@ -1719,8 +1511,6 @@ subroutine d_fix_coo_impl(a,info,idir) end if return - - end subroutine d_fix_coo_impl diff --git a/base/newserial/psbn_d_csr_impl.f03 b/base/newserial/psbn_d_csr_impl.f03 index 32352941..fe0181a6 100644 --- a/base/newserial/psbn_d_csr_impl.f03 +++ b/base/newserial/psbn_d_csr_impl.f03 @@ -1,4 +1,17 @@ +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + ! + ! + ! + ! Computational routines + ! + ! + ! + ! + ! + ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + subroutine d_csr_csmv_impl(alpha,a,x,beta,y,info,trans) use psb_error_mod use psbn_d_csr_mat_mod, psb_protect_name => d_csr_csmv_impl @@ -940,6 +953,45 @@ contains end subroutine d_csr_cssm_impl +function d_csr_csnmi_impl(a) result(res) + use psb_error_mod + use psbn_d_csr_mat_mod, psb_protect_name => d_csr_csnmi_impl + implicit none + class(psbn_d_csr_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + + integer :: i,j,k,m,n, nr, ir, jc, nc + real(psb_dpk_) :: acc + logical :: tra + Integer :: err_act + character(len=20) :: name='d_csnmi' + logical, parameter :: debug=.false. + + + res = dzero + + do i = 1, a%get_nrows() + acc = dzero + do j=a%irp(i),a%irp(i+1)-1 + acc = acc + abs(a%val(j)) + end do + res = max(res,acc) + end do + +end function d_csr_csnmi_impl + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + ! + ! + ! + ! Data management + ! + ! + ! + ! + ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + subroutine d_csr_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) diff --git a/base/newserial/psbn_d_csr_mat_mod.f03 b/base/newserial/psbn_d_csr_mat_mod.f03 index a155e8b8..0bd2f549 100644 --- a/base/newserial/psbn_d_csr_mat_mod.f03 +++ b/base/newserial/psbn_d_csr_mat_mod.f03 @@ -8,29 +8,30 @@ module psbn_d_csr_mat_mod real(psb_dpk_), allocatable :: val(:) contains - procedure, pass(a) :: get_nzeros => d_csr_get_nzeros - procedure, pass(a) :: d_base_csmm => d_csr_csmm - procedure, pass(a) :: d_base_csmv => d_csr_csmv - procedure, pass(a) :: d_base_cssm => d_csr_cssm - procedure, pass(a) :: d_base_cssv => d_csr_cssv - procedure, pass(a) :: reallocate_nz => d_csr_reallocate_nz - procedure, pass(a) :: csput => d_csr_csput - procedure, pass(a) :: allocate_mnnz => d_csr_allocate_mnnz - procedure, pass(a) :: cp_to_coo => d_cp_csr_to_coo - procedure, pass(a) :: cp_from_coo => d_cp_csr_from_coo - procedure, pass(a) :: cp_to_fmt => d_cp_csr_to_fmt - procedure, pass(a) :: cp_from_fmt => d_cp_csr_from_fmt - procedure, pass(a) :: mv_to_coo => d_mv_csr_to_coo - procedure, pass(a) :: mv_from_coo => d_mv_csr_from_coo - 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) :: free => d_csr_free - procedure, pass(a) :: print => d_csr_print - procedure, pass(a) :: get_fmt => d_csr_get_fmt + procedure, pass(a) :: get_nzeros => d_csr_get_nzeros + procedure, pass(a) :: get_fmt => d_csr_get_fmt + procedure, pass(a) :: d_base_csmm => d_csr_csmm + procedure, pass(a) :: d_base_csmv => d_csr_csmv + procedure, pass(a) :: d_base_cssm => d_csr_cssm + procedure, pass(a) :: d_base_cssv => d_csr_cssv + procedure, pass(a) :: csnmi => d_csr_csnmi + procedure, pass(a) :: reallocate_nz => d_csr_reallocate_nz + procedure, pass(a) :: csput => d_csr_csput + procedure, pass(a) :: allocate_mnnz => d_csr_allocate_mnnz + procedure, pass(a) :: cp_to_coo => d_cp_csr_to_coo + procedure, pass(a) :: cp_from_coo => d_cp_csr_from_coo + procedure, pass(a) :: cp_to_fmt => d_cp_csr_to_fmt + procedure, pass(a) :: cp_from_fmt => d_cp_csr_from_fmt + procedure, pass(a) :: mv_to_coo => d_mv_csr_to_coo + procedure, pass(a) :: mv_from_coo => d_mv_csr_from_coo + 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) :: free => d_csr_free + procedure, pass(a) :: print => d_csr_print end type psbn_d_csr_sparse_mat private :: d_csr_get_nzeros, d_csr_csmm, d_csr_csmv, d_csr_cssm, d_csr_cssv, & & d_csr_csput, d_csr_reallocate_nz, d_csr_allocate_mnnz, & - & d_csr_free, d_csr_print, d_csr_get_fmt, & + & d_csr_free, d_csr_print, d_csr_get_fmt, d_csr_csnmi, & & d_cp_csr_to_coo, d_cp_csr_from_coo, & & d_mv_csr_to_coo, d_mv_csr_from_coo, & & d_cp_csr_to_fmt, d_cp_csr_from_fmt, & @@ -181,10 +182,31 @@ module psbn_d_csr_mat_mod end subroutine d_csr_csmm_impl end interface + interface d_csr_csnmi_impl + function d_csr_csnmi_impl(a) result(res) + use psb_const_mod + import psbn_d_csr_sparse_mat + class(psbn_d_csr_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + end function d_csr_csnmi_impl + end interface + contains + !===================================== + ! + ! + ! + ! Getters + ! + ! + ! + ! + ! + !===================================== + function d_csr_get_fmt(a) result(res) implicit none class(psbn_d_csr_sparse_mat), intent(in) :: a @@ -192,6 +214,26 @@ contains res = 'CSR' end function d_csr_get_fmt + function d_csr_get_nzeros(a) result(res) + implicit none + class(psbn_d_csr_sparse_mat), intent(in) :: a + integer :: res + res = a%irp(a%m+1)-1 + end function d_csr_get_nzeros + + +!===================================== + ! + ! + ! + ! Data management + ! + ! + ! + ! + ! +!===================================== + subroutine d_csr_reallocate_nz(nz,a) use psb_error_mod @@ -227,14 +269,6 @@ contains end subroutine d_csr_reallocate_nz - function d_csr_get_nzeros(a) result(res) - implicit none - class(psbn_d_csr_sparse_mat), intent(in) :: a - integer :: res - res = a%irp(a%m+1)-1 - end function d_csr_get_nzeros - - subroutine d_csr_csput(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) use psb_const_mod use psb_error_mod @@ -299,195 +333,6 @@ contains end subroutine d_csr_csput - subroutine d_csr_csmv(alpha,a,x,beta,y,info,trans) - use psb_error_mod - implicit none - class(psbn_d_csr_sparse_mat), intent(in) :: a - real(psb_dpk_), intent(in) :: alpha, beta, x(:) - real(psb_dpk_), intent(inout) :: y(:) - integer, intent(out) :: info - character, optional, intent(in) :: trans - - character :: trans_ - integer :: i,j,k,m,n, nnz, ir, jc - real(psb_dpk_) :: acc - logical :: tra - Integer :: err_act - character(len=20) :: name='d_csr_csmv' - logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - - if (.not.a%is_asb()) then - info = 1121 - call psb_errpush(info,name) - goto 9999 - endif - - - call d_csr_csmm_impl(alpha,a,x,beta,y,info,trans) - - 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_csr_csmv - - subroutine d_csr_csmm(alpha,a,x,beta,y,info,trans) - use psb_error_mod - implicit none - class(psbn_d_csr_sparse_mat), intent(in) :: a - real(psb_dpk_), intent(in) :: alpha, beta, x(:,:) - real(psb_dpk_), intent(inout) :: y(:,:) - integer, intent(out) :: info - character, optional, intent(in) :: trans - - character :: trans_ - integer :: i,j,k,m,n, nnz, ir, jc, nc - real(psb_dpk_), allocatable :: acc(:) - logical :: tra - Integer :: err_act - character(len=20) :: name='d_csr_csmm' - logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - - - - call d_csr_csmm_impl(alpha,a,x,beta,y,info,trans) - - 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_csr_csmm - - - subroutine d_csr_cssv(alpha,a,x,beta,y,info,trans) - use psb_error_mod - implicit none - class(psbn_d_csr_sparse_mat), intent(in) :: a - real(psb_dpk_), intent(in) :: alpha, beta, x(:) - real(psb_dpk_), intent(inout) :: y(:) - integer, intent(out) :: info - character, optional, intent(in) :: trans - - character :: trans_ - integer :: i,j,k,m,n, nnz, ir, jc - real(psb_dpk_) :: acc - real(psb_dpk_), allocatable :: tmp(:) - logical :: tra - Integer :: err_act - character(len=20) :: name='d_csr_cssv' - logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - - if (.not.a%is_asb()) then - info = 1121 - call psb_errpush(info,name) - goto 9999 - endif - - - if (.not. (a%is_triangle())) then - write(0,*) 'Called SM on a non-triangular mat!' - info = 1121 - call psb_errpush(info,name) - goto 9999 - end if - - call d_csr_cssm_impl(alpha,a,x,beta,y,info,trans) - - 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_csr_cssv - - - - subroutine d_csr_cssm(alpha,a,x,beta,y,info,trans) - use psb_error_mod - implicit none - class(psbn_d_csr_sparse_mat), intent(in) :: a - real(psb_dpk_), intent(in) :: alpha, beta, x(:,:) - real(psb_dpk_), intent(inout) :: y(:,:) - integer, intent(out) :: info - character, optional, intent(in) :: trans - - character :: trans_ - integer :: i,j,k,m,n, nnz, ir, jc, nc - real(psb_dpk_) :: acc - real(psb_dpk_), allocatable :: tmp(:,:) - logical :: tra - Integer :: err_act - character(len=20) :: name='d_csr_csmm' - logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - - if (.not.a%is_asb()) then - info = 1121 - call psb_errpush(info,name) - goto 9999 - endif - - - if (.not. (a%is_triangle())) then - write(0,*) 'Called SM on a non-triangular mat!' - info = 1121 - call psb_errpush(info,name) - goto 9999 - end if - - call d_csr_cssm_impl(alpha,a,x,beta,y,info,trans) - 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_csr_cssm - - subroutine d_csr_free(a) implicit none @@ -907,4 +752,226 @@ contains end subroutine d_csr_print + !===================================== + ! + ! + ! + ! Computational routines + ! + ! + ! + ! + ! + ! + !===================================== + + + subroutine d_csr_csmv(alpha,a,x,beta,y,info,trans) + use psb_error_mod + implicit none + class(psbn_d_csr_sparse_mat), intent(in) :: a + real(psb_dpk_), intent(in) :: alpha, beta, x(:) + real(psb_dpk_), intent(inout) :: y(:) + integer, intent(out) :: info + character, optional, intent(in) :: trans + + character :: trans_ + integer :: i,j,k,m,n, nnz, ir, jc + real(psb_dpk_) :: acc + logical :: tra + Integer :: err_act + character(len=20) :: name='d_csr_csmv' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + + if (.not.a%is_asb()) then + info = 1121 + call psb_errpush(info,name) + goto 9999 + endif + + + call d_csr_csmm_impl(alpha,a,x,beta,y,info,trans) + + 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_csr_csmv + + subroutine d_csr_csmm(alpha,a,x,beta,y,info,trans) + use psb_error_mod + implicit none + class(psbn_d_csr_sparse_mat), intent(in) :: a + real(psb_dpk_), intent(in) :: alpha, beta, x(:,:) + real(psb_dpk_), intent(inout) :: y(:,:) + integer, intent(out) :: info + character, optional, intent(in) :: trans + + character :: trans_ + integer :: i,j,k,m,n, nnz, ir, jc, nc + real(psb_dpk_), allocatable :: acc(:) + logical :: tra + Integer :: err_act + character(len=20) :: name='d_csr_csmm' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + + + + call d_csr_csmm_impl(alpha,a,x,beta,y,info,trans) + + 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_csr_csmm + + + subroutine d_csr_cssv(alpha,a,x,beta,y,info,trans) + use psb_error_mod + implicit none + class(psbn_d_csr_sparse_mat), intent(in) :: a + real(psb_dpk_), intent(in) :: alpha, beta, x(:) + real(psb_dpk_), intent(inout) :: y(:) + integer, intent(out) :: info + character, optional, intent(in) :: trans + + character :: trans_ + integer :: i,j,k,m,n, nnz, ir, jc + real(psb_dpk_) :: acc + real(psb_dpk_), allocatable :: tmp(:) + logical :: tra + Integer :: err_act + character(len=20) :: name='d_csr_cssv' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + + if (.not.a%is_asb()) then + info = 1121 + call psb_errpush(info,name) + goto 9999 + endif + + + if (.not. (a%is_triangle())) then + write(0,*) 'Called SM on a non-triangular mat!' + info = 1121 + call psb_errpush(info,name) + goto 9999 + end if + + call d_csr_cssm_impl(alpha,a,x,beta,y,info,trans) + + 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_csr_cssv + + + + subroutine d_csr_cssm(alpha,a,x,beta,y,info,trans) + use psb_error_mod + implicit none + class(psbn_d_csr_sparse_mat), intent(in) :: a + real(psb_dpk_), intent(in) :: alpha, beta, x(:,:) + real(psb_dpk_), intent(inout) :: y(:,:) + integer, intent(out) :: info + character, optional, intent(in) :: trans + + character :: trans_ + integer :: i,j,k,m,n, nnz, ir, jc, nc + real(psb_dpk_) :: acc + real(psb_dpk_), allocatable :: tmp(:,:) + logical :: tra + Integer :: err_act + character(len=20) :: name='d_csr_csmm' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + + if (.not.a%is_asb()) then + info = 1121 + call psb_errpush(info,name) + goto 9999 + endif + + + if (.not. (a%is_triangle())) then + write(0,*) 'Called SM on a non-triangular mat!' + info = 1121 + call psb_errpush(info,name) + goto 9999 + end if + + call d_csr_cssm_impl(alpha,a,x,beta,y,info,trans) + 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_csr_cssm + + function d_csr_csnmi(a) result(res) + use psb_error_mod + use psb_const_mod + implicit none + class(psbn_d_csr_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + + Integer :: err_act + character(len=20) :: name='csnmi' + logical, parameter :: debug=.false. + + + res = d_csr_csnmi_impl(a) + + return + + end function d_csr_csnmi + + + end module psbn_d_csr_mat_mod diff --git a/base/newserial/psbn_mat_impl.f03 b/base/newserial/psbn_mat_impl.f03 index 4dc44a22..6d6555bf 100644 --- a/base/newserial/psbn_mat_impl.f03 +++ b/base/newserial/psbn_mat_impl.f03 @@ -125,7 +125,7 @@ subroutine psbn_d_spcnv(a,b,info,type,mold,upd,dupl) call psb_errpush(info,name) goto 9999 end if - + call altmp%cp_from_fmt(a%a, info) if (info /= 0) then @@ -208,15 +208,20 @@ subroutine psbn_d_spcnv_ip(a,info,type,mold,dupl) else allocate(psbn_d_csr_sparse_mat :: altmp, stat=info) end if + if (info /= 0) then info = 4000 call psb_errpush(info,name) goto 9999 end if - - call altmp%mv_from_fmt(a%a, info) + if (allocated(altmp)) then + call altmp%mv_from_fmt(a%a, info) + else + write(0,*) 'Unallocated altmp??' + info = -1 + end if if (info /= 0) then info = 4010 call psb_errpush(info,name,a_err="mv_from") diff --git a/base/newserial/psbn_mat_mod.f03 b/base/newserial/psbn_mat_mod.f03 index 43aea727..9b9a47e3 100644 --- a/base/newserial/psbn_mat_mod.f03 +++ b/base/newserial/psbn_mat_mod.f03 @@ -54,6 +54,7 @@ module psbn_d_mat_mod procedure, pass(a) :: print => sparse_print ! Computational routines + procedure, pass(a) :: csnmi procedure, pass(a) :: d_csmv procedure, pass(a) :: d_csmm generic, public :: csmm => d_csmm, d_csmv @@ -71,12 +72,12 @@ module psbn_d_mat_mod & reallocate_nz, free, d_csmv, d_csmm, d_cssv, d_cssm, sparse_print, & & set_nrows, set_ncols, set_dupl, set_state, set_null, set_bld, & & set_upd, set_asb, set_sorted, set_upper, set_lower, set_triangle, & - & set_unit + & set_unit, csnmi contains -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !===================================== ! ! ! @@ -86,7 +87,7 @@ contains ! ! ! -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !===================================== function sparse_get_fmt(a) result(res) @@ -342,7 +343,7 @@ contains -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !===================================== ! ! ! @@ -353,7 +354,7 @@ contains ! ! ! -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !===================================== subroutine set_nrows(m,a) @@ -762,7 +763,7 @@ contains end subroutine set_upper -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !===================================== ! ! ! @@ -772,7 +773,7 @@ contains ! ! ! -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !===================================== @@ -793,7 +794,7 @@ contains logical, parameter :: debug=.false. info = 0 - call psb_erractionsave(err_act) + call psb_get_erraction(err_act) if (.not.allocated(a%a)) then info = 1121 call psb_errpush(info,name) @@ -802,12 +803,9 @@ contains call a%a%print(iout,iv,eirs,eics,head,ivr,ivc) - - call psb_erractionrestore(err_act) return 9999 continue - call psb_erractionrestore(err_act) if (err_act == psb_act_abort_) then call psb_error() @@ -873,8 +871,8 @@ contains character(len=20) :: name='csall' logical, parameter :: debug=.false. + call psb_get_erraction(err_act) - call psb_erractionsave(err_act) info = 0 allocate(psbn_d_coo_sparse_mat :: a%a, stat=info) if (info /= 0) then @@ -885,11 +883,9 @@ contains call a%a%allocate(nr,nc,nz) call a%set_bld() - call psb_erractionrestore(err_act) return 9999 continue - call psb_erractionrestore(err_act) if (err_act == psb_act_abort_) then call psb_error() @@ -916,9 +912,6 @@ contains call a%a%reallocate(nz) - if (info /= 0) goto 9999 - - call psb_erractionrestore(err_act) return 9999 continue @@ -940,7 +933,6 @@ contains character(len=20) :: name='free' logical, parameter :: debug=.false. - call psb_erractionsave(err_act) if (.not.allocated(a%a)) then info = 1121 call psb_errpush(info,name) @@ -949,7 +941,6 @@ contains call a%a%free() - call psb_erractionrestore(err_act) return 9999 continue @@ -1113,6 +1104,9 @@ contains goto 9999 end if + if (debug) write(0,*) 'Converting from ',& + & a%get_fmt(),' to ',altmp%get_fmt() + call altmp%cp_from_fmt(a%a, info) if (info /= 0) then @@ -1202,6 +1196,9 @@ contains goto 9999 end if + if (debug) write(0,*) 'Converting in-place from ',& + & a%get_fmt(),' to ',altmp%get_fmt() + call altmp%mv_from_fmt(a%a, info) if (info /= 0) then @@ -1229,9 +1226,7 @@ contains - - -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !===================================== ! ! ! @@ -1242,7 +1237,7 @@ contains ! ! ! -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !===================================== subroutine d_csmm(alpha,a,x,beta,y,info,trans) @@ -1266,7 +1261,7 @@ contains endif call a%a%csmm(alpha,x,beta,y,info,trans) - + if (info /= 0) goto 9999 call psb_erractionrestore(err_act) return @@ -1302,7 +1297,7 @@ contains endif call a%a%csmm(alpha,x,beta,y,info,trans) - + if (info /= 0) goto 9999 call psb_erractionrestore(err_act) return @@ -1338,6 +1333,7 @@ contains endif call a%a%cssm(alpha,x,beta,y,info,trans) + if (info /= 0) goto 9999 call psb_erractionrestore(err_act) return @@ -1375,6 +1371,7 @@ contains call a%a%cssm(alpha,x,beta,y,info,trans) + if (info /= 0) goto 9999 call psb_erractionrestore(err_act) return @@ -1390,5 +1387,38 @@ contains end subroutine d_cssv + + function csnmi(a) result(res) + use psb_error_mod + use psb_const_mod + implicit none + class(psbn_d_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + + Integer :: err_act, info + character(len=20) :: name='csnmi' + logical, parameter :: debug=.false. + + if (.not.allocated(a%a)) then + info = 1121 + call psb_errpush(info,name) + goto 9999 + endif + + res = a%a%csnmi() + + + return + +9999 continue + + if (err_act == psb_act_abort_) then + call psb_error() + return + end if + return + + end function csnmi + end module psbn_d_mat_mod diff --git a/test/serial/Makefile b/test/serial/Makefile index a9f7e3ea..b51bffcc 100644 --- a/test/serial/Makefile +++ b/test/serial/Makefile @@ -38,7 +38,8 @@ spde: spde.o clean: - /bin/rm -f d_coo_matgen.o d_matgen.o tpg.o ppde.o spde.o $(EXEDIR)/ppde + /bin/rm -f d_coo_matgen.o d_matgen.o tpg.o ppde.o spde.o \ + psbn_d_cxx_mat_mod.o psbn_d_cxx_impl.o $(EXEDIR)/ppde verycleanlib: (cd ../..; make veryclean) lib: diff --git a/test/serial/d_matgen.f03 b/test/serial/d_matgen.f03 index a2a7110d..e8f9b6de 100644 --- a/test/serial/d_matgen.f03 +++ b/test/serial/d_matgen.f03 @@ -161,7 +161,7 @@ contains type(psbn_d_cxx_sparse_mat) :: acxx ! deltah dimension of each grid cell ! deltat discretization time - real(psb_dpk_) :: deltah + real(psb_dpk_) :: deltah, anorm real(psb_dpk_),parameter :: rhs=0.d0,one=1.d0,zero=0.d0 real(psb_dpk_) :: t0, t1, t2, t3, tasb, talc, ttot, tgen, tcpy, tmov real(psb_dpk_) :: a1, a2, a3, a4, b1, b2, b3 @@ -172,7 +172,7 @@ contains info = 0 name = 'create_matrix' - call psb_erractionsave(err_act) +!!$ call psb_erractionsave(err_act) call psb_info(ictxt, iam, np) @@ -380,6 +380,8 @@ contains end if tasb = psb_wtime()-t1 call a_n%print(20) + anorm = a_n%csnmi() + write(0,*) 'Nrm infinity ',anorm !!$ t1 = psb_wtime() call a_n%cscnv(info,mold=acxx) @@ -392,6 +394,8 @@ contains end if tmov = psb_wtime()-t1 call a_n%print(21) + anorm = a_n%csnmi() + write(0,*) 'Nrm infinity ',anorm !!$ if(iam == psb_root_) then @@ -404,7 +408,7 @@ contains !!$ write(*,'("-total time : ",es12.5)') ttot end if - call psb_erractionrestore(err_act) +!!$ call psb_erractionrestore(err_act) return 9999 continue diff --git a/test/serial/psbn_d_cxx_impl.f03 b/test/serial/psbn_d_cxx_impl.f03 index caee3d06..1bc1b7c0 100644 --- a/test/serial/psbn_d_cxx_impl.f03 +++ b/test/serial/psbn_d_cxx_impl.f03 @@ -1,4 +1,17 @@ +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + ! + ! + ! + ! Computational routines + ! + ! + ! + ! + ! + ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + subroutine d_cxx_csmv_impl(alpha,a,x,beta,y,info,trans) use psb_error_mod use psbn_d_cxx_mat_mod, psb_protect_name => d_cxx_csmv_impl @@ -940,6 +953,45 @@ contains end subroutine d_cxx_cssm_impl +function d_cxx_csnmi_impl(a) result(res) + use psb_error_mod + use psbn_d_cxx_mat_mod, psb_protect_name => d_cxx_csnmi_impl + implicit none + class(psbn_d_cxx_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + + integer :: i,j,k,m,n, nr, ir, jc, nc + real(psb_dpk_) :: acc + logical :: tra + Integer :: err_act + character(len=20) :: name='d_csnmi' + logical, parameter :: debug=.false. + + + res = dzero + + do i = 1, a%get_nrows() + acc = dzero + do j=a%irp(i),a%irp(i+1)-1 + acc = acc + abs(a%val(j)) + end do + res = max(res,acc) + end do + +end function d_cxx_csnmi_impl + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + ! + ! + ! + ! Data management + ! + ! + ! + ! + ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + subroutine d_cxx_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) @@ -1449,7 +1501,6 @@ subroutine d_mv_cxx_to_fmt_impl(a,b,info) end subroutine d_mv_cxx_to_fmt_impl - subroutine d_cp_cxx_to_fmt_impl(a,b,info) use psb_const_mod use psb_realloc_mod @@ -1471,6 +1522,7 @@ subroutine d_cp_cxx_to_fmt_impl(a,b,info) info = 0 + select type (b) class is (psbn_d_coo_sparse_mat) call a%cp_to_coo(b,info) @@ -1478,19 +1530,20 @@ subroutine d_cp_cxx_to_fmt_impl(a,b,info) call tmp%cp_from_fmt(a,info) if (info == 0) call b%mv_from_coo(tmp,info) end select + end subroutine d_cp_cxx_to_fmt_impl -subroutine d_cp_cxx_from_fmt_impl(a,b,info) +subroutine d_mv_cxx_from_fmt_impl(a,b,info) use psb_const_mod use psb_realloc_mod use psbn_d_base_mat_mod - use psbn_d_cxx_mat_mod, psb_protect_name => d_cp_cxx_from_fmt_impl + use psbn_d_cxx_mat_mod, psb_protect_name => d_mv_cxx_from_fmt_impl implicit none - class(psbn_d_cxx_sparse_mat), intent(inout) :: a - class(psbn_d_base_sparse_mat), intent(in) :: b - integer, intent(out) :: info + class(psbn_d_cxx_sparse_mat), intent(inout) :: a + class(psbn_d_base_sparse_mat), intent(inout) :: b + integer, intent(out) :: info !locals type(psbn_d_coo_sparse_mat) :: tmp @@ -1504,24 +1557,26 @@ subroutine d_cp_cxx_from_fmt_impl(a,b,info) select type (b) class is (psbn_d_coo_sparse_mat) - call a%cp_from_coo(b,info) + call a%mv_from_coo(b,info) class default - call tmp%cp_from_fmt(b,info) + call tmp%mv_from_fmt(b,info) if (info == 0) call a%mv_from_coo(tmp,info) end select -end subroutine d_cp_cxx_from_fmt_impl +end subroutine d_mv_cxx_from_fmt_impl -subroutine d_mv_cxx_from_fmt_impl(a,b,info) + + +subroutine d_cp_cxx_from_fmt_impl(a,b,info) use psb_const_mod use psb_realloc_mod use psbn_d_base_mat_mod - use psbn_d_cxx_mat_mod, psb_protect_name => d_mv_cxx_from_fmt_impl + use psbn_d_cxx_mat_mod, psb_protect_name => d_cp_cxx_from_fmt_impl implicit none - class(psbn_d_cxx_sparse_mat), intent(inout) :: a - class(psbn_d_base_sparse_mat), intent(inout) :: b - integer, intent(out) :: info + class(psbn_d_cxx_sparse_mat), intent(inout) :: a + class(psbn_d_base_sparse_mat), intent(in) :: b + integer, intent(out) :: info !locals type(psbn_d_coo_sparse_mat) :: tmp @@ -1535,10 +1590,9 @@ subroutine d_mv_cxx_from_fmt_impl(a,b,info) select type (b) class is (psbn_d_coo_sparse_mat) - call a%mv_from_coo(b,info) + call a%cp_from_coo(b,info) class default - call tmp%mv_from_fmt(b,info) + call tmp%cp_from_fmt(b,info) if (info == 0) call a%mv_from_coo(tmp,info) end select -end subroutine d_mv_cxx_from_fmt_impl - +end subroutine d_cp_cxx_from_fmt_impl diff --git a/test/serial/psbn_d_cxx_mat_mod.f03 b/test/serial/psbn_d_cxx_mat_mod.f03 index ebad3431..0d4414b6 100644 --- a/test/serial/psbn_d_cxx_mat_mod.f03 +++ b/test/serial/psbn_d_cxx_mat_mod.f03 @@ -8,29 +8,30 @@ module psbn_d_cxx_mat_mod real(psb_dpk_), allocatable :: val(:) contains - procedure, pass(a) :: get_nzeros => d_cxx_get_nzeros - procedure, pass(a) :: d_base_csmm => d_cxx_csmm - procedure, pass(a) :: d_base_csmv => d_cxx_csmv - procedure, pass(a) :: d_base_cssm => d_cxx_cssm - procedure, pass(a) :: d_base_cssv => d_cxx_cssv - procedure, pass(a) :: reallocate_nz => d_cxx_reallocate_nz - procedure, pass(a) :: csput => d_cxx_csput - procedure, pass(a) :: allocate_mnnz => d_cxx_allocate_mnnz - procedure, pass(a) :: cp_to_coo => d_cp_cxx_to_coo - procedure, pass(a) :: cp_from_coo => d_cp_cxx_from_coo - procedure, pass(a) :: cp_to_fmt => d_cp_cxx_to_fmt - procedure, pass(a) :: cp_from_fmt => d_cp_cxx_from_fmt - procedure, pass(a) :: mv_to_coo => d_mv_cxx_to_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_from_fmt => d_mv_cxx_from_fmt - procedure, pass(a) :: free => d_cxx_free - procedure, pass(a) :: print => d_cxx_print - procedure, pass(a) :: get_fmt => d_cxx_get_fmt + procedure, pass(a) :: get_nzeros => d_cxx_get_nzeros + procedure, pass(a) :: get_fmt => d_cxx_get_fmt + procedure, pass(a) :: d_base_csmm => d_cxx_csmm + procedure, pass(a) :: d_base_csmv => d_cxx_csmv + procedure, pass(a) :: d_base_cssm => d_cxx_cssm + procedure, pass(a) :: d_base_cssv => d_cxx_cssv + procedure, pass(a) :: csnmi => d_cxx_csnmi + procedure, pass(a) :: reallocate_nz => d_cxx_reallocate_nz + procedure, pass(a) :: csput => d_cxx_csput + procedure, pass(a) :: allocate_mnnz => d_cxx_allocate_mnnz + procedure, pass(a) :: cp_to_coo => d_cp_cxx_to_coo + procedure, pass(a) :: cp_from_coo => d_cp_cxx_from_coo + procedure, pass(a) :: cp_to_fmt => d_cp_cxx_to_fmt + procedure, pass(a) :: cp_from_fmt => d_cp_cxx_from_fmt + procedure, pass(a) :: mv_to_coo => d_mv_cxx_to_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_from_fmt => d_mv_cxx_from_fmt + procedure, pass(a) :: free => d_cxx_free + procedure, pass(a) :: print => d_cxx_print end type psbn_d_cxx_sparse_mat 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_free, d_cxx_print, d_cxx_get_fmt, & + & d_cxx_free, d_cxx_print, d_cxx_get_fmt, d_cxx_csnmi, & & d_cp_cxx_to_coo, d_cp_cxx_from_coo, & & d_mv_cxx_to_coo, d_mv_cxx_from_coo, & & d_cp_cxx_to_fmt, d_cp_cxx_from_fmt, & @@ -181,10 +182,31 @@ module psbn_d_cxx_mat_mod end subroutine d_cxx_csmm_impl end interface + interface d_cxx_csnmi_impl + function d_cxx_csnmi_impl(a) result(res) + use psb_const_mod + import psbn_d_cxx_sparse_mat + class(psbn_d_cxx_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + end function d_cxx_csnmi_impl + end interface + contains +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + ! + ! + ! + ! Getters + ! + ! + ! + ! + ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + function d_cxx_get_fmt(a) result(res) implicit none class(psbn_d_cxx_sparse_mat), intent(in) :: a @@ -192,6 +214,26 @@ contains res = 'CXX' end function d_cxx_get_fmt + function d_cxx_get_nzeros(a) result(res) + implicit none + class(psbn_d_cxx_sparse_mat), intent(in) :: a + integer :: res + res = a%irp(a%m+1)-1 + end function d_cxx_get_nzeros + + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + ! + ! + ! + ! Data management + ! + ! + ! + ! + ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + subroutine d_cxx_reallocate_nz(nz,a) use psb_error_mod @@ -227,14 +269,6 @@ contains end subroutine d_cxx_reallocate_nz - function d_cxx_get_nzeros(a) result(res) - implicit none - class(psbn_d_cxx_sparse_mat), intent(in) :: a - integer :: res - res = a%irp(a%m+1)-1 - end function d_cxx_get_nzeros - - subroutine d_cxx_csput(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl) use psb_const_mod use psb_error_mod @@ -299,195 +333,6 @@ contains end subroutine d_cxx_csput - subroutine d_cxx_csmv(alpha,a,x,beta,y,info,trans) - use psb_error_mod - implicit none - class(psbn_d_cxx_sparse_mat), intent(in) :: a - real(psb_dpk_), intent(in) :: alpha, beta, x(:) - real(psb_dpk_), intent(inout) :: y(:) - integer, intent(out) :: info - character, optional, intent(in) :: trans - - character :: trans_ - integer :: i,j,k,m,n, nnz, ir, jc - real(psb_dpk_) :: acc - logical :: tra - Integer :: err_act - character(len=20) :: name='d_cxx_csmv' - logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - - if (.not.a%is_asb()) then - info = 1121 - call psb_errpush(info,name) - goto 9999 - endif - - - call d_cxx_csmm_impl(alpha,a,x,beta,y,info,trans) - - 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_csmv - - subroutine d_cxx_csmm(alpha,a,x,beta,y,info,trans) - use psb_error_mod - implicit none - class(psbn_d_cxx_sparse_mat), intent(in) :: a - real(psb_dpk_), intent(in) :: alpha, beta, x(:,:) - real(psb_dpk_), intent(inout) :: y(:,:) - integer, intent(out) :: info - character, optional, intent(in) :: trans - - character :: trans_ - integer :: i,j,k,m,n, nnz, ir, jc, nc - real(psb_dpk_), allocatable :: acc(:) - logical :: tra - Integer :: err_act - character(len=20) :: name='d_cxx_csmm' - logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - - - - call d_cxx_csmm_impl(alpha,a,x,beta,y,info,trans) - - 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_csmm - - - subroutine d_cxx_cssv(alpha,a,x,beta,y,info,trans) - use psb_error_mod - implicit none - class(psbn_d_cxx_sparse_mat), intent(in) :: a - real(psb_dpk_), intent(in) :: alpha, beta, x(:) - real(psb_dpk_), intent(inout) :: y(:) - integer, intent(out) :: info - character, optional, intent(in) :: trans - - character :: trans_ - integer :: i,j,k,m,n, nnz, ir, jc - real(psb_dpk_) :: acc - real(psb_dpk_), allocatable :: tmp(:) - logical :: tra - Integer :: err_act - character(len=20) :: name='d_cxx_cssv' - logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - - if (.not.a%is_asb()) then - info = 1121 - call psb_errpush(info,name) - goto 9999 - endif - - - if (.not. (a%is_triangle())) then - write(0,*) 'Called SM on a non-triangular mat!' - info = 1121 - call psb_errpush(info,name) - goto 9999 - end if - - call d_cxx_cssm_impl(alpha,a,x,beta,y,info,trans) - - 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_cssv - - - - subroutine d_cxx_cssm(alpha,a,x,beta,y,info,trans) - use psb_error_mod - implicit none - class(psbn_d_cxx_sparse_mat), intent(in) :: a - real(psb_dpk_), intent(in) :: alpha, beta, x(:,:) - real(psb_dpk_), intent(inout) :: y(:,:) - integer, intent(out) :: info - character, optional, intent(in) :: trans - - character :: trans_ - integer :: i,j,k,m,n, nnz, ir, jc, nc - real(psb_dpk_) :: acc - real(psb_dpk_), allocatable :: tmp(:,:) - logical :: tra - Integer :: err_act - character(len=20) :: name='d_cxx_csmm' - logical, parameter :: debug=.false. - - call psb_erractionsave(err_act) - - if (.not.a%is_asb()) then - info = 1121 - call psb_errpush(info,name) - goto 9999 - endif - - - if (.not. (a%is_triangle())) then - write(0,*) 'Called SM on a non-triangular mat!' - info = 1121 - call psb_errpush(info,name) - goto 9999 - end if - - call d_cxx_cssm_impl(alpha,a,x,beta,y,info,trans) - 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_cssm - - subroutine d_cxx_free(a) implicit none @@ -601,22 +446,22 @@ contains return end subroutine d_cp_cxx_to_fmt - - subroutine d_mv_cxx_to_coo(a,b,info) + + subroutine d_cp_cxx_from_fmt(a,b,info) use psb_error_mod use psb_realloc_mod implicit none class(psbn_d_cxx_sparse_mat), intent(inout) :: a - class(psbn_d_coo_sparse_mat), intent(out) :: b - integer, intent(out) :: info + class(psbn_d_base_sparse_mat), intent(in) :: b + integer, intent(out) :: info Integer :: err_act - character(len=20) :: name='to_coo' + character(len=20) :: name='from_fmt' logical, parameter :: debug=.false. call psb_erractionsave(err_act) info = 0 - call d_mv_cxx_to_coo_impl(a,b,info) + call d_cp_cxx_from_fmt_impl(a,b,info) if (info /= 0) goto 9999 call psb_erractionrestore(err_act) @@ -632,24 +477,24 @@ contains end if return - end subroutine d_mv_cxx_to_coo - - - subroutine d_cp_cxx_from_fmt(a,b,info) + end subroutine d_cp_cxx_from_fmt + + + subroutine d_mv_cxx_to_coo(a,b,info) use psb_error_mod use psb_realloc_mod implicit none class(psbn_d_cxx_sparse_mat), intent(inout) :: a - class(psbn_d_base_sparse_mat), intent(in) :: b - integer, intent(out) :: info + class(psbn_d_coo_sparse_mat), intent(out) :: b + integer, intent(out) :: info Integer :: err_act - character(len=20) :: name='from_fmt' + character(len=20) :: name='to_coo' logical, parameter :: debug=.false. call psb_erractionsave(err_act) info = 0 - call d_cp_cxx_from_fmt_impl(a,b,info) + call d_mv_cxx_to_coo_impl(a,b,info) if (info /= 0) goto 9999 call psb_erractionrestore(err_act) @@ -665,8 +510,8 @@ contains end if return - end subroutine d_cp_cxx_from_fmt - + end subroutine d_mv_cxx_to_coo + subroutine d_mv_cxx_from_coo(a,b,info) use psb_error_mod use psb_realloc_mod @@ -907,4 +752,226 @@ contains end subroutine d_cxx_print +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + ! + ! + ! + ! Computational routines + ! + ! + ! + ! + ! + ! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + + subroutine d_cxx_csmv(alpha,a,x,beta,y,info,trans) + use psb_error_mod + implicit none + class(psbn_d_cxx_sparse_mat), intent(in) :: a + real(psb_dpk_), intent(in) :: alpha, beta, x(:) + real(psb_dpk_), intent(inout) :: y(:) + integer, intent(out) :: info + character, optional, intent(in) :: trans + + character :: trans_ + integer :: i,j,k,m,n, nnz, ir, jc + real(psb_dpk_) :: acc + logical :: tra + Integer :: err_act + character(len=20) :: name='d_cxx_csmv' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + + if (.not.a%is_asb()) then + info = 1121 + call psb_errpush(info,name) + goto 9999 + endif + + + call d_cxx_csmm_impl(alpha,a,x,beta,y,info,trans) + + 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_csmv + + subroutine d_cxx_csmm(alpha,a,x,beta,y,info,trans) + use psb_error_mod + implicit none + class(psbn_d_cxx_sparse_mat), intent(in) :: a + real(psb_dpk_), intent(in) :: alpha, beta, x(:,:) + real(psb_dpk_), intent(inout) :: y(:,:) + integer, intent(out) :: info + character, optional, intent(in) :: trans + + character :: trans_ + integer :: i,j,k,m,n, nnz, ir, jc, nc + real(psb_dpk_), allocatable :: acc(:) + logical :: tra + Integer :: err_act + character(len=20) :: name='d_cxx_csmm' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + + + + call d_cxx_csmm_impl(alpha,a,x,beta,y,info,trans) + + 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_csmm + + + subroutine d_cxx_cssv(alpha,a,x,beta,y,info,trans) + use psb_error_mod + implicit none + class(psbn_d_cxx_sparse_mat), intent(in) :: a + real(psb_dpk_), intent(in) :: alpha, beta, x(:) + real(psb_dpk_), intent(inout) :: y(:) + integer, intent(out) :: info + character, optional, intent(in) :: trans + + character :: trans_ + integer :: i,j,k,m,n, nnz, ir, jc + real(psb_dpk_) :: acc + real(psb_dpk_), allocatable :: tmp(:) + logical :: tra + Integer :: err_act + character(len=20) :: name='d_cxx_cssv' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + + if (.not.a%is_asb()) then + info = 1121 + call psb_errpush(info,name) + goto 9999 + endif + + + if (.not. (a%is_triangle())) then + write(0,*) 'Called SM on a non-triangular mat!' + info = 1121 + call psb_errpush(info,name) + goto 9999 + end if + + call d_cxx_cssm_impl(alpha,a,x,beta,y,info,trans) + + 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_cssv + + + + subroutine d_cxx_cssm(alpha,a,x,beta,y,info,trans) + use psb_error_mod + implicit none + class(psbn_d_cxx_sparse_mat), intent(in) :: a + real(psb_dpk_), intent(in) :: alpha, beta, x(:,:) + real(psb_dpk_), intent(inout) :: y(:,:) + integer, intent(out) :: info + character, optional, intent(in) :: trans + + character :: trans_ + integer :: i,j,k,m,n, nnz, ir, jc, nc + real(psb_dpk_) :: acc + real(psb_dpk_), allocatable :: tmp(:,:) + logical :: tra + Integer :: err_act + character(len=20) :: name='d_cxx_csmm' + logical, parameter :: debug=.false. + + call psb_erractionsave(err_act) + + if (.not.a%is_asb()) then + info = 1121 + call psb_errpush(info,name) + goto 9999 + endif + + + if (.not. (a%is_triangle())) then + write(0,*) 'Called SM on a non-triangular mat!' + info = 1121 + call psb_errpush(info,name) + goto 9999 + end if + + call d_cxx_cssm_impl(alpha,a,x,beta,y,info,trans) + 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_cssm + + function d_cxx_csnmi(a) result(res) + use psb_error_mod + use psb_const_mod + implicit none + class(psbn_d_cxx_sparse_mat), intent(in) :: a + real(psb_dpk_) :: res + + Integer :: err_act + character(len=20) :: name='csnmi' + logical, parameter :: debug=.false. + + + res = d_cxx_csnmi_impl(a) + + return + + end function d_cxx_csnmi + + + end module psbn_d_cxx_mat_mod