New method par_csr_spspmm.

merge-paraggr
Salvatore Filippone 5 years ago
parent 21aaf5cc62
commit 40755d6279

@ -32,7 +32,7 @@
Module psb_c_tools_mod
use psb_desc_mod, only : psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_
use psb_c_vect_mod, only : psb_c_base_vect_type, psb_c_vect_type
use psb_c_mat_mod, only : psb_cspmat_type, psb_lcspmat_type, psb_c_base_sparse_mat
use psb_c_mat_mod, only : psb_cspmat_type, psb_lcspmat_type, psb_c_base_sparse_mat, psb_lc_csr_sparse_mat
use psb_l_vect_mod, only : psb_l_vect_type
use psb_c_multivect_mod, only : psb_c_base_multivect_type, psb_c_multivect_type
@ -320,4 +320,18 @@ Module psb_c_tools_mod
end subroutine psb_csprn
end interface
interface psb_par_spspmm
subroutine psb_lc_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
import :: psb_lc_csr_sparse_mat, psb_desc_type, psb_ipk_
Implicit None
type(psb_lc_csr_sparse_mat),intent(in) :: acsr
type(psb_lc_csr_sparse_mat),intent(inout) :: bcsr
type(psb_lc_csr_sparse_mat),intent(out) :: ccsr
type(psb_desc_type),intent(in) :: desc_a
type(psb_desc_type),intent(inout) :: desc_c
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: data
End Subroutine psb_lc_par_csr_spspmm
end interface psb_par_spspmm
end module psb_c_tools_mod

@ -32,7 +32,7 @@
Module psb_d_tools_mod
use psb_desc_mod, only : psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_
use psb_d_vect_mod, only : psb_d_base_vect_type, psb_d_vect_type
use psb_d_mat_mod, only : psb_dspmat_type, psb_ldspmat_type, psb_d_base_sparse_mat
use psb_d_mat_mod, only : psb_dspmat_type, psb_ldspmat_type, psb_d_base_sparse_mat, psb_ld_csr_sparse_mat
use psb_l_vect_mod, only : psb_l_vect_type
use psb_d_multivect_mod, only : psb_d_base_multivect_type, psb_d_multivect_type
@ -320,4 +320,18 @@ Module psb_d_tools_mod
end subroutine psb_dsprn
end interface
interface psb_par_spspmm
subroutine psb_ld_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
import :: psb_ld_csr_sparse_mat, psb_desc_type, psb_ipk_
Implicit None
type(psb_ld_csr_sparse_mat),intent(in) :: acsr
type(psb_ld_csr_sparse_mat),intent(inout) :: bcsr
type(psb_ld_csr_sparse_mat),intent(out) :: ccsr
type(psb_desc_type),intent(in) :: desc_a
type(psb_desc_type),intent(inout) :: desc_c
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: data
End Subroutine psb_ld_par_csr_spspmm
end interface psb_par_spspmm
end module psb_d_tools_mod

@ -32,7 +32,7 @@
Module psb_s_tools_mod
use psb_desc_mod, only : psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_
use psb_s_vect_mod, only : psb_s_base_vect_type, psb_s_vect_type
use psb_s_mat_mod, only : psb_sspmat_type, psb_lsspmat_type, psb_s_base_sparse_mat
use psb_s_mat_mod, only : psb_sspmat_type, psb_lsspmat_type, psb_s_base_sparse_mat, psb_ls_csr_sparse_mat
use psb_l_vect_mod, only : psb_l_vect_type
use psb_s_multivect_mod, only : psb_s_base_multivect_type, psb_s_multivect_type
@ -320,4 +320,18 @@ Module psb_s_tools_mod
end subroutine psb_ssprn
end interface
interface psb_par_spspmm
subroutine psb_ls_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
import :: psb_ls_csr_sparse_mat, psb_desc_type, psb_ipk_
Implicit None
type(psb_ls_csr_sparse_mat),intent(in) :: acsr
type(psb_ls_csr_sparse_mat),intent(inout) :: bcsr
type(psb_ls_csr_sparse_mat),intent(out) :: ccsr
type(psb_desc_type),intent(in) :: desc_a
type(psb_desc_type),intent(inout) :: desc_c
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: data
End Subroutine psb_ls_par_csr_spspmm
end interface psb_par_spspmm
end module psb_s_tools_mod

@ -32,7 +32,7 @@
Module psb_z_tools_mod
use psb_desc_mod, only : psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_
use psb_z_vect_mod, only : psb_z_base_vect_type, psb_z_vect_type
use psb_z_mat_mod, only : psb_zspmat_type, psb_lzspmat_type, psb_z_base_sparse_mat
use psb_z_mat_mod, only : psb_zspmat_type, psb_lzspmat_type, psb_z_base_sparse_mat, psb_lz_csr_sparse_mat
use psb_l_vect_mod, only : psb_l_vect_type
use psb_z_multivect_mod, only : psb_z_base_multivect_type, psb_z_multivect_type
@ -320,4 +320,18 @@ Module psb_z_tools_mod
end subroutine psb_zsprn
end interface
interface psb_par_spspmm
subroutine psb_lz_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
import :: psb_lz_csr_sparse_mat, psb_desc_type, psb_ipk_
Implicit None
type(psb_lz_csr_sparse_mat),intent(in) :: acsr
type(psb_lz_csr_sparse_mat),intent(inout) :: bcsr
type(psb_lz_csr_sparse_mat),intent(out) :: ccsr
type(psb_desc_type),intent(in) :: desc_a
type(psb_desc_type),intent(inout) :: desc_c
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: data
End Subroutine psb_lz_par_csr_spspmm
end interface psb_par_spspmm
end module psb_z_tools_mod

