Fixed conversion bug, changed SPASB interface

non-diag
Salvatore Filippone 2 years ago
parent f09e25524e
commit 86b8a261ef

@ -250,7 +250,7 @@ Module psb_c_tools_mod
end interface end interface
interface psb_spasb interface psb_spasb
subroutine psb_cspasb(a,desc_a, info, afmt, upd, mold) subroutine psb_cspasb(a,desc_a, info, afmt, upd, mold, bld_and)
import import
implicit none implicit none
type(psb_cspmat_type), intent (inout) :: a type(psb_cspmat_type), intent (inout) :: a
@ -259,6 +259,7 @@ Module psb_c_tools_mod
integer(psb_ipk_),optional, intent(in) :: upd integer(psb_ipk_),optional, intent(in) :: upd
character(len=*), optional, intent(in) :: afmt character(len=*), optional, intent(in) :: afmt
class(psb_c_base_sparse_mat), intent(in), optional :: mold class(psb_c_base_sparse_mat), intent(in), optional :: mold
logical, intent(in), optional :: bld_and
end subroutine psb_cspasb end subroutine psb_cspasb
end interface end interface

@ -250,7 +250,7 @@ Module psb_d_tools_mod
end interface end interface
interface psb_spasb interface psb_spasb
subroutine psb_dspasb(a,desc_a, info, afmt, upd, mold) subroutine psb_dspasb(a,desc_a, info, afmt, upd, mold, bld_and)
import import
implicit none implicit none
type(psb_dspmat_type), intent (inout) :: a type(psb_dspmat_type), intent (inout) :: a
@ -259,6 +259,7 @@ Module psb_d_tools_mod
integer(psb_ipk_),optional, intent(in) :: upd integer(psb_ipk_),optional, intent(in) :: upd
character(len=*), optional, intent(in) :: afmt character(len=*), optional, intent(in) :: afmt
class(psb_d_base_sparse_mat), intent(in), optional :: mold class(psb_d_base_sparse_mat), intent(in), optional :: mold
logical, intent(in), optional :: bld_and
end subroutine psb_dspasb end subroutine psb_dspasb
end interface end interface

@ -250,7 +250,7 @@ Module psb_s_tools_mod
end interface end interface
interface psb_spasb interface psb_spasb
subroutine psb_sspasb(a,desc_a, info, afmt, upd, mold) subroutine psb_sspasb(a,desc_a, info, afmt, upd, mold, bld_and)
import import
implicit none implicit none
type(psb_sspmat_type), intent (inout) :: a type(psb_sspmat_type), intent (inout) :: a
@ -259,6 +259,7 @@ Module psb_s_tools_mod
integer(psb_ipk_),optional, intent(in) :: upd integer(psb_ipk_),optional, intent(in) :: upd
character(len=*), optional, intent(in) :: afmt character(len=*), optional, intent(in) :: afmt
class(psb_s_base_sparse_mat), intent(in), optional :: mold class(psb_s_base_sparse_mat), intent(in), optional :: mold
logical, intent(in), optional :: bld_and
end subroutine psb_sspasb end subroutine psb_sspasb
end interface end interface

@ -250,7 +250,7 @@ Module psb_z_tools_mod
end interface end interface
interface psb_spasb interface psb_spasb
subroutine psb_zspasb(a,desc_a, info, afmt, upd, mold) subroutine psb_zspasb(a,desc_a, info, afmt, upd, mold, bld_and)
import import
implicit none implicit none
type(psb_zspmat_type), intent (inout) :: a type(psb_zspmat_type), intent (inout) :: a
@ -259,6 +259,7 @@ Module psb_z_tools_mod
integer(psb_ipk_),optional, intent(in) :: upd integer(psb_ipk_),optional, intent(in) :: upd
character(len=*), optional, intent(in) :: afmt character(len=*), optional, intent(in) :: afmt
class(psb_z_base_sparse_mat), intent(in), optional :: mold class(psb_z_base_sparse_mat), intent(in), optional :: mold
logical, intent(in), optional :: bld_and
end subroutine psb_zspasb end subroutine psb_zspasb
end interface end interface

