fdurastante
/
amg4psblas
mirror of https://github.com/sfilippone/amg4psblas.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

mld2p4-extaggr:

Browse Source 
mlprec/impl/mld_daggrmat_asb.f90
 mlprec/impl/mld_dcoarse_bld.f90
 mlprec/mld_d_inner_mod.f90

Refactoring steps: work split between coarse_bld and internals.
stopcriterion
Salvatore Filippone 9 years ago
parent 26119298bd
commit 37943c7b98
3 changed files with 249 additions and 132 deletions
Show Stats Download Patch File Download Diff File

231
mlprec/impl/mld_daggrmat_asb.f90
Unescape Escape View File

@ -98,7 +98,7 @@
! info - integer, output. ! info - integer, output.
! Error code. ! Error code.
! !
subroutine mld_daggrmat_asb(a,desc_a,ilaggr,nlaggr,p,info) subroutine mld_daggrmat_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info)
use psb_base_mod use psb_base_mod
use mld_d_inner_mod, mld_protect_name => mld_daggrmat_asb use mld_d_inner_mod, mld_protect_name => mld_daggrmat_asb
@ -109,11 +109,12 @@ subroutine mld_daggrmat_asb(a,desc_a,ilaggr,nlaggr,p,info)
type(psb_dspmat_type), intent(in) :: a type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a type(psb_desc_type), intent(in) :: desc_a
integer(psb_ipk_), intent(inout) :: ilaggr(:), nlaggr(:) integer(psb_ipk_), intent(inout) :: ilaggr(:), nlaggr(:)
type(mld_d_onelev_type), intent(inout), target :: p type(mld_dml_parms), intent(inout) :: parms
!!$ type(mld_d_onelev_type), intent(inout), target :: p
integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), intent(out) :: info
type(psb_dspmat_type), intent(inout) :: ac, op_prol,op_restr
! Local variables ! Local variables
type(psb_dspmat_type) :: ac, op_prol,op_restr
type(psb_d_coo_sparse_mat) :: acoo, bcoo type(psb_d_coo_sparse_mat) :: acoo, bcoo
type(psb_d_csr_sparse_mat) :: acsr1 type(psb_d_csr_sparse_mat) :: acsr1
integer(psb_ipk_) :: nzl,ntaggr, err_act integer(psb_ipk_) :: nzl,ntaggr, err_act
@ -133,26 +134,26 @@ subroutine mld_daggrmat_asb(a,desc_a,ilaggr,nlaggr,p,info)
call psb_info(ictxt, me, np) call psb_info(ictxt, me, np)
select case (p%parms%aggr_kind) select case (parms%aggr_kind)
case (mld_no_smooth_) case (mld_no_smooth_)
call mld_daggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,& call mld_daggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,&
& p%parms,ac,op_prol,op_restr,info) & parms,ac,op_prol,op_restr,info)
case(mld_smooth_prol_) case(mld_smooth_prol_)
call mld_daggrmat_smth_asb(a,desc_a,ilaggr,nlaggr, & call mld_daggrmat_smth_asb(a,desc_a,ilaggr,nlaggr, &
& p%parms,ac,op_prol,op_restr,info) & parms,ac,op_prol,op_restr,info)
case(mld_biz_prol_) case(mld_biz_prol_)
call mld_daggrmat_biz_asb(a,desc_a,ilaggr,nlaggr, & call mld_daggrmat_biz_asb(a,desc_a,ilaggr,nlaggr, &
& p%parms,ac,op_prol,op_restr,info) & parms,ac,op_prol,op_restr,info)
case(mld_min_energy_) case(mld_min_energy_)
call mld_daggrmat_minnrg_asb(a,desc_a,ilaggr,nlaggr, & call mld_daggrmat_minnrg_asb(a,desc_a,ilaggr,nlaggr, &
& p%parms,ac,op_prol,op_restr,info) & parms,ac,op_prol,op_restr,info)
case default case default
info = psb_err_internal_error_ info = psb_err_internal_error_
@ -166,113 +167,113 @@ subroutine mld_daggrmat_asb(a,desc_a,ilaggr,nlaggr,p,info)
end if end if
!!$
ntaggr = sum(nlaggr) !!$ ntaggr = sum(nlaggr)
!!$
select case(p%parms%coarse_mat) !!$ select case(p%parms%coarse_mat)
!!$
case(mld_distr_mat_) !!$ case(mld_distr_mat_)
!!$
call ac%mv_to(bcoo) !!$ call ac%mv_to(bcoo)
if (p%parms%clean_zeros) call bcoo%clean_zeros(info) !!$ if (p%parms%clean_zeros) call bcoo%clean_zeros(info)
nzl = bcoo%get_nzeros() !!$ nzl = bcoo%get_nzeros()
!!$
if (info == psb_success_) call psb_cdall(ictxt,p%desc_ac,info,nl=nlaggr(me+1)) !!$ if (info == psb_success_) call psb_cdall(ictxt,p%desc_ac,info,nl=nlaggr(me+1))
if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,p%desc_ac,info) !!$ if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,p%desc_ac,info)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info) !!$ if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),p%desc_ac,info,iact='I') !!$ if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),p%desc_ac,info,iact='I')
if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),p%desc_ac,info,iact='I') !!$ if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),p%desc_ac,info,iact='I')
if (info /= psb_success_) then !!$ if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,& !!$ call psb_errpush(psb_err_internal_error_,name,&
& a_err='Creating p%desc_ac and converting ac') !!$ & a_err='Creating p%desc_ac and converting ac')
goto 9999 !!$ goto 9999
end if !!$ end if
if (debug_level >= psb_debug_outer_) & !!$ if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),& !!$ & write(debug_unit,*) me,' ',trim(name),&
& 'Assembld aux descr. distr.' !!$ & 'Assembld aux descr. distr.'
call p%ac%mv_from(bcoo) !!$ call p%ac%mv_from(bcoo)
!!$
call p%ac%set_nrows(p%desc_ac%get_local_rows()) !!$ call p%ac%set_nrows(p%desc_ac%get_local_rows())
call p%ac%set_ncols(p%desc_ac%get_local_cols()) !!$ call p%ac%set_ncols(p%desc_ac%get_local_cols())
call p%ac%set_asb() !!$ call p%ac%set_asb()
!!$
if (info /= psb_success_) then !!$ if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free') !!$ call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free')
goto 9999 !!$ goto 9999
end if !!$ end if
!!$
if (np>1) then !!$ if (np>1) then
call op_prol%mv_to(acsr1) !!$ call op_prol%mv_to(acsr1)
nzl = acsr1%get_nzeros() !!$ nzl = acsr1%get_nzeros()
call psb_glob_to_loc(acsr1%ja(1:nzl),p%desc_ac,info,'I') !!$ call psb_glob_to_loc(acsr1%ja(1:nzl),p%desc_ac,info,'I')
if(info /= psb_success_) then !!$ if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc') !!$ call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc')
goto 9999 !!$ goto 9999
end if !!$ end if
call op_prol%mv_from(acsr1) !!$ call op_prol%mv_from(acsr1)
endif !!$ endif
call op_prol%set_ncols(p%desc_ac%get_local_cols()) !!$ call op_prol%set_ncols(p%desc_ac%get_local_cols())
!!$
if (np>1) then !!$ if (np>1) then
call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_) !!$ call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_)
call op_restr%mv_to(acoo) !!$ call op_restr%mv_to(acoo)
nzl = acoo%get_nzeros() !!$ nzl = acoo%get_nzeros()
if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),p%desc_ac,info,'I') !!$ if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),p%desc_ac,info,'I')
call acoo%set_dupl(psb_dupl_add_) !!$ call acoo%set_dupl(psb_dupl_add_)
if (info == psb_success_) call op_restr%mv_from(acoo) !!$ if (info == psb_success_) call op_restr%mv_from(acoo)
if (info == psb_success_) call op_restr%cscnv(info,type='csr') !!$ if (info == psb_success_) call op_restr%cscnv(info,type='csr')
if(info /= psb_success_) then !!$ if(info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,& !!$ call psb_errpush(psb_err_internal_error_,name,&
& a_err='Converting op_restr to local') !!$ & a_err='Converting op_restr to local')
goto 9999 !!$ goto 9999
end if !!$ end if
end if !!$ end if
! !!$ !
! Clip to local rows. !!$ ! Clip to local rows.
! !!$ !
call op_restr%set_nrows(p%desc_ac%get_local_rows()) !!$ call op_restr%set_nrows(p%desc_ac%get_local_rows())
!!$
if (debug_level >= psb_debug_outer_) & !!$ if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),& !!$ & write(debug_unit,*) me,' ',trim(name),&
& 'Done ac ' !!$ & 'Done ac '
!!$
case(mld_repl_mat_) !!$ case(mld_repl_mat_)
! !!$ !
! !!$ !
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) !!$ call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info) !!$ if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if ((info == psb_success_).and.p%parms%clean_zeros) call ac%clean_zeros(info) !!$ if ((info == psb_success_).and.p%parms%clean_zeros) call ac%clean_zeros(info)
if (info == psb_success_) & !!$ if (info == psb_success_) &
& call psb_gather(p%ac,ac,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) !!$ & call psb_gather(p%ac,ac,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.)
!!$
if (info /= psb_success_) goto 9999 !!$ if (info /= psb_success_) goto 9999
!!$
case default !!$ case default
info = psb_err_internal_error_ !!$ info = psb_err_internal_error_
call psb_errpush(info,name,a_err='invalid mld_coarse_mat_') !!$ call psb_errpush(info,name,a_err='invalid mld_coarse_mat_')
goto 9999 !!$ goto 9999
end select !!$ end select
!!$
call p%ac%cscnv(info,type='csr',dupl=psb_dupl_add_) !!$ call p%ac%cscnv(info,type='csr',dupl=psb_dupl_add_)
if(info /= psb_success_) then !!$ if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='spcnv') !!$ call psb_errpush(psb_err_from_subroutine_,name,a_err='spcnv')
goto 9999 !!$ goto 9999
end if !!$ end if
!!$
! !!$ !
! Copy the prolongation/restriction matrices into the descriptor map. !!$ ! Copy the prolongation/restriction matrices into the descriptor map.
! op_restr => PR^T i.e. restriction operator !!$ ! op_restr => PR^T i.e. restriction operator
! op_prol => PR i.e. prolongation operator !!$ ! op_prol => PR i.e. prolongation operator
! !!$ !
!!$
p%map = psb_linmap(psb_map_aggr_,desc_a,& !!$ p%map = psb_linmap(psb_map_aggr_,desc_a,&
& p%desc_ac,op_restr,op_prol,ilaggr,nlaggr) !!$ & p%desc_ac,op_restr,op_prol,ilaggr,nlaggr)
if (info == psb_success_) call op_prol%free() !!$ if (info == psb_success_) call op_prol%free()
if (info == psb_success_) call op_restr%free() !!$ if (info == psb_success_) call op_restr%free()
if(info /= psb_success_) then !!$ if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='sp_Free') !!$ call psb_errpush(psb_err_from_subroutine_,name,a_err='sp_Free')
goto 9999 !!$ goto 9999
end if !!$ end if
call psb_erractionrestore(err_act) call psb_erractionrestore(err_act)