@ -24,7 +24,8 @@ FOBJS = psb_cdall.o psb_cdals.o psb_cdalv.o psb_cd_inloc.o psb_cdins.o psb_cdprt
psb_zspins.o psb_zsprn.o \
psb_cspalloc.o psb_cspasb.o psb_cspfree.o\
psb_cspins.o psb_csprn.o psb_cd_set_bld.o \
psb_s_map.o psb_d_map.o psb_c_map.o psb_z_map.o
psb_s_map.o psb_d_map.o psb_c_map.o psb_z_map.o \
psb_s_par_csr_spspmm.o psb_d_par_csr_spspmm.o psb_c_par_csr_spspmm.o psb_z_par_csr_spspmm.o
# psb_lallc.o psb_lasb.o psb_lfree.o psb_lins.o \
MPFOBJS = psb_icdasb.o psb_ssphalo.o psb_dsphalo.o psb_csphalo.o psb_zsphalo.o \

@ -0,0 +1,249 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: psb_c_par_csr_spspmm.f90
!
! Subroutine: psb_c_par_csr_spspmm
! Version: complex
!
! This routine computes a parallel product of two sparse matrices
!
! C = A * B
!
! where all the matrices are stored in CSR. On input and output the matrices
! are stored with column indices in local numbering, but inermediate quantities
! are in global numbering because gathering the halo of B to multiply it
! by A implies a potential enlargement of the support.
! Also, B may have a column index space different from its row index space,
! which is obviously the same as the column space of A.
!
!
! Arguments:
! acsr - type(psb_c_csr_sparse_mat), input.
! The sparse matrix structure A
! desc_a - type(psb_desc_type), input.
! The communication descriptor of the column space of A
! bcsr - type(psb_c_csr_sparse_mat), input/output.
! The sparse matrix structure B, gets row-extended on output
! ccsr - type(psb_c_csr_sparse_mat), output
! The sparse matrix structure C
! desc_c - type(psb_desc_type), input/output.
! The communication descriptor of the column space of B
!
! info - integer, output.
! Error code.
!
!!$Subroutine psb_c_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
!!$ use psb_base_mod, psb_protect_name => psb_c_par_csr_spspmm
!!$ Implicit None
!!$
!!$ type(psb_c_csr_sparse_mat),intent(in) :: acsr
!!$ type(psb_c_csr_sparse_mat),intent(inout) :: bcsr
!!$ type(psb_c_csr_sparse_mat),intent(out) :: ccsr
!!$ type(psb_desc_type),intent(in) :: desc_a
!!$ type(psb_desc_type),intent(inout) :: desc_c
!!$ integer(psb_ipk_), intent(out) :: info
!!$ integer(psb_ipk_), intent(in), optional :: data
!!$ ! ...local scalars....
!!$ integer(psb_ipk_) :: ictxt, np,me
!!$ integer(psb_ipk_) :: ncol, nnz
!!$ type(psb_c_csr_sparse_mat) :: tcsr1
!!$ logical :: update_desc_c
!!$ integer(psb_ipk_) :: debug_level, debug_unit, err_act
!!$ character(len=20) :: name, ch_err
!!$
!!$ if(psb_get_errstatus() /= 0) return
!!$ info=psb_success_
!!$ name='psb_c_p_csr_spspmm'
!!$ call psb_erractionsave(err_act)
!!$ if (psb_errstatus_fatal()) then
!!$ info = psb_err_internal_error_ ; goto 9999
!!$ end if
!!$ debug_unit = psb_get_debug_unit()
!!$ debug_level = psb_get_debug_level()
!!$
!!$ ictxt = desc_a%get_context()
!!$
!!$ call psb_info(ictxt, me, np)
!!$
!!$ if (debug_level >= psb_debug_outer_) &
!!$ & write(debug_unit,*) me,' ',trim(name),': Start'
!!$
!!$ update_desc_c = desc_c%is_bld()
!!$
!!$ !
!!$ ! This is a bit tricky.
!!$ ! DESC_A is the descriptor of (the columns of) A, and therefore
!!$ ! of the rows of B; the columns of B, in the intended usage, span
!!$ ! a different space for which we have DESC_C.
!!$ ! We are gathering the halo rows of B to multiply by A;
!!$ ! now, the columns of B would ideally be kept in
!!$ ! global numbering, so that we can call this repeatedly to accumulate
!!$ ! the product of multiple operators, and convert to local numbering
!!$ ! at the last possible moment. However, this would imply calling
!!$ ! the serial SPSPMM with a matrix B with the GLOBAL number of columns
!!$ ! and this could be very expensive in memory. The solution is to keep B
!!$ ! in local numbering, so that only columns really appearing count, but to
!!$ ! expand the descriptor when gathering the halo, because by performing
!!$ ! the products we are extending the support of the operator; hence
!!$ ! this routine is intended to be called with a temporary descriptor
!!$ ! DESC_C which is in the BUILD state, to allow for such expansion
!!$ ! across multiple products.
!!$ ! The caller will at some later point finalize the descriptor DESC_C.
!!$ !
!!$
!!$ ncol = desc_a%get_local_cols()
!!$ call psb_sphalo(bcsr,desc_a,tcsr1,info,&
!!$ & colcnv=.true.,rowscale=.true.,outcol_glob=.true.,col_desc=desc_c,data=data)
!!$ nnz = tcsr1%get_nzeros()
!!$ if (update_desc_c) then
!!$ call desc_c%indxmap%g2lip_ins(tcsr1%ja(1:nnz),info)
!!$ else
!!$ call desc_c%indxmap%g2lip(tcsr1%ja(1:nnz),info)
!!$ end if
!!$ if (info == psb_success_) call psb_rwextd(ncol,bcsr,info,b=tcsr1)
!!$ if (info == psb_success_) call tcsr1%free()
!!$ if(info /= psb_success_) then
!!$ call psb_errpush(psb_err_internal_error_,name,a_err='Extend am3')
!!$ goto 9999
!!$ end if
!!$ call bcsr%set_ncols(desc_c%get_local_cols())
!!$
!!$
!!$ if (debug_level >= psb_debug_outer_) &
!!$ & write(debug_unit,*) me,' ',trim(name),&
!!$ & 'starting spspmm 3'
!!$ if (debug_level >= psb_debug_outer_) write(debug_unit,*) me,' ',trim(name),&
!!$ & 'starting spspmm ',acsr%get_nrows(),acsr%get_ncols(),bcsr%get_nrows(),bcsr%get_ncols()
!!$ call psb_spspmm(acsr,bcsr,ccsr,info)
!!$
!!$ call psb_erractionrestore(err_act)
!!$ return
!!$
!!$9999 call psb_error_handler(ictxt,err_act)
!!$
!!$ return
!!$
!!$End Subroutine psb_c_par_csr_spspmm
Subroutine psb_lc_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
use psb_base_mod, psb_protect_name => psb_lc_par_csr_spspmm
Implicit None
type(psb_lc_csr_sparse_mat),intent(in) :: acsr
type(psb_lc_csr_sparse_mat),intent(inout) :: bcsr
type(psb_lc_csr_sparse_mat),intent(out) :: ccsr
type(psb_desc_type),intent(in) :: desc_a
type(psb_desc_type),intent(inout) :: desc_c
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: data
! ...local scalars....
integer(psb_ipk_) :: ictxt, np,me
integer(psb_lpk_) :: nacol, nccol, nnz
type(psb_lc_csr_sparse_mat) :: tcsr1
logical :: update_desc_c
integer(psb_ipk_) :: debug_level, debug_unit, err_act
character(len=20) :: name, ch_err
if(psb_get_errstatus() /= 0) return
info=psb_success_
name='psb_lc_p_csr_spspmm'
call psb_erractionsave(err_act)
if (psb_errstatus_fatal()) then
info = psb_err_internal_error_ ; goto 9999
end if
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
ictxt = desc_a%get_context()
call psb_info(ictxt, me, np)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),': Start'
update_desc_c = desc_c%is_bld()
!
! This is a bit tricky.
! DESC_A is the descriptor of (the columns of) A, and therefore
! of the rows of B; the columns of B, in the intended usage, span
! a different space for which we have DESC_C.
! We are gathering the halo rows of B to multiply by A;
! now, the columns of B would ideally be kept in
! global numbering, so that we can call this repeatedly to accumulate
! the product of multiple operators, and convert to local numbering
! at the last possible moment. However, this would imply calling
! the serial SPSPMM with a matrix B with the GLOBAL number of columns
! and this could be very expensive in memory. The solution is to keep B
! in local numbering, so that only columns really appearing count, but to
! expand the descriptor when gathering the halo, because by performing
! the products we are extending the support of the operator; hence
! this routine is intended to be called with a temporary descriptor
! DESC_C which is in the BUILD state, to allow for such expansion
! across multiple products.
! The caller will at some later point finalize the descriptor DESC_C.
!
nacol = desc_a%get_local_cols()
call psb_sphalo(bcsr,desc_a,tcsr1,info,&
& colcnv=.true.,rowscale=.true.,outcol_glob=.true.,col_desc=desc_c,data=data)
nnz = tcsr1%get_nzeros()
if (update_desc_c) then
call desc_c%indxmap%g2lip_ins(tcsr1%ja(1:nnz),info)
else
call desc_c%indxmap%g2lip(tcsr1%ja(1:nnz),info)
end if
if (info == psb_success_) call psb_rwextd(nacol,bcsr,info,b=tcsr1)
if (info == psb_success_) call tcsr1%free()
if(info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,a_err='Extend am3')
goto 9999
end if
nccol = desc_c%get_local_cols()
call bcsr%set_ncols(nccol)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'starting spspmm 3'
if (debug_level >= psb_debug_outer_) write(debug_unit,*) me,' ',trim(name),&
& 'starting spspmm ',acsr%get_nrows(),acsr%get_ncols(),bcsr%get_nrows(),bcsr%get_ncols()
call psb_spspmm(acsr,bcsr,ccsr,info)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(ictxt,err_act)
return
End Subroutine psb_lc_par_csr_spspmm