@ -179,11 +179,11 @@ subroutine psb_cspmv_vect(alpha,a,x,beta,y,desc_a,info,&
if (trans_ == 'N') then if (trans_ == 'N') then
! Matrix is not transposed ! Matrix is not transposed
if (.true.) then if (allocated(a%ad)) then
call psi_swapdata(psb_swap_send_,& if (doswap_) call psi_swapdata(psb_swap_send_,&
& czero,x%v,desc_a,iwork,info,data=psb_comm_halo_) & czero,x%v,desc_a,iwork,info,data=psb_comm_halo_)
call a%ad%spmm(alpha,x%v,beta,y%v,info) call a%ad%spmm(alpha,x%v,beta,y%v,info)
call psi_swapdata(psb_swap_recv_,& if (doswap_) call psi_swapdata(psb_swap_recv_,&
& czero,x%v,desc_a,iwork,info,data=psb_comm_halo_) & czero,x%v,desc_a,iwork,info,data=psb_comm_halo_)
call a%and%spmm(alpha,x%v,cone,y%v,info) call a%and%spmm(alpha,x%v,cone,y%v,info)

@ -179,11 +179,11 @@ subroutine psb_dspmv_vect(alpha,a,x,beta,y,desc_a,info,&
if (trans_ == 'N') then if (trans_ == 'N') then
! Matrix is not transposed ! Matrix is not transposed
if (.true.) then if (allocated(a%ad)) then
call psi_swapdata(psb_swap_send_,& if (doswap_) call psi_swapdata(psb_swap_send_,&
& dzero,x%v,desc_a,iwork,info,data=psb_comm_halo_) & dzero,x%v,desc_a,iwork,info,data=psb_comm_halo_)
call a%ad%spmm(alpha,x%v,beta,y%v,info) call a%ad%spmm(alpha,x%v,beta,y%v,info)
call psi_swapdata(psb_swap_recv_,& if (doswap_) call psi_swapdata(psb_swap_recv_,&
& dzero,x%v,desc_a,iwork,info,data=psb_comm_halo_) & dzero,x%v,desc_a,iwork,info,data=psb_comm_halo_)
call a%and%spmm(alpha,x%v,done,y%v,info) call a%and%spmm(alpha,x%v,done,y%v,info)

@ -179,11 +179,11 @@ subroutine psb_sspmv_vect(alpha,a,x,beta,y,desc_a,info,&
if (trans_ == 'N') then if (trans_ == 'N') then
! Matrix is not transposed ! Matrix is not transposed
if (.true.) then if (allocated(a%ad)) then
call psi_swapdata(psb_swap_send_,& if (doswap_) call psi_swapdata(psb_swap_send_,&
& szero,x%v,desc_a,iwork,info,data=psb_comm_halo_) & szero,x%v,desc_a,iwork,info,data=psb_comm_halo_)
call a%ad%spmm(alpha,x%v,beta,y%v,info) call a%ad%spmm(alpha,x%v,beta,y%v,info)
call psi_swapdata(psb_swap_recv_,& if (doswap_) call psi_swapdata(psb_swap_recv_,&
& szero,x%v,desc_a,iwork,info,data=psb_comm_halo_) & szero,x%v,desc_a,iwork,info,data=psb_comm_halo_)
call a%and%spmm(alpha,x%v,sone,y%v,info) call a%and%spmm(alpha,x%v,sone,y%v,info)

@ -179,11 +179,11 @@ subroutine psb_zspmv_vect(alpha,a,x,beta,y,desc_a,info,&
if (trans_ == 'N') then if (trans_ == 'N') then
! Matrix is not transposed ! Matrix is not transposed
if (.true.) then if (allocated(a%ad)) then
call psi_swapdata(psb_swap_send_,& if (doswap_) call psi_swapdata(psb_swap_send_,&
& zzero,x%v,desc_a,iwork,info,data=psb_comm_halo_) & zzero,x%v,desc_a,iwork,info,data=psb_comm_halo_)
call a%ad%spmm(alpha,x%v,beta,y%v,info) call a%ad%spmm(alpha,x%v,beta,y%v,info)
call psi_swapdata(psb_swap_recv_,& if (doswap_) call psi_swapdata(psb_swap_recv_,&
& zzero,x%v,desc_a,iwork,info,data=psb_comm_halo_) & zzero,x%v,desc_a,iwork,info,data=psb_comm_halo_)
call a%and%spmm(alpha,x%v,zone,y%v,info) call a%and%spmm(alpha,x%v,zone,y%v,info)

@ -3643,9 +3643,8 @@ subroutine psb_c_ecsr_csmv(alpha,a,x,beta,y,info,trans)
goto 9999 goto 9999
end if end if
if (((beta == cone).and..not.(tra.or.ctra))& if ((beta == cone).and.&
& .or.(a%is_triangle()).or.(a%is_unit())) then & .not.(tra.or.ctra.or.(a%is_triangle()).or.(a%is_unit()))) then
call psb_c_ecsr_csmv_inner(m,n,alpha,a%irp,a%ja,a%val,& call psb_c_ecsr_csmv_inner(m,n,alpha,a%irp,a%ja,a%val,&
& a%nnerws,a%nerwp,x,y) & a%nnerws,a%nerwp,x,y)
else else
@ -3672,9 +3671,6 @@ contains
if (alpha == czero) return if (alpha == czero) return
if (alpha == cone) then if (alpha == cone) then
!$omp parallel do private(ir,i,j,acc) !$omp parallel do private(ir,i,j,acc)
do ir=1,nnerws do ir=1,nnerws
@ -3740,6 +3736,7 @@ subroutine psb_c_ecsr_cmp_nerwp(a,info)
end if end if
end do end do
call psb_realloc(nnerws,a%nerwp,info) call psb_realloc(nnerws,a%nerwp,info)
a%nnerws = nnerws
end subroutine psb_c_ecsr_cmp_nerwp end subroutine psb_c_ecsr_cmp_nerwp
subroutine psb_c_cp_ecsr_from_coo(a,b,info) subroutine psb_c_cp_ecsr_from_coo(a,b,info)

@ -3643,9 +3643,8 @@ subroutine psb_d_ecsr_csmv(alpha,a,x,beta,y,info,trans)
goto 9999 goto 9999
end if end if
if (((beta == done).and..not.(tra.or.ctra))& if ((beta == done).and.&
& .or.(a%is_triangle()).or.(a%is_unit())) then & .not.(tra.or.ctra.or.(a%is_triangle()).or.(a%is_unit()))) then
call psb_d_ecsr_csmv_inner(m,n,alpha,a%irp,a%ja,a%val,& call psb_d_ecsr_csmv_inner(m,n,alpha,a%irp,a%ja,a%val,&
& a%nnerws,a%nerwp,x,y) & a%nnerws,a%nerwp,x,y)
else else
@ -3672,9 +3671,6 @@ contains
if (alpha == dzero) return if (alpha == dzero) return
if (alpha == done) then if (alpha == done) then
!$omp parallel do private(ir,i,j,acc) !$omp parallel do private(ir,i,j,acc)
do ir=1,nnerws do ir=1,nnerws
@ -3740,6 +3736,7 @@ subroutine psb_d_ecsr_cmp_nerwp(a,info)
end if end if
end do end do
call psb_realloc(nnerws,a%nerwp,info) call psb_realloc(nnerws,a%nerwp,info)
a%nnerws = nnerws
end subroutine psb_d_ecsr_cmp_nerwp end subroutine psb_d_ecsr_cmp_nerwp
subroutine psb_d_cp_ecsr_from_coo(a,b,info) subroutine psb_d_cp_ecsr_from_coo(a,b,info)

@ -3643,9 +3643,8 @@ subroutine psb_s_ecsr_csmv(alpha,a,x,beta,y,info,trans)
goto 9999 goto 9999
end if end if
if (((beta == sone).and..not.(tra.or.ctra))& if ((beta == sone).and.&
& .or.(a%is_triangle()).or.(a%is_unit())) then & .not.(tra.or.ctra.or.(a%is_triangle()).or.(a%is_unit()))) then
call psb_s_ecsr_csmv_inner(m,n,alpha,a%irp,a%ja,a%val,& call psb_s_ecsr_csmv_inner(m,n,alpha,a%irp,a%ja,a%val,&
& a%nnerws,a%nerwp,x,y) & a%nnerws,a%nerwp,x,y)
else else
@ -3672,9 +3671,6 @@ contains
if (alpha == szero) return if (alpha == szero) return
if (alpha == sone) then if (alpha == sone) then
!$omp parallel do private(ir,i,j,acc) !$omp parallel do private(ir,i,j,acc)
do ir=1,nnerws do ir=1,nnerws
@ -3740,6 +3736,7 @@ subroutine psb_s_ecsr_cmp_nerwp(a,info)
end if end if
end do end do
call psb_realloc(nnerws,a%nerwp,info) call psb_realloc(nnerws,a%nerwp,info)
a%nnerws = nnerws
end subroutine psb_s_ecsr_cmp_nerwp end subroutine psb_s_ecsr_cmp_nerwp
subroutine psb_s_cp_ecsr_from_coo(a,b,info) subroutine psb_s_cp_ecsr_from_coo(a,b,info)

@ -3643,9 +3643,8 @@ subroutine psb_z_ecsr_csmv(alpha,a,x,beta,y,info,trans)
goto 9999 goto 9999
end if end if
if (((beta == zone).and..not.(tra.or.ctra))& if ((beta == zone).and.&
& .or.(a%is_triangle()).or.(a%is_unit())) then & .not.(tra.or.ctra.or.(a%is_triangle()).or.(a%is_unit()))) then
call psb_z_ecsr_csmv_inner(m,n,alpha,a%irp,a%ja,a%val,& call psb_z_ecsr_csmv_inner(m,n,alpha,a%irp,a%ja,a%val,&
& a%nnerws,a%nerwp,x,y) & a%nnerws,a%nerwp,x,y)
else else
@ -3672,9 +3671,6 @@ contains
if (alpha == zzero) return if (alpha == zzero) return
if (alpha == zone) then if (alpha == zone) then
!$omp parallel do private(ir,i,j,acc) !$omp parallel do private(ir,i,j,acc)
do ir=1,nnerws do ir=1,nnerws
@ -3740,6 +3736,7 @@ subroutine psb_z_ecsr_cmp_nerwp(a,info)
end if end if
end do end do
call psb_realloc(nnerws,a%nerwp,info) call psb_realloc(nnerws,a%nerwp,info)
a%nnerws = nnerws
end subroutine psb_z_ecsr_cmp_nerwp end subroutine psb_z_ecsr_cmp_nerwp
subroutine psb_z_cp_ecsr_from_coo(a,b,info) subroutine psb_z_cp_ecsr_from_coo(a,b,info)

@ -44,7 +44,7 @@
! psb_upd_perm_ Permutation(more memory) ! psb_upd_perm_ Permutation(more memory)
! !
! !
subroutine psb_cspasb(a,desc_a, info, afmt, upd, mold) subroutine psb_cspasb(a,desc_a, info, afmt, upd, mold, bld_and)
use psb_base_mod, psb_protect_name => psb_cspasb use psb_base_mod, psb_protect_name => psb_cspasb
use psb_sort_mod use psb_sort_mod
use psi_mod use psi_mod
@ -58,6 +58,7 @@ subroutine psb_cspasb(a,desc_a, info, afmt, upd, mold)
integer(psb_ipk_), optional, intent(in) :: upd integer(psb_ipk_), optional, intent(in) :: upd
character(len=*), optional, intent(in) :: afmt character(len=*), optional, intent(in) :: afmt
class(psb_c_base_sparse_mat), intent(in), optional :: mold class(psb_c_base_sparse_mat), intent(in), optional :: mold
logical, intent(in), optional :: bld_and
!....Locals.... !....Locals....
type(psb_ctxt_type) :: ctxt type(psb_ctxt_type) :: ctxt
integer(psb_ipk_) :: np,me, err_act integer(psb_ipk_) :: np,me, err_act
@ -65,6 +66,7 @@ subroutine psb_cspasb(a,desc_a, info, afmt, upd, mold)
integer(psb_ipk_) :: debug_level, debug_unit integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err character(len=20) :: name, ch_err
class(psb_i_base_vect_type), allocatable :: ivm class(psb_i_base_vect_type), allocatable :: ivm
logical :: bld_and_
info = psb_success_ info = psb_success_
name = 'psb_spasb' name = 'psb_spasb'
@ -93,7 +95,11 @@ subroutine psb_cspasb(a,desc_a, info, afmt, upd, mold)
if (debug_level >= psb_debug_ext_)& if (debug_level >= psb_debug_ext_)&
& write(debug_unit, *) me,' ',trim(name),& & write(debug_unit, *) me,' ',trim(name),&
& ' Begin matrix assembly...' & ' Begin matrix assembly...'
if (present(bld_and)) then
bld_and_ = bld_and
else
bld_and_ = .false.
end if
!check on errors encountered in psdspins !check on errors encountered in psdspins
if (a%is_bld()) then if (a%is_bld()) then
@ -171,19 +177,26 @@ subroutine psb_cspasb(a,desc_a, info, afmt, upd, mold)
end if end if
if (.true.) then if (bld_and_) then
block block
character(len=1024) :: fname character(len=1024) :: fname
type(psb_c_coo_sparse_mat) :: acoo type(psb_c_coo_sparse_mat) :: acoo
type(psb_c_csr_sparse_mat), allocatable :: aclip type(psb_c_csr_sparse_mat), allocatable :: aclip
type(psb_c_ecsr_sparse_mat), allocatable :: andclip type(psb_c_ecsr_sparse_mat), allocatable :: andclip
allocate(aclip,andclip) logical, parameter :: use_ecsr=.false.
allocate(aclip)
call a%a%csclip(acoo,info,jmax=n_row,rscale=.false.,cscale=.false.) call a%a%csclip(acoo,info,jmax=n_row,rscale=.false.,cscale=.false.)
allocate(a%ad,mold=a%a) allocate(a%ad,mold=a%a)
call a%ad%mv_from_coo(acoo,info) call a%ad%mv_from_coo(acoo,info)
call a%a%csclip(acoo,info,jmin=n_row+1,jmax=n_col,rscale=.false.,cscale=.false.) call a%a%csclip(acoo,info,jmin=n_row+1,jmax=n_col,rscale=.false.,cscale=.false.)
call andclip%mv_from_coo(acoo,info) if (use_ecsr) then
call move_alloc(andclip,a%and) allocate(andclip)
call andclip%mv_from_coo(acoo,info)
call move_alloc(andclip,a%and)
else
allocate(a%and,mold=a%a)
call a%and%mv_from_coo(acoo,info)
end if
if (.false.) then if (.false.) then
write(fname,'(a,i2.2,a)') 'adclip_',me,'.mtx' write(fname,'(a,i2.2,a)') 'adclip_',me,'.mtx'
open(25,file=fname) open(25,file=fname)
@ -200,6 +213,9 @@ subroutine psb_cspasb(a,desc_a, info, afmt, upd, mold)
&a%and%get_nrows(),a%and%get_ncols(),n_row,n_col &a%and%get_nrows(),a%and%get_ncols(),n_row,n_col
end if end if
end block end block
else
if (allocated(a%ad)) deallocate(a%ad)
if (allocated(a%and)) deallocate(a%and)
end if end if
if (debug_level >= psb_debug_ext_) then if (debug_level >= psb_debug_ext_) then
ch_err=a%get_fmt() ch_err=a%get_fmt()

@ -44,7 +44,7 @@
! psb_upd_perm_ Permutation(more memory) ! psb_upd_perm_ Permutation(more memory)
! !
! !
subroutine psb_dspasb(a,desc_a, info, afmt, upd, mold) subroutine psb_dspasb(a,desc_a, info, afmt, upd, mold, bld_and)
use psb_base_mod, psb_protect_name => psb_dspasb use psb_base_mod, psb_protect_name => psb_dspasb
use psb_sort_mod use psb_sort_mod
use psi_mod use psi_mod
@ -58,6 +58,7 @@ subroutine psb_dspasb(a,desc_a, info, afmt, upd, mold)
integer(psb_ipk_), optional, intent(in) :: upd integer(psb_ipk_), optional, intent(in) :: upd
character(len=*), optional, intent(in) :: afmt character(len=*), optional, intent(in) :: afmt
class(psb_d_base_sparse_mat), intent(in), optional :: mold class(psb_d_base_sparse_mat), intent(in), optional :: mold
logical, intent(in), optional :: bld_and
!....Locals.... !....Locals....
type(psb_ctxt_type) :: ctxt type(psb_ctxt_type) :: ctxt
integer(psb_ipk_) :: np,me, err_act integer(psb_ipk_) :: np,me, err_act
@ -65,6 +66,7 @@ subroutine psb_dspasb(a,desc_a, info, afmt, upd, mold)
integer(psb_ipk_) :: debug_level, debug_unit integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err character(len=20) :: name, ch_err
class(psb_i_base_vect_type), allocatable :: ivm class(psb_i_base_vect_type), allocatable :: ivm
logical :: bld_and_
info = psb_success_ info = psb_success_
name = 'psb_spasb' name = 'psb_spasb'
@ -93,7 +95,11 @@ subroutine psb_dspasb(a,desc_a, info, afmt, upd, mold)
if (debug_level >= psb_debug_ext_)& if (debug_level >= psb_debug_ext_)&
& write(debug_unit, *) me,' ',trim(name),& & write(debug_unit, *) me,' ',trim(name),&
& ' Begin matrix assembly...' & ' Begin matrix assembly...'
if (present(bld_and)) then
bld_and_ = bld_and
else
bld_and_ = .false.
end if
!check on errors encountered in psdspins !check on errors encountered in psdspins
if (a%is_bld()) then if (a%is_bld()) then
@ -171,19 +177,26 @@ subroutine psb_dspasb(a,desc_a, info, afmt, upd, mold)
end if end if
if (.true.) then if (bld_and_) then
block block
character(len=1024) :: fname character(len=1024) :: fname
type(psb_d_coo_sparse_mat) :: acoo type(psb_d_coo_sparse_mat) :: acoo
type(psb_d_csr_sparse_mat), allocatable :: aclip type(psb_d_csr_sparse_mat), allocatable :: aclip
type(psb_d_ecsr_sparse_mat), allocatable :: andclip type(psb_d_ecsr_sparse_mat), allocatable :: andclip
allocate(aclip,andclip) logical, parameter :: use_ecsr=.true.
allocate(aclip)
call a%a%csclip(acoo,info,jmax=n_row,rscale=.false.,cscale=.false.) call a%a%csclip(acoo,info,jmax=n_row,rscale=.false.,cscale=.false.)
allocate(a%ad,mold=a%a) allocate(a%ad,mold=a%a)
call a%ad%mv_from_coo(acoo,info) call a%ad%mv_from_coo(acoo,info)
call a%a%csclip(acoo,info,jmin=n_row+1,jmax=n_col,rscale=.false.,cscale=.false.) call a%a%csclip(acoo,info,jmin=n_row+1,jmax=n_col,rscale=.false.,cscale=.false.)
call andclip%mv_from_coo(acoo,info) if (use_ecsr) then
call move_alloc(andclip,a%and) allocate(andclip)
call andclip%mv_from_coo(acoo,info)
call move_alloc(andclip,a%and)
else
allocate(a%and,mold=a%a)
call a%and%mv_from_coo(acoo,info)
end if
if (.false.) then if (.false.) then
write(fname,'(a,i2.2,a)') 'adclip_',me,'.mtx' write(fname,'(a,i2.2,a)') 'adclip_',me,'.mtx'
open(25,file=fname) open(25,file=fname)
@ -200,6 +213,9 @@ subroutine psb_dspasb(a,desc_a, info, afmt, upd, mold)
&a%and%get_nrows(),a%and%get_ncols(),n_row,n_col &a%and%get_nrows(),a%and%get_ncols(),n_row,n_col
end if end if
end block end block
else
if (allocated(a%ad)) deallocate(a%ad)
if (allocated(a%and)) deallocate(a%and)
end if end if
if (debug_level >= psb_debug_ext_) then if (debug_level >= psb_debug_ext_) then
ch_err=a%get_fmt() ch_err=a%get_fmt()

@ -44,7 +44,7 @@
! psb_upd_perm_ Permutation(more memory) ! psb_upd_perm_ Permutation(more memory)
! !
! !
subroutine psb_sspasb(a,desc_a, info, afmt, upd, mold) subroutine psb_sspasb(a,desc_a, info, afmt, upd, mold, bld_and)
use psb_base_mod, psb_protect_name => psb_sspasb use psb_base_mod, psb_protect_name => psb_sspasb
use psb_sort_mod use psb_sort_mod
use psi_mod use psi_mod
@ -58,6 +58,7 @@ subroutine psb_sspasb(a,desc_a, info, afmt, upd, mold)
integer(psb_ipk_), optional, intent(in) :: upd integer(psb_ipk_), optional, intent(in) :: upd
character(len=*), optional, intent(in) :: afmt character(len=*), optional, intent(in) :: afmt
class(psb_s_base_sparse_mat), intent(in), optional :: mold class(psb_s_base_sparse_mat), intent(in), optional :: mold
logical, intent(in), optional :: bld_and
!....Locals.... !....Locals....
type(psb_ctxt_type) :: ctxt type(psb_ctxt_type) :: ctxt
integer(psb_ipk_) :: np,me, err_act integer(psb_ipk_) :: np,me, err_act
@ -65,6 +66,7 @@ subroutine psb_sspasb(a,desc_a, info, afmt, upd, mold)
integer(psb_ipk_) :: debug_level, debug_unit integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err character(len=20) :: name, ch_err
class(psb_i_base_vect_type), allocatable :: ivm class(psb_i_base_vect_type), allocatable :: ivm
logical :: bld_and_
info = psb_success_ info = psb_success_
name = 'psb_spasb' name = 'psb_spasb'
@ -93,7 +95,11 @@ subroutine psb_sspasb(a,desc_a, info, afmt, upd, mold)
if (debug_level >= psb_debug_ext_)& if (debug_level >= psb_debug_ext_)&
& write(debug_unit, *) me,' ',trim(name),& & write(debug_unit, *) me,' ',trim(name),&
& ' Begin matrix assembly...' & ' Begin matrix assembly...'
if (present(bld_and)) then
bld_and_ = bld_and
else
bld_and_ = .false.
end if
!check on errors encountered in psdspins !check on errors encountered in psdspins
if (a%is_bld()) then if (a%is_bld()) then
@ -171,19 +177,26 @@ subroutine psb_sspasb(a,desc_a, info, afmt, upd, mold)
end if end if
if (.true.) then if (bld_and_) then
block block
character(len=1024) :: fname character(len=1024) :: fname
type(psb_s_coo_sparse_mat) :: acoo type(psb_s_coo_sparse_mat) :: acoo
type(psb_s_csr_sparse_mat), allocatable :: aclip type(psb_s_csr_sparse_mat), allocatable :: aclip
type(psb_s_ecsr_sparse_mat), allocatable :: andclip type(psb_s_ecsr_sparse_mat), allocatable :: andclip
allocate(aclip,andclip) logical, parameter :: use_ecsr=.false.
allocate(aclip)
call a%a%csclip(acoo,info,jmax=n_row,rscale=.false.,cscale=.false.) call a%a%csclip(acoo,info,jmax=n_row,rscale=.false.,cscale=.false.)
allocate(a%ad,mold=a%a) allocate(a%ad,mold=a%a)
call a%ad%mv_from_coo(acoo,info) call a%ad%mv_from_coo(acoo,info)
call a%a%csclip(acoo,info,jmin=n_row+1,jmax=n_col,rscale=.false.,cscale=.false.) call a%a%csclip(acoo,info,jmin=n_row+1,jmax=n_col,rscale=.false.,cscale=.false.)
call andclip%mv_from_coo(acoo,info) if (use_ecsr) then
call move_alloc(andclip,a%and) allocate(andclip)
call andclip%mv_from_coo(acoo,info)
call move_alloc(andclip,a%and)
else
allocate(a%and,mold=a%a)
call a%and%mv_from_coo(acoo,info)
end if
if (.false.) then if (.false.) then
write(fname,'(a,i2.2,a)') 'adclip_',me,'.mtx' write(fname,'(a,i2.2,a)') 'adclip_',me,'.mtx'
open(25,file=fname) open(25,file=fname)
@ -200,6 +213,9 @@ subroutine psb_sspasb(a,desc_a, info, afmt, upd, mold)
&a%and%get_nrows(),a%and%get_ncols(),n_row,n_col &a%and%get_nrows(),a%and%get_ncols(),n_row,n_col
end if end if
end block end block
else
if (allocated(a%ad)) deallocate(a%ad)
if (allocated(a%and)) deallocate(a%and)
end if end if
if (debug_level >= psb_debug_ext_) then if (debug_level >= psb_debug_ext_) then
ch_err=a%get_fmt() ch_err=a%get_fmt()

@ -44,7 +44,7 @@
! psb_upd_perm_ Permutation(more memory) ! psb_upd_perm_ Permutation(more memory)
! !
! !
subroutine psb_zspasb(a,desc_a, info, afmt, upd, mold) subroutine psb_zspasb(a,desc_a, info, afmt, upd, mold, bld_and)
use psb_base_mod, psb_protect_name => psb_zspasb use psb_base_mod, psb_protect_name => psb_zspasb
use psb_sort_mod use psb_sort_mod
use psi_mod use psi_mod
@ -58,6 +58,7 @@ subroutine psb_zspasb(a,desc_a, info, afmt, upd, mold)
integer(psb_ipk_), optional, intent(in) :: upd integer(psb_ipk_), optional, intent(in) :: upd
character(len=*), optional, intent(in) :: afmt character(len=*), optional, intent(in) :: afmt
class(psb_z_base_sparse_mat), intent(in), optional :: mold class(psb_z_base_sparse_mat), intent(in), optional :: mold
logical, intent(in), optional :: bld_and
!....Locals.... !....Locals....
type(psb_ctxt_type) :: ctxt type(psb_ctxt_type) :: ctxt
integer(psb_ipk_) :: np,me, err_act integer(psb_ipk_) :: np,me, err_act
@ -65,6 +66,7 @@ subroutine psb_zspasb(a,desc_a, info, afmt, upd, mold)
integer(psb_ipk_) :: debug_level, debug_unit integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err character(len=20) :: name, ch_err
class(psb_i_base_vect_type), allocatable :: ivm class(psb_i_base_vect_type), allocatable :: ivm
logical :: bld_and_
info = psb_success_ info = psb_success_
name = 'psb_spasb' name = 'psb_spasb'
@ -93,7 +95,11 @@ subroutine psb_zspasb(a,desc_a, info, afmt, upd, mold)
if (debug_level >= psb_debug_ext_)& if (debug_level >= psb_debug_ext_)&
& write(debug_unit, *) me,' ',trim(name),& & write(debug_unit, *) me,' ',trim(name),&
& ' Begin matrix assembly...' & ' Begin matrix assembly...'
if (present(bld_and)) then
bld_and_ = bld_and
else
bld_and_ = .false.
end if
!check on errors encountered in psdspins !check on errors encountered in psdspins
if (a%is_bld()) then if (a%is_bld()) then
@ -171,19 +177,26 @@ subroutine psb_zspasb(a,desc_a, info, afmt, upd, mold)
end if end if
if (.true.) then if (bld_and_) then
block block
character(len=1024) :: fname character(len=1024) :: fname
type(psb_z_coo_sparse_mat) :: acoo type(psb_z_coo_sparse_mat) :: acoo
type(psb_z_csr_sparse_mat), allocatable :: aclip type(psb_z_csr_sparse_mat), allocatable :: aclip
type(psb_z_ecsr_sparse_mat), allocatable :: andclip type(psb_z_ecsr_sparse_mat), allocatable :: andclip
allocate(aclip,andclip) logical, parameter :: use_ecsr=.false.
allocate(aclip)
call a%a%csclip(acoo,info,jmax=n_row,rscale=.false.,cscale=.false.) call a%a%csclip(acoo,info,jmax=n_row,rscale=.false.,cscale=.false.)
allocate(a%ad,mold=a%a) allocate(a%ad,mold=a%a)
call a%ad%mv_from_coo(acoo,info) call a%ad%mv_from_coo(acoo,info)
call a%a%csclip(acoo,info,jmin=n_row+1,jmax=n_col,rscale=.false.,cscale=.false.) call a%a%csclip(acoo,info,jmin=n_row+1,jmax=n_col,rscale=.false.,cscale=.false.)
call andclip%mv_from_coo(acoo,info) if (use_ecsr) then
call move_alloc(andclip,a%and) allocate(andclip)
call andclip%mv_from_coo(acoo,info)
call move_alloc(andclip,a%and)
else
allocate(a%and,mold=a%a)
call a%and%mv_from_coo(acoo,info)
end if
if (.false.) then if (.false.) then
write(fname,'(a,i2.2,a)') 'adclip_',me,'.mtx' write(fname,'(a,i2.2,a)') 'adclip_',me,'.mtx'
open(25,file=fname) open(25,file=fname)
@ -200,6 +213,9 @@ subroutine psb_zspasb(a,desc_a, info, afmt, upd, mold)
&a%and%get_nrows(),a%and%get_ncols(),n_row,n_col &a%and%get_nrows(),a%and%get_ncols(),n_row,n_col
end if end if
end block end block
else
if (allocated(a%ad)) deallocate(a%ad)
if (allocated(a%and)) deallocate(a%and)
end if end if
if (debug_level >= psb_debug_ext_) then if (debug_level >= psb_debug_ext_) then
ch_err=a%get_fmt() ch_err=a%get_fmt()

@ -680,9 +680,9 @@ contains
t1 = psb_wtime() t1 = psb_wtime()
if (info == psb_success_) then if (info == psb_success_) then
if (present(amold)) then if (present(amold)) then
call psb_spasb(a,desc_a,info,mold=amold) call psb_spasb(a,desc_a,info,mold=amold,bld_and=.true.)
else else
call psb_spasb(a,desc_a,info,afmt=afmt) call psb_spasb(a,desc_a,info,afmt=afmt,bld_and=.true.)
end if end if
end if end if
call psb_barrier(ctxt) call psb_barrier(ctxt)

@ -2,11 +2,11 @@
BICGSTAB Iterative method BICGSTAB CGS BICG BICGSTABL RGMRES FCG CGR BICGSTAB Iterative method BICGSTAB CGS BICG BICGSTABL RGMRES FCG CGR
BJAC Preconditioner NONE DIAG BJAC BJAC Preconditioner NONE DIAG BJAC
CSR Storage format for matrix A: CSR COO CSR Storage format for matrix A: CSR COO
100 Domain size (acutal system is this**3 (pde3d) or **2 (pde2d) ) 200 Domain size (acutal system is this**3 (pde3d) or **2 (pde2d) )
3 Partition: 1 BLOCK 3 3D 3 Partition: 1 BLOCK 3 3D
2 Stopping criterion 1 2 2 Stopping criterion 1 2
0100 MAXIT 0300 MAXIT
05 ITRACE 10 ITRACE
002 IRST restart for RGMRES and BiCGSTABL 002 IRST restart for RGMRES and BiCGSTABL
ILU Block Solver ILU,ILUT,INVK,AINVT,AORTH ILU Block Solver ILU,ILUT,INVK,AINVT,AORTH
NONE If ILU : MILU or NONE othewise ignored NONE If ILU : MILU or NONE othewise ignored

Loading…
Cancel
Save