120
mlprec/impl/mld_dcoarse_bld.f90
Unescape Escape View File

@ -83,11 +83,18 @@ subroutine mld_dcoarse_bld(a,desc_a,p,info)
integer(psb_mpik_) :: ictxt, np, me integer(psb_mpik_) :: ictxt, np, me
integer(psb_ipk_) :: err_act integer(psb_ipk_) :: err_act
integer(psb_ipk_), allocatable :: ilaggr(:), nlaggr(:) integer(psb_ipk_), allocatable :: ilaggr(:), nlaggr(:)
type(psb_dspmat_type) :: ac, op_prol,op_restr
type(psb_d_coo_sparse_mat) :: acoo, bcoo
type(psb_d_csr_sparse_mat) :: acsr1
integer(psb_ipk_) :: nzl, ntaggr
integer(psb_ipk_) :: debug_level, debug_unit
name='mld_dcoarse_bld' name='mld_dcoarse_bld'
if (psb_get_errstatus().ne.0) return if (psb_get_errstatus().ne.0) return
call psb_erractionsave(err_act) call psb_erractionsave(err_act)
info = psb_success_ debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
info = psb_success_
ictxt = desc_a%get_context() ictxt = desc_a%get_context()
call psb_info(ictxt,me,np) call psb_info(ictxt,me,np)
@ -134,7 +141,7 @@ subroutine mld_dcoarse_bld(a,desc_a,p,info)
! the mapping defined by mld_aggrmap_bld and applying the aggregation ! the mapping defined by mld_aggrmap_bld and applying the aggregation
! algorithm specified by p%iprcparm(mld_aggr_kind_) ! algorithm specified by p%iprcparm(mld_aggr_kind_)
! !
call mld_aggrmat_asb(a,desc_a,ilaggr,nlaggr,p,info) call mld_daggrmat_asb(a,desc_a,ilaggr,nlaggr,p%parms,ac,op_prol,op_restr,info)
if(info /= psb_success_) then if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_aggrmat_asb') call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_aggrmat_asb')
@ -157,6 +164,115 @@ subroutine mld_dcoarse_bld(a,desc_a,p,info)
end select end select
! Common code refactored here.
ntaggr = sum(nlaggr)
select case(p%parms%coarse_mat)
case(mld_distr_mat_)
call ac%mv_to(bcoo)
if (p%parms%clean_zeros) call bcoo%clean_zeros(info)
nzl = bcoo%get_nzeros()
if (info == psb_success_) call psb_cdall(ictxt,p%desc_ac,info,nl=nlaggr(me+1))
if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,p%desc_ac,info)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),p%desc_ac,info,iact='I')
if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),p%desc_ac,info,iact='I')
if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,&
& a_err='Creating p%desc_ac and converting ac')
goto 9999
end if
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'Assembld aux descr. distr.'
call p%ac%mv_from(bcoo)
call p%ac%set_nrows(p%desc_ac%get_local_rows())
call p%ac%set_ncols(p%desc_ac%get_local_cols())
call p%ac%set_asb()
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free')
goto 9999
end if
if (np>1) then
call op_prol%mv_to(acsr1)
nzl = acsr1%get_nzeros()
call psb_glob_to_loc(acsr1%ja(1:nzl),p%desc_ac,info,'I')
if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc')
goto 9999
end if
call op_prol%mv_from(acsr1)
endif
call op_prol%set_ncols(p%desc_ac%get_local_cols())
if (np>1) then
call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_)
call op_restr%mv_to(acoo)
nzl = acoo%get_nzeros()
if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),p%desc_ac,info,'I')
call acoo%set_dupl(psb_dupl_add_)
if (info == psb_success_) call op_restr%mv_from(acoo)
if (info == psb_success_) call op_restr%cscnv(info,type='csr')
if(info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,&
& a_err='Converting op_restr to local')
goto 9999
end if
end if
!
! Clip to local rows.
!
call op_restr%set_nrows(p%desc_ac%get_local_rows())
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& 'Done ac '
case(mld_repl_mat_)
!
!
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if ((info == psb_success_).and.p%parms%clean_zeros) call ac%clean_zeros(info)
if (info == psb_success_) &
& call psb_gather(p%ac,ac,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.)
if (info /= psb_success_) goto 9999
case default
info = psb_err_internal_error_
call psb_errpush(info,name,a_err='invalid mld_coarse_mat_')
goto 9999
end select
call p%ac%cscnv(info,type='csr',dupl=psb_dupl_add_)
if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='spcnv')
goto 9999
end if
!
! Copy the prolongation/restriction matrices into the descriptor map.
! op_restr => PR^T i.e. restriction operator
! op_prol => PR i.e. prolongation operator
!
p%map = psb_linmap(psb_map_aggr_,desc_a,&
& p%desc_ac,op_restr,op_prol,ilaggr,nlaggr)
if (info == psb_success_) call op_prol%free()
if (info == psb_success_) call op_restr%free()
if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='sp_Free')
goto 9999
end if
! !
! Fix the base_a and base_desc pointers for handling of residuals. ! Fix the base_a and base_desc pointers for handling of residuals.
! This is correct because this routine is only called at levels >=2. ! This is correct because this routine is only called at levels >=2.