@ -0,0 +1,249 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: psb_d_par_csr_spspmm.f90
!
! Subroutine: psb_d_par_csr_spspmm
! Version: real
!
! This routine computes a parallel product of two sparse matrices
!
! C = A * B
!
! where all the matrices are stored in CSR. On input and output the matrices
! are stored with column indices in local numbering, but inermediate quantities
! are in global numbering because gathering the halo of B to multiply it
! by A implies a potential enlargement of the support.
! Also, B may have a column index space different from its row index space,
! which is obviously the same as the column space of A.
!
!
! Arguments:
! acsr - type(psb_d_csr_sparse_mat), input.
! The sparse matrix structure A
! desc_a - type(psb_desc_type), input.
! The communication descriptor of the column space of A
! bcsr - type(psb_d_csr_sparse_mat), input/output.
! The sparse matrix structure B, gets row-extended on output
! ccsr - type(psb_d_csr_sparse_mat), output
! The sparse matrix structure C
! desc_c - type(psb_desc_type), input/output.
! The communication descriptor of the column space of B
!
! info - integer, output.
! Error code.
!
!!$Subroutine psb_d_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
!!$ use psb_base_mod, psb_protect_name => psb_d_par_csr_spspmm
!!$ Implicit None
!!$
!!$ type(psb_d_csr_sparse_mat),intent(in) :: acsr
!!$ type(psb_d_csr_sparse_mat),intent(inout) :: bcsr
!!$ type(psb_d_csr_sparse_mat),intent(out) :: ccsr
!!$ type(psb_desc_type),intent(in) :: desc_a
!!$ type(psb_desc_type),intent(inout) :: desc_c
!!$ integer(psb_ipk_), intent(out) :: info
!!$ integer(psb_ipk_), intent(in), optional :: data
!!$ ! ...local scalars....
!!$ integer(psb_ipk_) :: ictxt, np,me
!!$ integer(psb_ipk_) :: ncol, nnz
!!$ type(psb_d_csr_sparse_mat) :: tcsr1
!!$ logical :: update_desc_c
!!$ integer(psb_ipk_) :: debug_level, debug_unit, err_act
!!$ character(len=20) :: name, ch_err
!!$
!!$ if(psb_get_errstatus() /= 0) return
!!$ info=psb_success_
!!$ name='psb_d_p_csr_spspmm'
!!$ call psb_erractionsave(err_act)
!!$ if (psb_errstatus_fatal()) then
!!$ info = psb_err_internal_error_ ; goto 9999
!!$ end if
!!$ debug_unit = psb_get_debug_unit()
!!$ debug_level = psb_get_debug_level()
!!$
!!$ ictxt = desc_a%get_context()
!!$
!!$ call psb_info(ictxt, me, np)
!!$
!!$ if (debug_level >= psb_debug_outer_) &
!!$ & write(debug_unit,*) me,' ',trim(name),': Start'
!!$
!!$ update_desc_c = desc_c%is_bld()
!!$
!!$ !
!!$ ! This is a bit tricky.
!!$ ! DESC_A is the descriptor of (the columns of) A, and therefore
!!$ ! of the rows of B; the columns of B, in the intended usage, span
!!$ ! a different space for which we have DESC_C.
!!$ ! We are gathering the halo rows of B to multiply by A;
!!$ ! now, the columns of B would ideally be kept in
!!$ ! global numbering, so that we can call this repeatedly to accumulate
!!$ ! the product of multiple operators, and convert to local numbering
!!$ ! at the last possible moment. However, this would imply calling
!!$ ! the serial SPSPMM with a matrix B with the GLOBAL number of columns
!!$ ! and this could be very expensive in memory. The solution is to keep B
!!$ ! in local numbering, so that only columns really appearing count, but to
!!$ ! expand the descriptor when gathering the halo, because by performing
!!$ ! the products we are extending the support of the operator; hence
!!$ ! this routine is intended to be called with a temporary descriptor
!!$ ! DESC_C which is in the BUILD state, to allow for such expansion
!!$ ! across multiple products.
!!$ ! The caller will at some later point finalize the descriptor DESC_C.
!!$ !
!!$
!!$ ncol = desc_a%get_local_cols()
!!$ call psb_sphalo(bcsr,desc_a,tcsr1,info,&
!!$ & colcnv=.true.,rowscale=.true.,outcol_glob=.true.,col_desc=desc_c,data=data)
!!$ nnz = tcsr1%get_nzeros()
!!$ if (update_desc_c) then
!!$ call desc_c%indxmap%g2lip_ins(tcsr1%ja(1:nnz),info)
!!$ else
!!$ call desc_c%indxmap%g2lip(tcsr1%ja(1:nnz),info)
!!$ end if
!!$ if (info == psb_success_) call psb_rwextd(ncol,bcsr,info,b=tcsr1)
!!$ if (info == psb_success_) call tcsr1%free()
!!$ if(info /= psb_success_) then
!!$ call psb_errpush(psb_err_internal_error_,name,a_err='Extend am3')
!!$ goto 9999
!!$ end if
!!$ call bcsr%set_ncols(desc_c%get_local_cols())
!!$
!!$
!!$ if (debug_level >= psb_debug_outer_) &
!!$ & write(debug_unit,*) me,' ',trim(name),&
!!$ & 'starting spspmm 3'
!!$ if (debug_level >= psb_debug_outer_) write(debug_unit,*) me,' ',trim(name),&
!!$ & 'starting spspmm ',acsr%get_nrows(),acsr%get_ncols(),bcsr%get_nrows(),bcsr%get_ncols()
!!$ call psb_spspmm(acsr,bcsr,ccsr,info)
!!$
!!$ call psb_erractionrestore(err_act)
!!$ return
!!$
!!$9999 call psb_error_handler(ictxt,err_act)
!!$
!!$ return
!!$
!!$End Subroutine psb_d_par_csr_spspmm
Subroutine psb_ld_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
use psb_base_mod, psb_protect_name => psb_ld_par_csr_spspmm
Implicit None
type(psb_ld_csr_sparse_mat),intent(in) :: acsr
type(psb_ld_csr_sparse_mat),intent(inout) :: bcsr
type(psb_ld_csr_sparse_mat),intent(out) :: ccsr
type(psb_desc_type),intent(in) :: desc_a
type(psb_desc_type),intent(inout) :: desc_c
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: data
! ...local scalars....
integer(psb_ipk_) :: ictxt, np,me
integer(psb_lpk_) :: nacol, nccol, nnz
type(psb_ld_csr_sparse_mat) :: tcsr1
logical :: update_desc_c
integer(psb_ipk_) :: debug_level, debug_unit, err_act
character(len=20) :: name, ch_err
if(psb_get_errstatus() /= 0) return
info=psb_success_
name='psb_ld_p_csr_spspmm'
call psb_erractionsave(err_act)
if (psb_errstatus_fatal()) then
info = psb_err_internal_error_ ; goto 9999
end if
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
ictxt = desc_a%get_context()
call psb_info(ictxt, me, np)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),': Start'
update_desc_c = desc_c%is_bld()
!
! This is a bit tricky.
! DESC_A is the descriptor of (the columns of) A, and therefore
! of the rows of B; the columns of B, in the intended usage, span
! a different space for which we have DESC_C.
! We are gathering the halo rows of B to multiply by A;
! now, the columns of B would ideally be kept in
! global numbering, so that we can call this repeatedly to accumulate
! the product of multiple operators, and convert to local numbering
! at the last possible moment. However, this would imply calling
! the serial SPSPMM with a matrix B with the GLOBAL number of columns
! and this could be very expensive in memory. The solution is to keep B
! in local numbering, so that only columns really appearing count, but to
! expand the descriptor when gathering the halo, because by performing
! the products we are extending the support of the operator; hence
! this routine is intended to be called with a temporary descriptor
! DESC_C which is in the BUILD state, to allow for such expansion
! across multiple products.
! The caller will at some later point finalize the descriptor DESC_C.
!
nacol = desc_a%get_local_cols()
call psb_sphalo(bcsr,desc_a,tcsr1,info,&
& colcnv=.true.,rowscale=.true.,outcol_glob=.true.,col_desc=desc_c,data=data)
nnz = tcsr1%get_nzeros()
if (update_desc_c) then
call desc_c%indxmap%g2lip_ins(tcsr1%ja(1:nnz),info)
else
call desc_c%indxmap%g2lip(tcsr1%ja(1:nnz),info)
end if
if (info == psb_success_) call psb_rwextd(nacol,bcsr,info,b=tcsr1)
if (info == psb_success_) call tcsr1%free()
if(info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,a_err='Extend am3')
goto 9999
end if
nccol = desc_c%get_local_cols()
call bcsr%set_ncols(nccol)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'starting spspmm 3'
if (debug_level >= psb_debug_outer_) write(debug_unit,*) me,' ',trim(name),&
& 'starting spspmm ',acsr%get_nrows(),acsr%get_ncols(),bcsr%get_nrows(),bcsr%get_ncols()
call psb_spspmm(acsr,bcsr,ccsr,info)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(ictxt,err_act)
return
End Subroutine psb_ld_par_csr_spspmm

@ -0,0 +1,249 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: psb_s_par_csr_spspmm.f90
!
! Subroutine: psb_s_par_csr_spspmm
! Version: real
!
! This routine computes a parallel product of two sparse matrices
!
! C = A * B
!
! where all the matrices are stored in CSR. On input and output the matrices
! are stored with column indices in local numbering, but inermediate quantities
! are in global numbering because gathering the halo of B to multiply it
! by A implies a potential enlargement of the support.
! Also, B may have a column index space different from its row index space,
! which is obviously the same as the column space of A.
!
!
! Arguments:
! acsr - type(psb_s_csr_sparse_mat), input.
! The sparse matrix structure A
! desc_a - type(psb_desc_type), input.
! The communication descriptor of the column space of A
! bcsr - type(psb_s_csr_sparse_mat), input/output.
! The sparse matrix structure B, gets row-extended on output
! ccsr - type(psb_s_csr_sparse_mat), output
! The sparse matrix structure C
! desc_c - type(psb_desc_type), input/output.
! The communication descriptor of the column space of B
!
! info - integer, output.
! Error code.
!
!!$Subroutine psb_s_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
!!$ use psb_base_mod, psb_protect_name => psb_s_par_csr_spspmm
!!$ Implicit None
!!$
!!$ type(psb_s_csr_sparse_mat),intent(in) :: acsr
!!$ type(psb_s_csr_sparse_mat),intent(inout) :: bcsr
!!$ type(psb_s_csr_sparse_mat),intent(out) :: ccsr
!!$ type(psb_desc_type),intent(in) :: desc_a
!!$ type(psb_desc_type),intent(inout) :: desc_c
!!$ integer(psb_ipk_), intent(out) :: info
!!$ integer(psb_ipk_), intent(in), optional :: data
!!$ ! ...local scalars....
!!$ integer(psb_ipk_) :: ictxt, np,me
!!$ integer(psb_ipk_) :: ncol, nnz
!!$ type(psb_s_csr_sparse_mat) :: tcsr1
!!$ logical :: update_desc_c
!!$ integer(psb_ipk_) :: debug_level, debug_unit, err_act
!!$ character(len=20) :: name, ch_err
!!$
!!$ if(psb_get_errstatus() /= 0) return
!!$ info=psb_success_
!!$ name='psb_s_p_csr_spspmm'
!!$ call psb_erractionsave(err_act)
!!$ if (psb_errstatus_fatal()) then
!!$ info = psb_err_internal_error_ ; goto 9999
!!$ end if
!!$ debug_unit = psb_get_debug_unit()
!!$ debug_level = psb_get_debug_level()
!!$
!!$ ictxt = desc_a%get_context()
!!$
!!$ call psb_info(ictxt, me, np)
!!$
!!$ if (debug_level >= psb_debug_outer_) &
!!$ & write(debug_unit,*) me,' ',trim(name),': Start'
!!$
!!$ update_desc_c = desc_c%is_bld()
!!$
!!$ !
!!$ ! This is a bit tricky.
!!$ ! DESC_A is the descriptor of (the columns of) A, and therefore
!!$ ! of the rows of B; the columns of B, in the intended usage, span
!!$ ! a different space for which we have DESC_C.
!!$ ! We are gathering the halo rows of B to multiply by A;
!!$ ! now, the columns of B would ideally be kept in
!!$ ! global numbering, so that we can call this repeatedly to accumulate
!!$ ! the product of multiple operators, and convert to local numbering
!!$ ! at the last possible moment. However, this would imply calling
!!$ ! the serial SPSPMM with a matrix B with the GLOBAL number of columns
!!$ ! and this could be very expensive in memory. The solution is to keep B
!!$ ! in local numbering, so that only columns really appearing count, but to
!!$ ! expand the descriptor when gathering the halo, because by performing
!!$ ! the products we are extending the support of the operator; hence
!!$ ! this routine is intended to be called with a temporary descriptor
!!$ ! DESC_C which is in the BUILD state, to allow for such expansion
!!$ ! across multiple products.
!!$ ! The caller will at some later point finalize the descriptor DESC_C.
!!$ !
!!$
!!$ ncol = desc_a%get_local_cols()
!!$ call psb_sphalo(bcsr,desc_a,tcsr1,info,&
!!$ & colcnv=.true.,rowscale=.true.,outcol_glob=.true.,col_desc=desc_c,data=data)
!!$ nnz = tcsr1%get_nzeros()
!!$ if (update_desc_c) then
!!$ call desc_c%indxmap%g2lip_ins(tcsr1%ja(1:nnz),info)
!!$ else
!!$ call desc_c%indxmap%g2lip(tcsr1%ja(1:nnz),info)
!!$ end if
!!$ if (info == psb_success_) call psb_rwextd(ncol,bcsr,info,b=tcsr1)
!!$ if (info == psb_success_) call tcsr1%free()
!!$ if(info /= psb_success_) then
!!$ call psb_errpush(psb_err_internal_error_,name,a_err='Extend am3')
!!$ goto 9999
!!$ end if
!!$ call bcsr%set_ncols(desc_c%get_local_cols())
!!$
!!$
!!$ if (debug_level >= psb_debug_outer_) &
!!$ & write(debug_unit,*) me,' ',trim(name),&
!!$ & 'starting spspmm 3'
!!$ if (debug_level >= psb_debug_outer_) write(debug_unit,*) me,' ',trim(name),&
!!$ & 'starting spspmm ',acsr%get_nrows(),acsr%get_ncols(),bcsr%get_nrows(),bcsr%get_ncols()
!!$ call psb_spspmm(acsr,bcsr,ccsr,info)
!!$
!!$ call psb_erractionrestore(err_act)
!!$ return
!!$
!!$9999 call psb_error_handler(ictxt,err_act)
!!$
!!$ return
!!$
!!$End Subroutine psb_s_par_csr_spspmm
Subroutine psb_ls_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
use psb_base_mod, psb_protect_name => psb_ls_par_csr_spspmm
Implicit None
type(psb_ls_csr_sparse_mat),intent(in) :: acsr
type(psb_ls_csr_sparse_mat),intent(inout) :: bcsr
type(psb_ls_csr_sparse_mat),intent(out) :: ccsr
type(psb_desc_type),intent(in) :: desc_a
type(psb_desc_type),intent(inout) :: desc_c
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: data
! ...local scalars....
integer(psb_ipk_) :: ictxt, np,me
integer(psb_lpk_) :: nacol, nccol, nnz
type(psb_ls_csr_sparse_mat) :: tcsr1
logical :: update_desc_c
integer(psb_ipk_) :: debug_level, debug_unit, err_act
character(len=20) :: name, ch_err
if(psb_get_errstatus() /= 0) return
info=psb_success_
name='psb_ls_p_csr_spspmm'
call psb_erractionsave(err_act)
if (psb_errstatus_fatal()) then
info = psb_err_internal_error_ ; goto 9999
end if
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
ictxt = desc_a%get_context()
call psb_info(ictxt, me, np)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),': Start'
update_desc_c = desc_c%is_bld()
!
! This is a bit tricky.
! DESC_A is the descriptor of (the columns of) A, and therefore
! of the rows of B; the columns of B, in the intended usage, span
! a different space for which we have DESC_C.
! We are gathering the halo rows of B to multiply by A;
! now, the columns of B would ideally be kept in
! global numbering, so that we can call this repeatedly to accumulate
! the product of multiple operators, and convert to local numbering
! at the last possible moment. However, this would imply calling
! the serial SPSPMM with a matrix B with the GLOBAL number of columns
! and this could be very expensive in memory. The solution is to keep B
! in local numbering, so that only columns really appearing count, but to
! expand the descriptor when gathering the halo, because by performing
! the products we are extending the support of the operator; hence
! this routine is intended to be called with a temporary descriptor
! DESC_C which is in the BUILD state, to allow for such expansion
! across multiple products.
! The caller will at some later point finalize the descriptor DESC_C.
!
nacol = desc_a%get_local_cols()
call psb_sphalo(bcsr,desc_a,tcsr1,info,&
& colcnv=.true.,rowscale=.true.,outcol_glob=.true.,col_desc=desc_c,data=data)
nnz = tcsr1%get_nzeros()
if (update_desc_c) then
call desc_c%indxmap%g2lip_ins(tcsr1%ja(1:nnz),info)
else
call desc_c%indxmap%g2lip(tcsr1%ja(1:nnz),info)
end if
if (info == psb_success_) call psb_rwextd(nacol,bcsr,info,b=tcsr1)
if (info == psb_success_) call tcsr1%free()
if(info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,a_err='Extend am3')
goto 9999
end if
nccol = desc_c%get_local_cols()
call bcsr%set_ncols(nccol)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'starting spspmm 3'
if (debug_level >= psb_debug_outer_) write(debug_unit,*) me,' ',trim(name),&
& 'starting spspmm ',acsr%get_nrows(),acsr%get_ncols(),bcsr%get_nrows(),bcsr%get_ncols()
call psb_spspmm(acsr,bcsr,ccsr,info)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(ictxt,err_act)
return
End Subroutine psb_ls_par_csr_spspmm