26
mlprec/mld_d_inner_mod.f90
Unescape Escape View File

@ -165,18 +165,18 @@ module mld_d_inner_mod
end interface mld_dec_map_bld end interface mld_dec_map_bld
interface mld_aggrmat_asb !!$ interface mld_aggrmat_asb
subroutine mld_daggrmat_asb(a,desc_a,ilaggr,nlaggr,p,info) !!$ subroutine mld_daggrmat_asb(a,desc_a,ilaggr,nlaggr,p,info)
use psb_base_mod, only : psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_ !!$ use psb_base_mod, only : psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_
use mld_d_prec_type, only : mld_d_onelev_type !!$ use mld_d_prec_type, only : mld_d_onelev_type
implicit none !!$ implicit none
type(psb_dspmat_type), intent(in) :: a !!$ type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a !!$ type(psb_desc_type), intent(in) :: desc_a
integer(psb_ipk_), intent(inout) :: ilaggr(:), nlaggr(:) !!$ integer(psb_ipk_), intent(inout) :: ilaggr(:), nlaggr(:)
type(mld_d_onelev_type), intent(inout), target :: p !!$ type(mld_d_onelev_type), intent(inout), target :: p
integer(psb_ipk_), intent(out) :: info !!$ integer(psb_ipk_), intent(out) :: info
end subroutine mld_daggrmat_asb !!$ end subroutine mld_daggrmat_asb
end interface mld_aggrmat_asb !!$ end interface mld_aggrmat_asb
@ -197,7 +197,7 @@ module mld_d_inner_mod
procedure(mld_daggrmat_var_asb) :: mld_daggrmat_nosmth_asb, & procedure(mld_daggrmat_var_asb) :: mld_daggrmat_nosmth_asb, &
& mld_daggrmat_smth_asb, mld_daggrmat_minnrg_asb, & & mld_daggrmat_smth_asb, mld_daggrmat_minnrg_asb, &
& mld_daggrmat_biz_asb & mld_daggrmat_biz_asb, mld_daggrmat_asb
end module mld_d_inner_mod end module mld_d_inner_mod

Write Preview
Loading…
Cancel
Save