@ -0,0 +1,249 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: psb_z_par_csr_spspmm.f90
!
! Subroutine: psb_z_par_csr_spspmm
! Version: complex
!
! This routine computes a parallel product of two sparse matrices
!
! C = A * B
!
! where all the matrices are stored in CSR. On input and output the matrices
! are stored with column indices in local numbering, but inermediate quantities
! are in global numbering because gathering the halo of B to multiply it
! by A implies a potential enlargement of the support.
! Also, B may have a column index space different from its row index space,
! which is obviously the same as the column space of A.
!
!
! Arguments:
! acsr - type(psb_z_csr_sparse_mat), input.
! The sparse matrix structure A
! desc_a - type(psb_desc_type), input.
! The communication descriptor of the column space of A
! bcsr - type(psb_z_csr_sparse_mat), input/output.
! The sparse matrix structure B, gets row-extended on output
! ccsr - type(psb_z_csr_sparse_mat), output
! The sparse matrix structure C
! desc_c - type(psb_desc_type), input/output.
! The communication descriptor of the column space of B
!
! info - integer, output.
! Error code.
!
!!$Subroutine psb_z_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
!!$ use psb_base_mod, psb_protect_name => psb_z_par_csr_spspmm
!!$ Implicit None
!!$
!!$ type(psb_z_csr_sparse_mat),intent(in) :: acsr
!!$ type(psb_z_csr_sparse_mat),intent(inout) :: bcsr
!!$ type(psb_z_csr_sparse_mat),intent(out) :: ccsr
!!$ type(psb_desc_type),intent(in) :: desc_a
!!$ type(psb_desc_type),intent(inout) :: desc_c
!!$ integer(psb_ipk_), intent(out) :: info
!!$ integer(psb_ipk_), intent(in), optional :: data
!!$ ! ...local scalars....
!!$ integer(psb_ipk_) :: ictxt, np,me
!!$ integer(psb_ipk_) :: ncol, nnz
!!$ type(psb_z_csr_sparse_mat) :: tcsr1
!!$ logical :: update_desc_c
!!$ integer(psb_ipk_) :: debug_level, debug_unit, err_act
!!$ character(len=20) :: name, ch_err
!!$
!!$ if(psb_get_errstatus() /= 0) return
!!$ info=psb_success_
!!$ name='psb_z_p_csr_spspmm'
!!$ call psb_erractionsave(err_act)
!!$ if (psb_errstatus_fatal()) then
!!$ info = psb_err_internal_error_ ; goto 9999
!!$ end if
!!$ debug_unit = psb_get_debug_unit()
!!$ debug_level = psb_get_debug_level()
!!$
!!$ ictxt = desc_a%get_context()
!!$
!!$ call psb_info(ictxt, me, np)
!!$
!!$ if (debug_level >= psb_debug_outer_) &
!!$ & write(debug_unit,*) me,' ',trim(name),': Start'
!!$
!!$ update_desc_c = desc_c%is_bld()
!!$
!!$ !
!!$ ! This is a bit tricky.
!!$ ! DESC_A is the descriptor of (the columns of) A, and therefore
!!$ ! of the rows of B; the columns of B, in the intended usage, span
!!$ ! a different space for which we have DESC_C.
!!$ ! We are gathering the halo rows of B to multiply by A;
!!$ ! now, the columns of B would ideally be kept in
!!$ ! global numbering, so that we can call this repeatedly to accumulate
!!$ ! the product of multiple operators, and convert to local numbering
!!$ ! at the last possible moment. However, this would imply calling
!!$ ! the serial SPSPMM with a matrix B with the GLOBAL number of columns
!!$ ! and this could be very expensive in memory. The solution is to keep B
!!$ ! in local numbering, so that only columns really appearing count, but to
!!$ ! expand the descriptor when gathering the halo, because by performing
!!$ ! the products we are extending the support of the operator; hence
!!$ ! this routine is intended to be called with a temporary descriptor
!!$ ! DESC_C which is in the BUILD state, to allow for such expansion
!!$ ! across multiple products.
!!$ ! The caller will at some later point finalize the descriptor DESC_C.
!!$ !
!!$
!!$ ncol = desc_a%get_local_cols()
!!$ call psb_sphalo(bcsr,desc_a,tcsr1,info,&
!!$ & colcnv=.true.,rowscale=.true.,outcol_glob=.true.,col_desc=desc_c,data=data)
!!$ nnz = tcsr1%get_nzeros()
!!$ if (update_desc_c) then
!!$ call desc_c%indxmap%g2lip_ins(tcsr1%ja(1:nnz),info)
!!$ else
!!$ call desc_c%indxmap%g2lip(tcsr1%ja(1:nnz),info)
!!$ end if
!!$ if (info == psb_success_) call psb_rwextd(ncol,bcsr,info,b=tcsr1)
!!$ if (info == psb_success_) call tcsr1%free()
!!$ if(info /= psb_success_) then
!!$ call psb_errpush(psb_err_internal_error_,name,a_err='Extend am3')
!!$ goto 9999
!!$ end if
!!$ call bcsr%set_ncols(desc_c%get_local_cols())
!!$
!!$
!!$ if (debug_level >= psb_debug_outer_) &
!!$ & write(debug_unit,*) me,' ',trim(name),&
!!$ & 'starting spspmm 3'
!!$ if (debug_level >= psb_debug_outer_) write(debug_unit,*) me,' ',trim(name),&
!!$ & 'starting spspmm ',acsr%get_nrows(),acsr%get_ncols(),bcsr%get_nrows(),bcsr%get_ncols()
!!$ call psb_spspmm(acsr,bcsr,ccsr,info)
!!$
!!$ call psb_erractionrestore(err_act)
!!$ return
!!$
!!$9999 call psb_error_handler(ictxt,err_act)
!!$
!!$ return
!!$
!!$End Subroutine psb_z_par_csr_spspmm
Subroutine psb_lz_par_csr_spspmm(acsr,desc_a,bcsr,ccsr,desc_c,info,data)
use psb_base_mod, psb_protect_name => psb_lz_par_csr_spspmm
Implicit None
type(psb_lz_csr_sparse_mat),intent(in) :: acsr
type(psb_lz_csr_sparse_mat),intent(inout) :: bcsr
type(psb_lz_csr_sparse_mat),intent(out) :: ccsr
type(psb_desc_type),intent(in) :: desc_a
type(psb_desc_type),intent(inout) :: desc_c
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: data
! ...local scalars....
integer(psb_ipk_) :: ictxt, np,me
integer(psb_lpk_) :: nacol, nccol, nnz
type(psb_lz_csr_sparse_mat) :: tcsr1
logical :: update_desc_c
integer(psb_ipk_) :: debug_level, debug_unit, err_act
character(len=20) :: name, ch_err
if(psb_get_errstatus() /= 0) return
info=psb_success_
name='psb_lz_p_csr_spspmm'
call psb_erractionsave(err_act)
if (psb_errstatus_fatal()) then
info = psb_err_internal_error_ ; goto 9999
end if
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
ictxt = desc_a%get_context()
call psb_info(ictxt, me, np)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),': Start'
update_desc_c = desc_c%is_bld()
!
! This is a bit tricky.
! DESC_A is the descriptor of (the columns of) A, and therefore
! of the rows of B; the columns of B, in the intended usage, span
! a different space for which we have DESC_C.
! We are gathering the halo rows of B to multiply by A;
! now, the columns of B would ideally be kept in
! global numbering, so that we can call this repeatedly to accumulate
! the product of multiple operators, and convert to local numbering
! at the last possible moment. However, this would imply calling
! the serial SPSPMM with a matrix B with the GLOBAL number of columns
! and this could be very expensive in memory. The solution is to keep B
! in local numbering, so that only columns really appearing count, but to
! expand the descriptor when gathering the halo, because by performing
! the products we are extending the support of the operator; hence
! this routine is intended to be called with a temporary descriptor
! DESC_C which is in the BUILD state, to allow for such expansion
! across multiple products.
! The caller will at some later point finalize the descriptor DESC_C.
!
nacol = desc_a%get_local_cols()
call psb_sphalo(bcsr,desc_a,tcsr1,info,&
& colcnv=.true.,rowscale=.true.,outcol_glob=.true.,col_desc=desc_c,data=data)
nnz = tcsr1%get_nzeros()
if (update_desc_c) then
call desc_c%indxmap%g2lip_ins(tcsr1%ja(1:nnz),info)
else
call desc_c%indxmap%g2lip(tcsr1%ja(1:nnz),info)
end if
if (info == psb_success_) call psb_rwextd(nacol,bcsr,info,b=tcsr1)
if (info == psb_success_) call tcsr1%free()
if(info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,a_err='Extend am3')
goto 9999
end if
nccol = desc_c%get_local_cols()
call bcsr%set_ncols(nccol)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'starting spspmm 3'
if (debug_level >= psb_debug_outer_) write(debug_unit,*) me,' ',trim(name),&
& 'starting spspmm ',acsr%get_nrows(),acsr%get_ncols(),bcsr%get_nrows(),bcsr%get_ncols()
call psb_spspmm(acsr,bcsr,ccsr,info)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(ictxt,err_act)
return
End Subroutine psb_lz_par_csr_spspmm

@ -479,11 +479,6 @@ dnl FDEFINES="$psblas_cv_define_prepend-DMPI_MOD_F08 $FDEFINES";
fi
fi
dnl PAC_ARG_LONG_INTEGERS
dnl if test x"$pac_cv_long_integers" == x"yes" ; then
dnl FDEFINES="$psblas_cv_define_prepend-DLONG_INTEGERS $FDEFINES";
dnl CDEFINES="-DLONG_INTEGERS_ $CDEFINES";
dnl fi
PAC_ARG_WITH_IPK
PAC_ARG_WITH_LPK
# Defaults for IPK/LPK

Loading…
Cancel
Save