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

Changes to aggregation with new aggregation object. Fix CSETC.

Browse Source
stopcriterion
Salvatore Filippone 6 years ago
parent 5a1dace4e3
commit f3edc6a467
69 changed files with 715 additions and 593 deletions
Show Stats Download Patch File Download Diff File

19
mlprec/impl/aggregator/mld_c_dec_aggregator_tprol.f90
Unescape Escape View File

@ -40,7 +40,7 @@
! Subroutine: mld_c_dec_aggregator_tprol
! Version: complex
!
! This routine is mainly an interface to map_bld where the real work is performed.
! This routine is mainly an interface to soc_map_bld where the real work is performed.
! It takes care of some consistency checking, and calls map_to_tprol, which is
! refactored and shared among all the aggregation methods that produce a simple
! integer mapping.
@ -50,6 +50,8 @@
! ag - type(mld_c_dec_aggregator_type), input/output.
! The aggregator object, carrying with itself the mapping algorithm.
! parms - The auxiliary parameters object
! ag_data - Auxiliary global aggregation parameters object
!
!
! a - type(psb_cspmat_type).
! The sparse matrix structure containing the local part of the
@ -72,13 +74,15 @@
! info - integer, output.
! Error code.
!
subroutine mld_c_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_c_dec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
use mld_c_prec_type, mld_protect_name => mld_c_dec_aggregator_build_tprol
use mld_c_inner_mod
implicit none
class(mld_c_dec_aggregator_type), target, intent(inout) :: ag
type(mld_sml_parms), intent(inout) :: parms
type(mld_saggr_data), intent(in) :: ag_data
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -91,6 +95,7 @@ subroutine mld_c_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_
integer(psb_ipk_) :: err_act
integer(psb_lpk_) :: ntaggr
integer(psb_ipk_) :: debug_level, debug_unit
logical :: clean_zeros
name='mld_c_dec_aggregator_tprol'
call psb_erractionsave(err_act)
@ -111,13 +116,17 @@ subroutine mld_c_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_
& mld_aggr_ord_nat_,is_legal_ml_aggr_ord)
call mld_check_def(parms%aggr_thresh,'Aggr_Thresh',szero,is_legal_s_aggr_thrs)
call ag%map_bld(parms%aggr_ord,parms%aggr_thresh,a,desc_a,nlaggr,ilaggr,info)
!
! The decoupled aggregator based on SOC measures ignores
! ag_data except for clean_zeros; soc_map_bld is a procedure pointer.
!
clean_zeros = ag%do_clean_zeros
call ag%soc_map_bld(parms%aggr_ord,parms%aggr_thresh,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
if (info==psb_success_) call mld_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='map_bld/map_to_tprol')
call psb_errpush(info,name,a_err='soc_map_bld/map_to_tprol')
goto 9999
endif

30
mlprec/impl/aggregator/mld_c_soc1_map_bld.f90
Unescape Escape View File

@ -67,7 +67,7 @@
!
!
!
subroutine mld_c_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_c_soc1_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_base_prec_type
@ -77,6 +77,7 @@ subroutine mld_c_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Arguments
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_spk_), intent(in) :: theta
@ -109,8 +110,8 @@ subroutine mld_c_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
!
ictxt=desc_a%get_context()
call psb_info(ictxt,me,np)
nrow = desc_a%get_local_rows()
ncol = desc_a%get_local_cols()
nrow = desc_a%get_local_rows()
ncol = desc_a%get_local_cols()
nrglob = desc_a%get_global_rows()
nr = a%get_nrows()
@ -132,6 +133,7 @@ subroutine mld_c_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end if
call a%cp_to(acsr)
if (clean_zeros) call acsr%clean_zeros(info)
if (iorder == mld_aggr_ord_nat_) then
do i=1, nr
ilaggr(i) = -(nr+1)
@ -169,12 +171,6 @@ subroutine mld_c_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!!$ call a%csget(i,i,nz,irow,icol,val,info,chksz=.false.)
!!$ if (info /= psb_success_) then
!!$ info=psb_err_from_subroutine_
!!$ call psb_errpush(info,name,a_err='csget')
!!$ goto 9999
!!$ end if
!
! Build the set of all strongly coupled nodes
@ -194,7 +190,7 @@ subroutine mld_c_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! If the whole strongly coupled neighborhood of I is
! as yet unconnected, turn it into the next aggregate.
! Same if ip==0 (in which case, neighborhood only
! contains I even if it does not look from matrix)
! contains I even if it does not look like it from matrix)
!
disjoint = all(ilaggr(icol(1:ip)) == -(nr+1)).or.(ip==0)
if (disjoint) then
@ -222,14 +218,10 @@ subroutine mld_c_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
if (ilaggr(i) == -(nr+1)) then
nz = (acsr%irp(i+1)-acsr%irp(i))
if (nz == 1) cycle step2
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!!$ call a%csget(i,i,nz,irow,icol,val,info,chksz=.false.)
!!$ if (info /= psb_success_) then
!!$ info=psb_err_from_subroutine_
!!$ call psb_errpush(info,name,a_err='psb_sp_getrow')
!!$ goto 9999
!!$ end if
!
! Find the most strongly connected neighbour that is
! already aggregated, if any, and join its aggregate
@ -261,6 +253,7 @@ subroutine mld_c_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
if (ilaggr(i) < 0) then
nz = (acsr%irp(i+1)-acsr%irp(i))
if (nz == 1) cycle step3
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!
@ -288,8 +281,8 @@ subroutine mld_c_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end do
else
!
! This should not happen: we did not even connect with ourselves.
! Create an isolate anyway.
! This should not happen: we did not even connect with ourselves,
! but it's not a singleton.
!
naggr = naggr + 1
ilaggr(i) = naggr
@ -297,7 +290,6 @@ subroutine mld_c_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end if
end do step3
! Any leftovers?
do i=1, nr
if (ilaggr(i) < 0) then

6
mlprec/impl/aggregator/mld_c_soc2_map_bld.f90
Unescape Escape View File

@ -66,7 +66,7 @@
!
!
!
subroutine mld_c_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_c_soc2_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_base_prec_type
@ -76,6 +76,7 @@ subroutine mld_c_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Arguments
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_spk_), intent(in) :: theta
@ -134,13 +135,14 @@ subroutine mld_c_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Phase zero: compute muij
!
call a%cp_to(muij)
if (clean_zeros) call muij%clean_zeros(info)
do i=1, nr
do k=muij%irp(i),muij%irp(i+1)-1
j = muij%ja(k)
if (j<= nr) muij%val(k) = abs(muij%val(k))/sqrt(abs(diag(i)*diag(j)))
end do
end do
!write(*,*) 'murows/cols ',maxval(mu_rows),maxval(mu_cols)
!
! Compute the 1-neigbour; mark strong links with +1, weak links with -1
!

18
mlprec/impl/aggregator/mld_c_symdec_aggregator_tprol.f90
Unescape Escape View File

@ -41,7 +41,7 @@
! Version: complex
!
!
! This routine is mainly an interface to map_bld where the real work is performed.
! This routine is mainly an interface to soc_map_bld where the real work is performed.
! It takes care of some consistency checking, and calls map_to_tprol, which is
! refactored and shared among all the aggregation methods that produce a simple
! integer mapping. It also symmetrizes the pattern of the local matrix A.
@ -53,6 +53,7 @@
! ag - type(mld_c_dec_aggregator_type), input/output.
! The aggregator object, carrying with itself the mapping algorithm.
! parms - The auxiliary parameters object
! ag_data - Auxiliary global aggregation parameters object
!
! a - type(psb_cspmat_type).
! The sparse matrix structure containing the local part of the
@ -75,7 +76,8 @@
! info - integer, output.
! Error code.
!
subroutine mld_c_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_c_symdec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
use mld_c_prec_type
use mld_c_symdec_aggregator_mod, mld_protect_name => mld_c_symdec_aggregator_build_tprol
@ -83,6 +85,7 @@ subroutine mld_c_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
implicit none
class(mld_c_symdec_aggregator_type), target, intent(inout) :: ag
type(mld_sml_parms), intent(inout) :: parms
type(mld_saggr_data), intent(in) :: ag_data
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -97,6 +100,7 @@ subroutine mld_c_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
integer(psb_ipk_) :: nr
integer(psb_lpk_) :: ntaggr
integer(psb_ipk_) :: debug_level, debug_unit
logical :: clean_zeros
name='mld_c_symdec_aggregator_tprol'
call psb_erractionsave(err_act)
@ -130,14 +134,20 @@ subroutine mld_c_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
if (info == psb_success_) call atrans%free()
if (info == psb_success_) call atmp%cscnv(info,type='CSR')
!
! The decoupled aggregator based on SOC measures ignores
! ag_data except for clean_zeros; soc_map_bld is a procedure pointer.
!
clean_zeros = ag%do_clean_zeros
if (info == psb_success_) &
& call ag%map_bld(parms%aggr_ord,parms%aggr_thresh,atmp,desc_a,nlaggr,ilaggr,info)
& call ag%soc_map_bld(parms%aggr_ord,parms%aggr_thresh,clean_zeros,atmp,&
& desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call atmp%free()
if (info == psb_success_) call mld_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='map_bld/map_to_tprol')
call psb_errpush(info,name,a_err='soc_map_bld/map_to_tprol')
goto 9999
endif

19
mlprec/impl/aggregator/mld_d_dec_aggregator_tprol.f90
Unescape Escape View File

@ -40,7 +40,7 @@
! Subroutine: mld_d_dec_aggregator_tprol
! Version: real
!
! This routine is mainly an interface to map_bld where the real work is performed.
! This routine is mainly an interface to soc_map_bld where the real work is performed.
! It takes care of some consistency checking, and calls map_to_tprol, which is
! refactored and shared among all the aggregation methods that produce a simple
! integer mapping.
@ -50,6 +50,8 @@
! ag - type(mld_d_dec_aggregator_type), input/output.
! The aggregator object, carrying with itself the mapping algorithm.
! parms - The auxiliary parameters object
! ag_data - Auxiliary global aggregation parameters object
!
!
! a - type(psb_dspmat_type).
! The sparse matrix structure containing the local part of the
@ -72,13 +74,15 @@
! info - integer, output.
! Error code.
!
subroutine mld_d_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_d_dec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
use mld_d_prec_type, mld_protect_name => mld_d_dec_aggregator_build_tprol
use mld_d_inner_mod
implicit none
class(mld_d_dec_aggregator_type), target, intent(inout) :: ag
type(mld_dml_parms), intent(inout) :: parms
type(mld_daggr_data), intent(in) :: ag_data
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -91,6 +95,7 @@ subroutine mld_d_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_
integer(psb_ipk_) :: err_act
integer(psb_lpk_) :: ntaggr
integer(psb_ipk_) :: debug_level, debug_unit
logical :: clean_zeros
name='mld_d_dec_aggregator_tprol'
call psb_erractionsave(err_act)
@ -111,13 +116,17 @@ subroutine mld_d_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_
& mld_aggr_ord_nat_,is_legal_ml_aggr_ord)
call mld_check_def(parms%aggr_thresh,'Aggr_Thresh',dzero,is_legal_d_aggr_thrs)
call ag%map_bld(parms%aggr_ord,parms%aggr_thresh,a,desc_a,nlaggr,ilaggr,info)
!
! The decoupled aggregator based on SOC measures ignores
! ag_data except for clean_zeros; soc_map_bld is a procedure pointer.
!
clean_zeros = ag%do_clean_zeros
call ag%soc_map_bld(parms%aggr_ord,parms%aggr_thresh,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
if (info==psb_success_) call mld_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='map_bld/map_to_tprol')
call psb_errpush(info,name,a_err='soc_map_bld/map_to_tprol')
goto 9999
endif

30
mlprec/impl/aggregator/mld_d_soc1_map_bld.f90
Unescape Escape View File

@ -67,7 +67,7 @@
!
!
!
subroutine mld_d_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_d_soc1_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_base_prec_type
@ -77,6 +77,7 @@ subroutine mld_d_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Arguments
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_dpk_), intent(in) :: theta
@ -109,8 +110,8 @@ subroutine mld_d_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
!
ictxt=desc_a%get_context()
call psb_info(ictxt,me,np)
nrow = desc_a%get_local_rows()
ncol = desc_a%get_local_cols()
nrow = desc_a%get_local_rows()
ncol = desc_a%get_local_cols()
nrglob = desc_a%get_global_rows()
nr = a%get_nrows()
@ -132,6 +133,7 @@ subroutine mld_d_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end if
call a%cp_to(acsr)
if (clean_zeros) call acsr%clean_zeros(info)
if (iorder == mld_aggr_ord_nat_) then
do i=1, nr
ilaggr(i) = -(nr+1)
@ -169,12 +171,6 @@ subroutine mld_d_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!!$ call a%csget(i,i,nz,irow,icol,val,info,chksz=.false.)
!!$ if (info /= psb_success_) then
!!$ info=psb_err_from_subroutine_
!!$ call psb_errpush(info,name,a_err='csget')
!!$ goto 9999
!!$ end if
!
! Build the set of all strongly coupled nodes
@ -194,7 +190,7 @@ subroutine mld_d_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! If the whole strongly coupled neighborhood of I is
! as yet unconnected, turn it into the next aggregate.
! Same if ip==0 (in which case, neighborhood only
! contains I even if it does not look from matrix)
! contains I even if it does not look like it from matrix)
!
disjoint = all(ilaggr(icol(1:ip)) == -(nr+1)).or.(ip==0)
if (disjoint) then
@ -222,14 +218,10 @@ subroutine mld_d_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
if (ilaggr(i) == -(nr+1)) then
nz = (acsr%irp(i+1)-acsr%irp(i))
if (nz == 1) cycle step2
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!!$ call a%csget(i,i,nz,irow,icol,val,info,chksz=.false.)
!!$ if (info /= psb_success_) then
!!$ info=psb_err_from_subroutine_
!!$ call psb_errpush(info,name,a_err='psb_sp_getrow')
!!$ goto 9999
!!$ end if
!
! Find the most strongly connected neighbour that is
! already aggregated, if any, and join its aggregate
@ -261,6 +253,7 @@ subroutine mld_d_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
if (ilaggr(i) < 0) then
nz = (acsr%irp(i+1)-acsr%irp(i))
if (nz == 1) cycle step3
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!
@ -288,8 +281,8 @@ subroutine mld_d_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end do
else
!
! This should not happen: we did not even connect with ourselves.
! Create an isolate anyway.
! This should not happen: we did not even connect with ourselves,
! but it's not a singleton.
!
naggr = naggr + 1
ilaggr(i) = naggr
@ -297,7 +290,6 @@ subroutine mld_d_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end if
end do step3
! Any leftovers?
do i=1, nr
if (ilaggr(i) < 0) then

6
mlprec/impl/aggregator/mld_d_soc2_map_bld.f90
Unescape Escape View File

@ -66,7 +66,7 @@
!
!
!
subroutine mld_d_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_d_soc2_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_base_prec_type
@ -76,6 +76,7 @@ subroutine mld_d_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Arguments
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_dpk_), intent(in) :: theta
@ -134,13 +135,14 @@ subroutine mld_d_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Phase zero: compute muij
!
call a%cp_to(muij)
if (clean_zeros) call muij%clean_zeros(info)
do i=1, nr
do k=muij%irp(i),muij%irp(i+1)-1
j = muij%ja(k)
if (j<= nr) muij%val(k) = abs(muij%val(k))/sqrt(abs(diag(i)*diag(j)))
end do
end do
!write(*,*) 'murows/cols ',maxval(mu_rows),maxval(mu_cols)
!
! Compute the 1-neigbour; mark strong links with +1, weak links with -1
!

18
mlprec/impl/aggregator/mld_d_symdec_aggregator_tprol.f90
Unescape Escape View File

@ -41,7 +41,7 @@
! Version: real
!
!
! This routine is mainly an interface to map_bld where the real work is performed.
! This routine is mainly an interface to soc_map_bld where the real work is performed.
! It takes care of some consistency checking, and calls map_to_tprol, which is
! refactored and shared among all the aggregation methods that produce a simple
! integer mapping. It also symmetrizes the pattern of the local matrix A.
@ -53,6 +53,7 @@
! ag - type(mld_d_dec_aggregator_type), input/output.
! The aggregator object, carrying with itself the mapping algorithm.
! parms - The auxiliary parameters object
! ag_data - Auxiliary global aggregation parameters object
!
! a - type(psb_dspmat_type).
! The sparse matrix structure containing the local part of the
@ -75,7 +76,8 @@
! info - integer, output.
! Error code.
!
subroutine mld_d_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_d_symdec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
use mld_d_prec_type
use mld_d_symdec_aggregator_mod, mld_protect_name => mld_d_symdec_aggregator_build_tprol
@ -83,6 +85,7 @@ subroutine mld_d_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
implicit none
class(mld_d_symdec_aggregator_type), target, intent(inout) :: ag
type(mld_dml_parms), intent(inout) :: parms
type(mld_daggr_data), intent(in) :: ag_data
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -97,6 +100,7 @@ subroutine mld_d_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
integer(psb_ipk_) :: nr
integer(psb_lpk_) :: ntaggr
integer(psb_ipk_) :: debug_level, debug_unit
logical :: clean_zeros
name='mld_d_symdec_aggregator_tprol'
call psb_erractionsave(err_act)
@ -130,14 +134,20 @@ subroutine mld_d_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
if (info == psb_success_) call atrans%free()
if (info == psb_success_) call atmp%cscnv(info,type='CSR')
!
! The decoupled aggregator based on SOC measures ignores
! ag_data except for clean_zeros; soc_map_bld is a procedure pointer.
!
clean_zeros = ag%do_clean_zeros
if (info == psb_success_) &
& call ag%map_bld(parms%aggr_ord,parms%aggr_thresh,atmp,desc_a,nlaggr,ilaggr,info)
& call ag%soc_map_bld(parms%aggr_ord,parms%aggr_thresh,clean_zeros,atmp,&
& desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call atmp%free()
if (info == psb_success_) call mld_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='map_bld/map_to_tprol')
call psb_errpush(info,name,a_err='soc_map_bld/map_to_tprol')
goto 9999
endif

19
mlprec/impl/aggregator/mld_s_dec_aggregator_tprol.f90
Unescape Escape View File

@ -40,7 +40,7 @@
! Subroutine: mld_s_dec_aggregator_tprol
! Version: real
!
! This routine is mainly an interface to map_bld where the real work is performed.
! This routine is mainly an interface to soc_map_bld where the real work is performed.
! It takes care of some consistency checking, and calls map_to_tprol, which is
! refactored and shared among all the aggregation methods that produce a simple
! integer mapping.
@ -50,6 +50,8 @@
! ag - type(mld_s_dec_aggregator_type), input/output.
! The aggregator object, carrying with itself the mapping algorithm.
! parms - The auxiliary parameters object
! ag_data - Auxiliary global aggregation parameters object
!
!
! a - type(psb_sspmat_type).
! The sparse matrix structure containing the local part of the
@ -72,13 +74,15 @@
! info - integer, output.
! Error code.
!
subroutine mld_s_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_s_dec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
use mld_s_prec_type, mld_protect_name => mld_s_dec_aggregator_build_tprol
use mld_s_inner_mod
implicit none
class(mld_s_dec_aggregator_type), target, intent(inout) :: ag
type(mld_sml_parms), intent(inout) :: parms
type(mld_saggr_data), intent(in) :: ag_data
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -91,6 +95,7 @@ subroutine mld_s_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_
integer(psb_ipk_) :: err_act
integer(psb_lpk_) :: ntaggr
integer(psb_ipk_) :: debug_level, debug_unit
logical :: clean_zeros
name='mld_s_dec_aggregator_tprol'
call psb_erractionsave(err_act)
@ -111,13 +116,17 @@ subroutine mld_s_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_
& mld_aggr_ord_nat_,is_legal_ml_aggr_ord)
call mld_check_def(parms%aggr_thresh,'Aggr_Thresh',szero,is_legal_s_aggr_thrs)
call ag%map_bld(parms%aggr_ord,parms%aggr_thresh,a,desc_a,nlaggr,ilaggr,info)
!
! The decoupled aggregator based on SOC measures ignores
! ag_data except for clean_zeros; soc_map_bld is a procedure pointer.
!
clean_zeros = ag%do_clean_zeros
call ag%soc_map_bld(parms%aggr_ord,parms%aggr_thresh,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
if (info==psb_success_) call mld_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='map_bld/map_to_tprol')
call psb_errpush(info,name,a_err='soc_map_bld/map_to_tprol')
goto 9999
endif

30
mlprec/impl/aggregator/mld_s_soc1_map_bld.f90
Unescape Escape View File

@ -67,7 +67,7 @@
!
!
!
subroutine mld_s_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_s_soc1_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_base_prec_type
@ -77,6 +77,7 @@ subroutine mld_s_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Arguments
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_spk_), intent(in) :: theta
@ -109,8 +110,8 @@ subroutine mld_s_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
!
ictxt=desc_a%get_context()
call psb_info(ictxt,me,np)
nrow = desc_a%get_local_rows()
ncol = desc_a%get_local_cols()
nrow = desc_a%get_local_rows()
ncol = desc_a%get_local_cols()
nrglob = desc_a%get_global_rows()
nr = a%get_nrows()
@ -132,6 +133,7 @@ subroutine mld_s_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end if
call a%cp_to(acsr)
if (clean_zeros) call acsr%clean_zeros(info)
if (iorder == mld_aggr_ord_nat_) then
do i=1, nr
ilaggr(i) = -(nr+1)
@ -169,12 +171,6 @@ subroutine mld_s_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!!$ call a%csget(i,i,nz,irow,icol,val,info,chksz=.false.)
!!$ if (info /= psb_success_) then
!!$ info=psb_err_from_subroutine_
!!$ call psb_errpush(info,name,a_err='csget')
!!$ goto 9999
!!$ end if
!
! Build the set of all strongly coupled nodes
@ -194,7 +190,7 @@ subroutine mld_s_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! If the whole strongly coupled neighborhood of I is
! as yet unconnected, turn it into the next aggregate.
! Same if ip==0 (in which case, neighborhood only
! contains I even if it does not look from matrix)
! contains I even if it does not look like it from matrix)
!
disjoint = all(ilaggr(icol(1:ip)) == -(nr+1)).or.(ip==0)
if (disjoint) then
@ -222,14 +218,10 @@ subroutine mld_s_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
if (ilaggr(i) == -(nr+1)) then
nz = (acsr%irp(i+1)-acsr%irp(i))
if (nz == 1) cycle step2
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!!$ call a%csget(i,i,nz,irow,icol,val,info,chksz=.false.)
!!$ if (info /= psb_success_) then
!!$ info=psb_err_from_subroutine_
!!$ call psb_errpush(info,name,a_err='psb_sp_getrow')
!!$ goto 9999
!!$ end if
!
! Find the most strongly connected neighbour that is
! already aggregated, if any, and join its aggregate
@ -261,6 +253,7 @@ subroutine mld_s_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
if (ilaggr(i) < 0) then
nz = (acsr%irp(i+1)-acsr%irp(i))
if (nz == 1) cycle step3
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!
@ -288,8 +281,8 @@ subroutine mld_s_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end do
else
!
! This should not happen: we did not even connect with ourselves.
! Create an isolate anyway.
! This should not happen: we did not even connect with ourselves,
! but it's not a singleton.
!
naggr = naggr + 1
ilaggr(i) = naggr
@ -297,7 +290,6 @@ subroutine mld_s_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end if
end do step3
! Any leftovers?
do i=1, nr
if (ilaggr(i) < 0) then

6
mlprec/impl/aggregator/mld_s_soc2_map_bld.f90
Unescape Escape View File

@ -66,7 +66,7 @@
!
!
!
subroutine mld_s_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_s_soc2_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_base_prec_type
@ -76,6 +76,7 @@ subroutine mld_s_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Arguments
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_spk_), intent(in) :: theta
@ -134,13 +135,14 @@ subroutine mld_s_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Phase zero: compute muij
!
call a%cp_to(muij)
if (clean_zeros) call muij%clean_zeros(info)
do i=1, nr
do k=muij%irp(i),muij%irp(i+1)-1
j = muij%ja(k)
if (j<= nr) muij%val(k) = abs(muij%val(k))/sqrt(abs(diag(i)*diag(j)))
end do
end do
!write(*,*) 'murows/cols ',maxval(mu_rows),maxval(mu_cols)
!
! Compute the 1-neigbour; mark strong links with +1, weak links with -1
!

18
mlprec/impl/aggregator/mld_s_symdec_aggregator_tprol.f90
Unescape Escape View File

@ -41,7 +41,7 @@
! Version: real
!
!
! This routine is mainly an interface to map_bld where the real work is performed.
! This routine is mainly an interface to soc_map_bld where the real work is performed.
! It takes care of some consistency checking, and calls map_to_tprol, which is
! refactored and shared among all the aggregation methods that produce a simple
! integer mapping. It also symmetrizes the pattern of the local matrix A.
@ -53,6 +53,7 @@
! ag - type(mld_s_dec_aggregator_type), input/output.
! The aggregator object, carrying with itself the mapping algorithm.
! parms - The auxiliary parameters object
! ag_data - Auxiliary global aggregation parameters object
!
! a - type(psb_sspmat_type).
! The sparse matrix structure containing the local part of the
@ -75,7 +76,8 @@
! info - integer, output.
! Error code.
!
subroutine mld_s_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_s_symdec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
use mld_s_prec_type
use mld_s_symdec_aggregator_mod, mld_protect_name => mld_s_symdec_aggregator_build_tprol
@ -83,6 +85,7 @@ subroutine mld_s_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
implicit none
class(mld_s_symdec_aggregator_type), target, intent(inout) :: ag
type(mld_sml_parms), intent(inout) :: parms
type(mld_saggr_data), intent(in) :: ag_data
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -97,6 +100,7 @@ subroutine mld_s_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
integer(psb_ipk_) :: nr
integer(psb_lpk_) :: ntaggr
integer(psb_ipk_) :: debug_level, debug_unit
logical :: clean_zeros
name='mld_s_symdec_aggregator_tprol'
call psb_erractionsave(err_act)
@ -130,14 +134,20 @@ subroutine mld_s_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
if (info == psb_success_) call atrans%free()
if (info == psb_success_) call atmp%cscnv(info,type='CSR')
!
! The decoupled aggregator based on SOC measures ignores
! ag_data except for clean_zeros; soc_map_bld is a procedure pointer.
!
clean_zeros = ag%do_clean_zeros
if (info == psb_success_) &
& call ag%map_bld(parms%aggr_ord,parms%aggr_thresh,atmp,desc_a,nlaggr,ilaggr,info)
& call ag%soc_map_bld(parms%aggr_ord,parms%aggr_thresh,clean_zeros,atmp,&
& desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call atmp%free()
if (info == psb_success_) call mld_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='map_bld/map_to_tprol')
call psb_errpush(info,name,a_err='soc_map_bld/map_to_tprol')
goto 9999
endif

19
mlprec/impl/aggregator/mld_z_dec_aggregator_tprol.f90
Unescape Escape View File

@ -40,7 +40,7 @@
! Subroutine: mld_z_dec_aggregator_tprol
! Version: complex
!
! This routine is mainly an interface to map_bld where the real work is performed.
! This routine is mainly an interface to soc_map_bld where the real work is performed.
! It takes care of some consistency checking, and calls map_to_tprol, which is
! refactored and shared among all the aggregation methods that produce a simple
! integer mapping.
@ -50,6 +50,8 @@
! ag - type(mld_z_dec_aggregator_type), input/output.
! The aggregator object, carrying with itself the mapping algorithm.
! parms - The auxiliary parameters object
! ag_data - Auxiliary global aggregation parameters object
!
!
! a - type(psb_zspmat_type).
! The sparse matrix structure containing the local part of the
@ -72,13 +74,15 @@
! info - integer, output.
! Error code.
!
subroutine mld_z_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_z_dec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
use mld_z_prec_type, mld_protect_name => mld_z_dec_aggregator_build_tprol
use mld_z_inner_mod
implicit none
class(mld_z_dec_aggregator_type), target, intent(inout) :: ag
type(mld_dml_parms), intent(inout) :: parms
type(mld_daggr_data), intent(in) :: ag_data
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -91,6 +95,7 @@ subroutine mld_z_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_
integer(psb_ipk_) :: err_act
integer(psb_lpk_) :: ntaggr
integer(psb_ipk_) :: debug_level, debug_unit
logical :: clean_zeros
name='mld_z_dec_aggregator_tprol'
call psb_erractionsave(err_act)
@ -111,13 +116,17 @@ subroutine mld_z_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_
& mld_aggr_ord_nat_,is_legal_ml_aggr_ord)
call mld_check_def(parms%aggr_thresh,'Aggr_Thresh',dzero,is_legal_d_aggr_thrs)
call ag%map_bld(parms%aggr_ord,parms%aggr_thresh,a,desc_a,nlaggr,ilaggr,info)
!
! The decoupled aggregator based on SOC measures ignores
! ag_data except for clean_zeros; soc_map_bld is a procedure pointer.
!
clean_zeros = ag%do_clean_zeros
call ag%soc_map_bld(parms%aggr_ord,parms%aggr_thresh,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
if (info==psb_success_) call mld_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='map_bld/map_to_tprol')
call psb_errpush(info,name,a_err='soc_map_bld/map_to_tprol')
goto 9999
endif

30
mlprec/impl/aggregator/mld_z_soc1_map_bld.f90
Unescape Escape View File

@ -67,7 +67,7 @@
!
!
!
subroutine mld_z_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_z_soc1_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_base_prec_type
@ -77,6 +77,7 @@ subroutine mld_z_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Arguments
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_dpk_), intent(in) :: theta
@ -109,8 +110,8 @@ subroutine mld_z_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
!
ictxt=desc_a%get_context()
call psb_info(ictxt,me,np)
nrow = desc_a%get_local_rows()
ncol = desc_a%get_local_cols()
nrow = desc_a%get_local_rows()
ncol = desc_a%get_local_cols()
nrglob = desc_a%get_global_rows()
nr = a%get_nrows()
@ -132,6 +133,7 @@ subroutine mld_z_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end if
call a%cp_to(acsr)
if (clean_zeros) call acsr%clean_zeros(info)
if (iorder == mld_aggr_ord_nat_) then
do i=1, nr
ilaggr(i) = -(nr+1)
@ -169,12 +171,6 @@ subroutine mld_z_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!!$ call a%csget(i,i,nz,irow,icol,val,info,chksz=.false.)
!!$ if (info /= psb_success_) then
!!$ info=psb_err_from_subroutine_
!!$ call psb_errpush(info,name,a_err='csget')
!!$ goto 9999
!!$ end if
!
! Build the set of all strongly coupled nodes
@ -194,7 +190,7 @@ subroutine mld_z_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! If the whole strongly coupled neighborhood of I is
! as yet unconnected, turn it into the next aggregate.
! Same if ip==0 (in which case, neighborhood only
! contains I even if it does not look from matrix)
! contains I even if it does not look like it from matrix)
!
disjoint = all(ilaggr(icol(1:ip)) == -(nr+1)).or.(ip==0)
if (disjoint) then
@ -222,14 +218,10 @@ subroutine mld_z_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
if (ilaggr(i) == -(nr+1)) then
nz = (acsr%irp(i+1)-acsr%irp(i))
if (nz == 1) cycle step2
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!!$ call a%csget(i,i,nz,irow,icol,val,info,chksz=.false.)
!!$ if (info /= psb_success_) then
!!$ info=psb_err_from_subroutine_
!!$ call psb_errpush(info,name,a_err='psb_sp_getrow')
!!$ goto 9999
!!$ end if
!
! Find the most strongly connected neighbour that is
! already aggregated, if any, and join its aggregate
@ -261,6 +253,7 @@ subroutine mld_z_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
if (ilaggr(i) < 0) then
nz = (acsr%irp(i+1)-acsr%irp(i))
if (nz == 1) cycle step3
icol(1:nz) = acsr%ja(acsr%irp(i):acsr%irp(i+1)-1)
val(1:nz) = acsr%val(acsr%irp(i):acsr%irp(i+1)-1)
!
@ -288,8 +281,8 @@ subroutine mld_z_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end do
else
!
! This should not happen: we did not even connect with ourselves.
! Create an isolate anyway.
! This should not happen: we did not even connect with ourselves,
! but it's not a singleton.
!
naggr = naggr + 1
ilaggr(i) = naggr
@ -297,7 +290,6 @@ subroutine mld_z_soc1_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
end if
end do step3
! Any leftovers?
do i=1, nr
if (ilaggr(i) < 0) then

6
mlprec/impl/aggregator/mld_z_soc2_map_bld.f90
Unescape Escape View File

@ -66,7 +66,7 @@
!
!
!
subroutine mld_z_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_z_soc2_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_base_prec_type
@ -76,6 +76,7 @@ subroutine mld_z_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Arguments
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_dpk_), intent(in) :: theta
@ -134,13 +135,14 @@ subroutine mld_z_soc2_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
! Phase zero: compute muij
!
call a%cp_to(muij)
if (clean_zeros) call muij%clean_zeros(info)
do i=1, nr
do k=muij%irp(i),muij%irp(i+1)-1
j = muij%ja(k)
if (j<= nr) muij%val(k) = abs(muij%val(k))/sqrt(abs(diag(i)*diag(j)))
end do
end do
!write(*,*) 'murows/cols ',maxval(mu_rows),maxval(mu_cols)
!
! Compute the 1-neigbour; mark strong links with +1, weak links with -1
!

18
mlprec/impl/aggregator/mld_z_symdec_aggregator_tprol.f90
Unescape Escape View File

@ -41,7 +41,7 @@
! Version: complex
!
!
! This routine is mainly an interface to map_bld where the real work is performed.
! This routine is mainly an interface to soc_map_bld where the real work is performed.
! It takes care of some consistency checking, and calls map_to_tprol, which is
! refactored and shared among all the aggregation methods that produce a simple
! integer mapping. It also symmetrizes the pattern of the local matrix A.
@ -53,6 +53,7 @@
! ag - type(mld_z_dec_aggregator_type), input/output.
! The aggregator object, carrying with itself the mapping algorithm.
! parms - The auxiliary parameters object
! ag_data - Auxiliary global aggregation parameters object
!
! a - type(psb_zspmat_type).
! The sparse matrix structure containing the local part of the
@ -75,7 +76,8 @@
! info - integer, output.
! Error code.
!
subroutine mld_z_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_z_symdec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
use mld_z_prec_type
use mld_z_symdec_aggregator_mod, mld_protect_name => mld_z_symdec_aggregator_build_tprol
@ -83,6 +85,7 @@ subroutine mld_z_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
implicit none
class(mld_z_symdec_aggregator_type), target, intent(inout) :: ag
type(mld_dml_parms), intent(inout) :: parms
type(mld_daggr_data), intent(in) :: ag_data
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -97,6 +100,7 @@ subroutine mld_z_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
integer(psb_ipk_) :: nr
integer(psb_lpk_) :: ntaggr
integer(psb_ipk_) :: debug_level, debug_unit
logical :: clean_zeros
name='mld_z_symdec_aggregator_tprol'
call psb_erractionsave(err_act)
@ -130,14 +134,20 @@ subroutine mld_z_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,
if (info == psb_success_) call atrans%free()
if (info == psb_success_) call atmp%cscnv(info,type='CSR')
!
! The decoupled aggregator based on SOC measures ignores
! ag_data except for clean_zeros; soc_map_bld is a procedure pointer.
!
clean_zeros = ag%do_clean_zeros
if (info == psb_success_) &
& call ag%map_bld(parms%aggr_ord,parms%aggr_thresh,atmp,desc_a,nlaggr,ilaggr,info)
& call ag%soc_map_bld(parms%aggr_ord,parms%aggr_thresh,clean_zeros,atmp,&
& desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call atmp%free()
if (info == psb_success_) call mld_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='map_bld/map_to_tprol')
call psb_errpush(info,name,a_err='soc_map_bld/map_to_tprol')
goto 9999
endif

2
mlprec/impl/level/mld_c_base_onelev_cnv.f90
Unescape Escape View File

@ -62,6 +62,6 @@ subroutine mld_c_base_onelev_cnv(lv,info,amold,vmold,imold)
& call lv%ac%cscnv(info,mold=amold)
if (info == psb_success_ .and. lv%desc_ac%is_ok() &
& .and. present(imold)) call lv%desc_ac%cnv(imold)
call lv%map%cnv(info,mold=amold,imold=imold)
if (info == psb_success_) call lv%map%cnv(info,mold=amold,imold=imold)
end if
end subroutine mld_c_base_onelev_cnv

26
mlprec/impl/level/mld_c_base_onelev_dump.f90
Unescape Escape View File

@ -36,7 +36,7 @@
!
!
subroutine mld_c_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
& smoother,solver,global_num)
& smoother,solver,tprol,global_num)
use psb_base_mod
use mld_c_onelev_mod, mld_protect_name => mld_c_base_onelev_dump
@ -45,13 +45,13 @@ subroutine mld_c_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
integer(psb_ipk_), intent(in) :: level
integer(psb_ipk_), intent(out) :: info
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: ac, rp, smoother, solver,global_num
logical, optional, intent(in) :: ac, rp, smoother, solver, tprol, global_num
! Local variables
integer(psb_ipk_) :: i, j, il1, iln, lname, lev
integer(psb_ipk_) :: icontxt,iam, np
character(len=80) :: prefix_
character(len=120) :: fname ! len should be at least 20 more than
logical :: ac_, rp_,global_num_
logical :: ac_, rp_, tprol_, global_num_
integer(psb_ipk_), allocatable :: ivr(:), ivc(:)
info = 0
@ -79,6 +79,11 @@ subroutine mld_c_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
else
rp_ = .false.
end if
if (present(tprol)) then
tprol_ = tprol
else
tprol_ = .false.
end if
if (present(global_num)) then
global_num_ = global_num
else
@ -104,6 +109,15 @@ subroutine mld_c_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_p.mtx'
call lv%map%map_Y2X%print(fname,head=head,ivr=ivr,ivc=ivc)
end if
if (tprol_) then
! Tentative prolongator is stored with column indices already
! in global numbering, so only IVR is needed.
ivr = lv%map%p_desc_X%get_global_indices(owned=.false.)
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_tprol.mtx'
!
! This is not implemented yet.
!call lv%tprol%print(fname,head=head,ivr=ivr)
end if
end if
else
if (level >= 2) then
@ -117,6 +131,12 @@ subroutine mld_c_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_p.mtx'
call lv%map%map_Y2X%print(fname,head=head)
end if
if (tprol_) then
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_tprol.mtx'
!
! This is not implemented yet.
!call lv%tprol%print(fname,head=head)
end if
end if
end if

1
mlprec/impl/level/mld_c_base_onelev_setag.f90
Unescape Escape View File

@ -75,6 +75,7 @@ subroutine mld_c_base_onelev_setag(lv,val,info,pos)
end if
lv%parms%par_aggr_alg = mld_ext_aggr_
lv%parms%aggr_type = mld_noalg_
call lv%aggr%default()
end if
end subroutine mld_c_base_onelev_setag

2
mlprec/impl/level/mld_d_base_onelev_cnv.f90
Unescape Escape View File

@ -62,6 +62,6 @@ subroutine mld_d_base_onelev_cnv(lv,info,amold,vmold,imold)
& call lv%ac%cscnv(info,mold=amold)
if (info == psb_success_ .and. lv%desc_ac%is_ok() &
& .and. present(imold)) call lv%desc_ac%cnv(imold)
call lv%map%cnv(info,mold=amold,imold=imold)
if (info == psb_success_) call lv%map%cnv(info,mold=amold,imold=imold)
end if
end subroutine mld_d_base_onelev_cnv

26
mlprec/impl/level/mld_d_base_onelev_dump.f90
Unescape Escape View File

@ -36,7 +36,7 @@
!
!
subroutine mld_d_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
& smoother,solver,global_num)
& smoother,solver,tprol,global_num)
use psb_base_mod
use mld_d_onelev_mod, mld_protect_name => mld_d_base_onelev_dump
@ -45,13 +45,13 @@ subroutine mld_d_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
integer(psb_ipk_), intent(in) :: level
integer(psb_ipk_), intent(out) :: info
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: ac, rp, smoother, solver,global_num
logical, optional, intent(in) :: ac, rp, smoother, solver, tprol, global_num
! Local variables
integer(psb_ipk_) :: i, j, il1, iln, lname, lev
integer(psb_ipk_) :: icontxt,iam, np
character(len=80) :: prefix_
character(len=120) :: fname ! len should be at least 20 more than
logical :: ac_, rp_,global_num_
logical :: ac_, rp_, tprol_, global_num_
integer(psb_ipk_), allocatable :: ivr(:), ivc(:)
info = 0
@ -79,6 +79,11 @@ subroutine mld_d_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
else
rp_ = .false.
end if
if (present(tprol)) then
tprol_ = tprol
else
tprol_ = .false.
end if
if (present(global_num)) then
global_num_ = global_num
else
@ -104,6 +109,15 @@ subroutine mld_d_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_p.mtx'
call lv%map%map_Y2X%print(fname,head=head,ivr=ivr,ivc=ivc)
end if
if (tprol_) then
! Tentative prolongator is stored with column indices already
! in global numbering, so only IVR is needed.
ivr = lv%map%p_desc_X%get_global_indices(owned=.false.)
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_tprol.mtx'
!
! This is not implemented yet.
!call lv%tprol%print(fname,head=head,ivr=ivr)
end if
end if
else
if (level >= 2) then
@ -117,6 +131,12 @@ subroutine mld_d_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_p.mtx'
call lv%map%map_Y2X%print(fname,head=head)
end if
if (tprol_) then
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_tprol.mtx'
!
! This is not implemented yet.
!call lv%tprol%print(fname,head=head)
end if
end if
end if

1
mlprec/impl/level/mld_d_base_onelev_setag.f90
Unescape Escape View File

@ -75,6 +75,7 @@ subroutine mld_d_base_onelev_setag(lv,val,info,pos)
end if
lv%parms%par_aggr_alg = mld_ext_aggr_
lv%parms%aggr_type = mld_noalg_
call lv%aggr%default()
end if
end subroutine mld_d_base_onelev_setag

2
mlprec/impl/level/mld_s_base_onelev_cnv.f90
Unescape Escape View File

@ -62,6 +62,6 @@ subroutine mld_s_base_onelev_cnv(lv,info,amold,vmold,imold)
& call lv%ac%cscnv(info,mold=amold)
if (info == psb_success_ .and. lv%desc_ac%is_ok() &
& .and. present(imold)) call lv%desc_ac%cnv(imold)
call lv%map%cnv(info,mold=amold,imold=imold)
if (info == psb_success_) call lv%map%cnv(info,mold=amold,imold=imold)
end if
end subroutine mld_s_base_onelev_cnv

26
mlprec/impl/level/mld_s_base_onelev_dump.f90
Unescape Escape View File

@ -36,7 +36,7 @@
!
!
subroutine mld_s_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
& smoother,solver,global_num)
& smoother,solver,tprol,global_num)
use psb_base_mod
use mld_s_onelev_mod, mld_protect_name => mld_s_base_onelev_dump
@ -45,13 +45,13 @@ subroutine mld_s_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
integer(psb_ipk_), intent(in) :: level
integer(psb_ipk_), intent(out) :: info
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: ac, rp, smoother, solver,global_num
logical, optional, intent(in) :: ac, rp, smoother, solver, tprol, global_num
! Local variables
integer(psb_ipk_) :: i, j, il1, iln, lname, lev
integer(psb_ipk_) :: icontxt,iam, np
character(len=80) :: prefix_
character(len=120) :: fname ! len should be at least 20 more than
logical :: ac_, rp_,global_num_
logical :: ac_, rp_, tprol_, global_num_
integer(psb_ipk_), allocatable :: ivr(:), ivc(:)
info = 0
@ -79,6 +79,11 @@ subroutine mld_s_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
else
rp_ = .false.
end if
if (present(tprol)) then
tprol_ = tprol
else
tprol_ = .false.
end if
if (present(global_num)) then
global_num_ = global_num
else
@ -104,6 +109,15 @@ subroutine mld_s_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_p.mtx'
call lv%map%map_Y2X%print(fname,head=head,ivr=ivr,ivc=ivc)
end if
if (tprol_) then
! Tentative prolongator is stored with column indices already
! in global numbering, so only IVR is needed.
ivr = lv%map%p_desc_X%get_global_indices(owned=.false.)
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_tprol.mtx'
!
! This is not implemented yet.
!call lv%tprol%print(fname,head=head,ivr=ivr)
end if
end if
else
if (level >= 2) then
@ -117,6 +131,12 @@ subroutine mld_s_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_p.mtx'
call lv%map%map_Y2X%print(fname,head=head)
end if
if (tprol_) then
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_tprol.mtx'
!
! This is not implemented yet.
!call lv%tprol%print(fname,head=head)
end if
end if
end if

1
mlprec/impl/level/mld_s_base_onelev_setag.f90
Unescape Escape View File

@ -75,6 +75,7 @@ subroutine mld_s_base_onelev_setag(lv,val,info,pos)
end if
lv%parms%par_aggr_alg = mld_ext_aggr_
lv%parms%aggr_type = mld_noalg_
call lv%aggr%default()
end if
end subroutine mld_s_base_onelev_setag

2
mlprec/impl/level/mld_z_base_onelev_cnv.f90
Unescape Escape View File

@ -62,6 +62,6 @@ subroutine mld_z_base_onelev_cnv(lv,info,amold,vmold,imold)
& call lv%ac%cscnv(info,mold=amold)
if (info == psb_success_ .and. lv%desc_ac%is_ok() &
& .and. present(imold)) call lv%desc_ac%cnv(imold)
call lv%map%cnv(info,mold=amold,imold=imold)
if (info == psb_success_) call lv%map%cnv(info,mold=amold,imold=imold)
end if
end subroutine mld_z_base_onelev_cnv

26
mlprec/impl/level/mld_z_base_onelev_dump.f90
Unescape Escape View File

@ -36,7 +36,7 @@
!
!
subroutine mld_z_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
& smoother,solver,global_num)
& smoother,solver,tprol,global_num)
use psb_base_mod
use mld_z_onelev_mod, mld_protect_name => mld_z_base_onelev_dump
@ -45,13 +45,13 @@ subroutine mld_z_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
integer(psb_ipk_), intent(in) :: level
integer(psb_ipk_), intent(out) :: info
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: ac, rp, smoother, solver,global_num
logical, optional, intent(in) :: ac, rp, smoother, solver, tprol, global_num
! Local variables
integer(psb_ipk_) :: i, j, il1, iln, lname, lev
integer(psb_ipk_) :: icontxt,iam, np
character(len=80) :: prefix_
character(len=120) :: fname ! len should be at least 20 more than
logical :: ac_, rp_,global_num_
logical :: ac_, rp_, tprol_, global_num_
integer(psb_ipk_), allocatable :: ivr(:), ivc(:)
info = 0
@ -79,6 +79,11 @@ subroutine mld_z_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
else
rp_ = .false.
end if
if (present(tprol)) then
tprol_ = tprol
else
tprol_ = .false.
end if
if (present(global_num)) then
global_num_ = global_num
else
@ -104,6 +109,15 @@ subroutine mld_z_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_p.mtx'
call lv%map%map_Y2X%print(fname,head=head,ivr=ivr,ivc=ivc)
end if
if (tprol_) then
! Tentative prolongator is stored with column indices already
! in global numbering, so only IVR is needed.
ivr = lv%map%p_desc_X%get_global_indices(owned=.false.)
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_tprol.mtx'
!
! This is not implemented yet.
!call lv%tprol%print(fname,head=head,ivr=ivr)
end if
end if
else
if (level >= 2) then
@ -117,6 +131,12 @@ subroutine mld_z_base_onelev_dump(lv,level,info,prefix,head,ac,rp,&
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_p.mtx'
call lv%map%map_Y2X%print(fname,head=head)
end if
if (tprol_) then
write(fname(lname+1:),'(a,i3.3,a)')'_l',level,'_tprol.mtx'
!
! This is not implemented yet.
!call lv%tprol%print(fname,head=head)
end if
end if
end if

1
mlprec/impl/level/mld_z_base_onelev_setag.f90
Unescape Escape View File

@ -75,6 +75,7 @@ subroutine mld_z_base_onelev_setag(lv,val,info,pos)
end if
lv%parms%par_aggr_alg = mld_ext_aggr_
lv%parms%aggr_type = mld_noalg_
call lv%aggr%default()
end if
end subroutine mld_z_base_onelev_setag

14
mlprec/impl/mld_c_extprol_bld.F90
Unescape Escape View File

@ -160,9 +160,9 @@ subroutine mld_c_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
! Check to ensure all procs have the same
!
newsz = -1
mxplevs = p%max_levs
mnaggratio = p%min_cr_ratio
casize = p%min_coarse_size
mxplevs = p%ag_data%max_levs
mnaggratio = p%ag_data%min_cr_ratio
casize = p%ag_data%min_coarse_size
iszv = size(p%precv)
nprolv = size(prolv)
nrestrv = size(restrv)
@ -172,17 +172,17 @@ subroutine mld_c_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
call psb_bcast(ictxt,mnaggratio)
call psb_bcast(ictxt,nprolv)
call psb_bcast(ictxt,nrestrv)
if (casize /= p%min_coarse_size) then
if (casize /= p%ag_data%min_coarse_size) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_coarse_size')
goto 9999
end if
if (mxplevs /= p%max_levs) then
if (mxplevs /= p%ag_data%max_levs) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent max_levs')
goto 9999
end if
if (mnaggratio /= p%min_cr_ratio) then
if (mnaggratio /= p%ag_data%min_cr_ratio) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_cr_ratio')
goto 9999
@ -226,7 +226,7 @@ subroutine mld_c_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
!
nplevs = nrestrv + 1
p%max_levs = nplevs
p%ag_data%max_levs = nplevs
!
! Fixed number of levels.

16
mlprec/impl/mld_c_hierarchy_bld.f90
Unescape Escape View File

@ -124,25 +124,25 @@ subroutine mld_c_hierarchy_bld(a,desc_a,prec,info)
! Check to ensure all procs have the same
!
newsz = -1
mxplevs = prec%max_levs
mnaggratio = prec%min_cr_ratio
casize = prec%min_coarse_size
mxplevs = prec%ag_data%max_levs
mnaggratio = prec%ag_data%min_cr_ratio
casize = prec%ag_data%min_coarse_size
iszv = size(prec%precv)
call psb_bcast(ictxt,iszv)
call psb_bcast(ictxt,casize)
call psb_bcast(ictxt,mxplevs)
call psb_bcast(ictxt,mnaggratio)
if (casize /= prec%min_coarse_size) then
if (casize /= prec%ag_data%min_coarse_size) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_coarse_size')
goto 9999
end if
if (mxplevs /= prec%max_levs) then
if (mxplevs /= prec%ag_data%max_levs) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent max_levs')
goto 9999
end if
if (mnaggratio /= prec%min_cr_ratio) then
if (mnaggratio /= prec%ag_data%min_cr_ratio) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_cr_ratio')
goto 9999
@ -196,7 +196,7 @@ subroutine mld_c_hierarchy_bld(a,desc_a,prec,info)
casize = max(casize,ione)
casize = casize*40_psb_ipk_
call psb_bcast(ictxt,casize)
prec%min_coarse_size = casize
prec%ag_data%min_coarse_size = casize
end if
nplevs = max(itwo,mxplevs)
@ -298,7 +298,7 @@ subroutine mld_c_hierarchy_bld(a,desc_a,prec,info)
if (info == psb_success_)&
& call prec%precv(i)%bld_tprol(prec%precv(i-1)%base_a,&
& prec%precv(i-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info)
& ilaggr,nlaggr,op_prol,prec%ag_data,info)
if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,&

12
mlprec/impl/mld_ccprecset.F90
Unescape Escape View File

@ -129,7 +129,7 @@ subroutine mld_ccprecseti(p,what,val,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1
@ -146,10 +146,10 @@ subroutine mld_ccprecseti(p,what,val,info,ilev,ilmax,pos,idx)
select case(psb_toupper(what))
case ('MIN_COARSE_SIZE')
p%min_coarse_size = max(val,-1)
p%ag_data%min_coarse_size = max(val,-1)
return
case('MAX_LEVS')
p%max_levs = max(val,1)
p%ag_data%max_levs = max(val,1)
return
case ('OUTER_SWEEPS')
p%outer_sweeps = max(val,1)
@ -502,7 +502,7 @@ subroutine mld_ccprecsetc(p,what,string,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1
@ -593,7 +593,7 @@ subroutine mld_ccprecsetr(p,what,val,info,ilev,ilmax,pos,idx)
select case(psb_toupper(what))
case ('MIN_CR_RATIO')
p%min_cr_ratio = max(sone,val)
p%ag_data%min_cr_ratio = max(sone,val)
return
end select
@ -615,7 +615,7 @@ subroutine mld_ccprecsetr(p,what,val,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1

4
mlprec/impl/mld_cprecinit.F90
Unescape Escape View File

@ -119,7 +119,7 @@ subroutine mld_cprecinit(ictxt,prec,ptype,info)
endif
endif
prec%ictxt = ictxt
prec%min_coarse_size = -1
prec%ag_data%min_coarse_size = -1
select case(psb_toupper(ptype(1:len_trim(ptype))))
case ('NOPREC','NONE')
@ -186,7 +186,7 @@ subroutine mld_cprecinit(ictxt,prec,ptype,info)
case ('ML')
nlev_ = prec%max_levs
nlev_ = prec%ag_data%max_levs
ilev_ = 1
allocate(prec%precv(nlev_),stat=info)

14
mlprec/impl/mld_d_extprol_bld.F90
Unescape Escape View File

@ -160,9 +160,9 @@ subroutine mld_d_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
! Check to ensure all procs have the same
!
newsz = -1
mxplevs = p%max_levs
mnaggratio = p%min_cr_ratio
casize = p%min_coarse_size
mxplevs = p%ag_data%max_levs
mnaggratio = p%ag_data%min_cr_ratio
casize = p%ag_data%min_coarse_size
iszv = size(p%precv)
nprolv = size(prolv)
nrestrv = size(restrv)
@ -172,17 +172,17 @@ subroutine mld_d_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
call psb_bcast(ictxt,mnaggratio)
call psb_bcast(ictxt,nprolv)
call psb_bcast(ictxt,nrestrv)
if (casize /= p%min_coarse_size) then
if (casize /= p%ag_data%min_coarse_size) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_coarse_size')
goto 9999
end if
if (mxplevs /= p%max_levs) then
if (mxplevs /= p%ag_data%max_levs) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent max_levs')
goto 9999
end if
if (mnaggratio /= p%min_cr_ratio) then
if (mnaggratio /= p%ag_data%min_cr_ratio) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_cr_ratio')
goto 9999
@ -226,7 +226,7 @@ subroutine mld_d_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
!
nplevs = nrestrv + 1
p%max_levs = nplevs
p%ag_data%max_levs = nplevs
!
! Fixed number of levels.

16
mlprec/impl/mld_d_hierarchy_bld.f90
Unescape Escape View File

@ -124,25 +124,25 @@ subroutine mld_d_hierarchy_bld(a,desc_a,prec,info)
! Check to ensure all procs have the same
!
newsz = -1
mxplevs = prec%max_levs
mnaggratio = prec%min_cr_ratio
casize = prec%min_coarse_size
mxplevs = prec%ag_data%max_levs
mnaggratio = prec%ag_data%min_cr_ratio
casize = prec%ag_data%min_coarse_size
iszv = size(prec%precv)
call psb_bcast(ictxt,iszv)
call psb_bcast(ictxt,casize)
call psb_bcast(ictxt,mxplevs)
call psb_bcast(ictxt,mnaggratio)
if (casize /= prec%min_coarse_size) then
if (casize /= prec%ag_data%min_coarse_size) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_coarse_size')
goto 9999
end if
if (mxplevs /= prec%max_levs) then
if (mxplevs /= prec%ag_data%max_levs) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent max_levs')
goto 9999
end if
if (mnaggratio /= prec%min_cr_ratio) then
if (mnaggratio /= prec%ag_data%min_cr_ratio) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_cr_ratio')
goto 9999
@ -196,7 +196,7 @@ subroutine mld_d_hierarchy_bld(a,desc_a,prec,info)
casize = max(casize,ione)
casize = casize*40_psb_ipk_
call psb_bcast(ictxt,casize)
prec%min_coarse_size = casize
prec%ag_data%min_coarse_size = casize
end if
nplevs = max(itwo,mxplevs)
@ -298,7 +298,7 @@ subroutine mld_d_hierarchy_bld(a,desc_a,prec,info)
if (info == psb_success_)&
& call prec%precv(i)%bld_tprol(prec%precv(i-1)%base_a,&
& prec%precv(i-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info)
& ilaggr,nlaggr,op_prol,prec%ag_data,info)
if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,&

12
mlprec/impl/mld_dcprecset.F90
Unescape Escape View File

@ -135,7 +135,7 @@ subroutine mld_dcprecseti(p,what,val,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1
@ -152,10 +152,10 @@ subroutine mld_dcprecseti(p,what,val,info,ilev,ilmax,pos,idx)
select case(psb_toupper(what))
case ('MIN_COARSE_SIZE')
p%min_coarse_size = max(val,-1)
p%ag_data%min_coarse_size = max(val,-1)
return
case('MAX_LEVS')
p%max_levs = max(val,1)
p%ag_data%max_levs = max(val,1)
return
case ('OUTER_SWEEPS')
p%outer_sweeps = max(val,1)
@ -536,7 +536,7 @@ subroutine mld_dcprecsetc(p,what,string,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1
@ -627,7 +627,7 @@ subroutine mld_dcprecsetr(p,what,val,info,ilev,ilmax,pos,idx)
select case(psb_toupper(what))
case ('MIN_CR_RATIO')
p%min_cr_ratio = max(done,val)
p%ag_data%min_cr_ratio = max(done,val)
return
end select
@ -649,7 +649,7 @@ subroutine mld_dcprecsetr(p,what,val,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1

4
mlprec/impl/mld_dprecinit.F90
Unescape Escape View File

@ -122,7 +122,7 @@ subroutine mld_dprecinit(ictxt,prec,ptype,info)
endif
endif
prec%ictxt = ictxt
prec%min_coarse_size = -1
prec%ag_data%min_coarse_size = -1
select case(psb_toupper(ptype(1:len_trim(ptype))))
case ('NOPREC','NONE')
@ -189,7 +189,7 @@ subroutine mld_dprecinit(ictxt,prec,ptype,info)
case ('ML')
nlev_ = prec%max_levs
nlev_ = prec%ag_data%max_levs
ilev_ = 1
allocate(prec%precv(nlev_),stat=info)

14
mlprec/impl/mld_s_extprol_bld.F90
Unescape Escape View File

@ -160,9 +160,9 @@ subroutine mld_s_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
! Check to ensure all procs have the same
!
newsz = -1
mxplevs = p%max_levs
mnaggratio = p%min_cr_ratio
casize = p%min_coarse_size
mxplevs = p%ag_data%max_levs
mnaggratio = p%ag_data%min_cr_ratio
casize = p%ag_data%min_coarse_size
iszv = size(p%precv)
nprolv = size(prolv)
nrestrv = size(restrv)
@ -172,17 +172,17 @@ subroutine mld_s_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
call psb_bcast(ictxt,mnaggratio)
call psb_bcast(ictxt,nprolv)
call psb_bcast(ictxt,nrestrv)
if (casize /= p%min_coarse_size) then
if (casize /= p%ag_data%min_coarse_size) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_coarse_size')
goto 9999
end if
if (mxplevs /= p%max_levs) then
if (mxplevs /= p%ag_data%max_levs) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent max_levs')
goto 9999
end if
if (mnaggratio /= p%min_cr_ratio) then
if (mnaggratio /= p%ag_data%min_cr_ratio) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_cr_ratio')
goto 9999
@ -226,7 +226,7 @@ subroutine mld_s_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
!
nplevs = nrestrv + 1
p%max_levs = nplevs
p%ag_data%max_levs = nplevs
!
! Fixed number of levels.

16
mlprec/impl/mld_s_hierarchy_bld.f90
Unescape Escape View File

@ -124,25 +124,25 @@ subroutine mld_s_hierarchy_bld(a,desc_a,prec,info)
! Check to ensure all procs have the same
!
newsz = -1
mxplevs = prec%max_levs
mnaggratio = prec%min_cr_ratio
casize = prec%min_coarse_size
mxplevs = prec%ag_data%max_levs
mnaggratio = prec%ag_data%min_cr_ratio
casize = prec%ag_data%min_coarse_size
iszv = size(prec%precv)
call psb_bcast(ictxt,iszv)
call psb_bcast(ictxt,casize)
call psb_bcast(ictxt,mxplevs)
call psb_bcast(ictxt,mnaggratio)
if (casize /= prec%min_coarse_size) then
if (casize /= prec%ag_data%min_coarse_size) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_coarse_size')
goto 9999
end if
if (mxplevs /= prec%max_levs) then
if (mxplevs /= prec%ag_data%max_levs) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent max_levs')
goto 9999
end if
if (mnaggratio /= prec%min_cr_ratio) then
if (mnaggratio /= prec%ag_data%min_cr_ratio) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_cr_ratio')
goto 9999
@ -196,7 +196,7 @@ subroutine mld_s_hierarchy_bld(a,desc_a,prec,info)
casize = max(casize,ione)
casize = casize*40_psb_ipk_
call psb_bcast(ictxt,casize)
prec%min_coarse_size = casize
prec%ag_data%min_coarse_size = casize
end if
nplevs = max(itwo,mxplevs)
@ -298,7 +298,7 @@ subroutine mld_s_hierarchy_bld(a,desc_a,prec,info)
if (info == psb_success_)&
& call prec%precv(i)%bld_tprol(prec%precv(i-1)%base_a,&
& prec%precv(i-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info)
& ilaggr,nlaggr,op_prol,prec%ag_data,info)
if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,&

12
mlprec/impl/mld_scprecset.F90
Unescape Escape View File

@ -129,7 +129,7 @@ subroutine mld_scprecseti(p,what,val,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1
@ -146,10 +146,10 @@ subroutine mld_scprecseti(p,what,val,info,ilev,ilmax,pos,idx)
select case(psb_toupper(what))
case ('MIN_COARSE_SIZE')
p%min_coarse_size = max(val,-1)
p%ag_data%min_coarse_size = max(val,-1)
return
case('MAX_LEVS')
p%max_levs = max(val,1)
p%ag_data%max_levs = max(val,1)
return
case ('OUTER_SWEEPS')
p%outer_sweeps = max(val,1)
@ -502,7 +502,7 @@ subroutine mld_scprecsetc(p,what,string,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1
@ -593,7 +593,7 @@ subroutine mld_scprecsetr(p,what,val,info,ilev,ilmax,pos,idx)
select case(psb_toupper(what))
case ('MIN_CR_RATIO')
p%min_cr_ratio = max(sone,val)
p%ag_data%min_cr_ratio = max(sone,val)
return
end select
@ -615,7 +615,7 @@ subroutine mld_scprecsetr(p,what,val,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1

4
mlprec/impl/mld_sprecinit.F90
Unescape Escape View File

@ -119,7 +119,7 @@ subroutine mld_sprecinit(ictxt,prec,ptype,info)
endif
endif
prec%ictxt = ictxt
prec%min_coarse_size = -1
prec%ag_data%min_coarse_size = -1
select case(psb_toupper(ptype(1:len_trim(ptype))))
case ('NOPREC','NONE')
@ -186,7 +186,7 @@ subroutine mld_sprecinit(ictxt,prec,ptype,info)
case ('ML')
nlev_ = prec%max_levs
nlev_ = prec%ag_data%max_levs
ilev_ = 1
allocate(prec%precv(nlev_),stat=info)

14
mlprec/impl/mld_z_extprol_bld.F90
Unescape Escape View File

@ -160,9 +160,9 @@ subroutine mld_z_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
! Check to ensure all procs have the same
!
newsz = -1
mxplevs = p%max_levs
mnaggratio = p%min_cr_ratio
casize = p%min_coarse_size
mxplevs = p%ag_data%max_levs
mnaggratio = p%ag_data%min_cr_ratio
casize = p%ag_data%min_coarse_size
iszv = size(p%precv)
nprolv = size(prolv)
nrestrv = size(restrv)
@ -172,17 +172,17 @@ subroutine mld_z_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
call psb_bcast(ictxt,mnaggratio)
call psb_bcast(ictxt,nprolv)
call psb_bcast(ictxt,nrestrv)
if (casize /= p%min_coarse_size) then
if (casize /= p%ag_data%min_coarse_size) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_coarse_size')
goto 9999
end if
if (mxplevs /= p%max_levs) then
if (mxplevs /= p%ag_data%max_levs) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent max_levs')
goto 9999
end if
if (mnaggratio /= p%min_cr_ratio) then
if (mnaggratio /= p%ag_data%min_cr_ratio) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_cr_ratio')
goto 9999
@ -226,7 +226,7 @@ subroutine mld_z_extprol_bld(a,desc_a,p,prolv,restrv,info,amold,vmold,imold)
!
nplevs = nrestrv + 1
p%max_levs = nplevs
p%ag_data%max_levs = nplevs
!
! Fixed number of levels.

16
mlprec/impl/mld_z_hierarchy_bld.f90
Unescape Escape View File

@ -124,25 +124,25 @@ subroutine mld_z_hierarchy_bld(a,desc_a,prec,info)
! Check to ensure all procs have the same
!
newsz = -1
mxplevs = prec%max_levs
mnaggratio = prec%min_cr_ratio
casize = prec%min_coarse_size
mxplevs = prec%ag_data%max_levs
mnaggratio = prec%ag_data%min_cr_ratio
casize = prec%ag_data%min_coarse_size
iszv = size(prec%precv)
call psb_bcast(ictxt,iszv)
call psb_bcast(ictxt,casize)
call psb_bcast(ictxt,mxplevs)
call psb_bcast(ictxt,mnaggratio)
if (casize /= prec%min_coarse_size) then
if (casize /= prec%ag_data%min_coarse_size) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_coarse_size')
goto 9999
end if
if (mxplevs /= prec%max_levs) then
if (mxplevs /= prec%ag_data%max_levs) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent max_levs')
goto 9999
end if
if (mnaggratio /= prec%min_cr_ratio) then
if (mnaggratio /= prec%ag_data%min_cr_ratio) then
info=psb_err_internal_error_
call psb_errpush(info,name,a_err='Inconsistent min_cr_ratio')
goto 9999
@ -196,7 +196,7 @@ subroutine mld_z_hierarchy_bld(a,desc_a,prec,info)
casize = max(casize,ione)
casize = casize*40_psb_ipk_
call psb_bcast(ictxt,casize)
prec%min_coarse_size = casize
prec%ag_data%min_coarse_size = casize
end if
nplevs = max(itwo,mxplevs)
@ -298,7 +298,7 @@ subroutine mld_z_hierarchy_bld(a,desc_a,prec,info)
if (info == psb_success_)&
& call prec%precv(i)%bld_tprol(prec%precv(i-1)%base_a,&
& prec%precv(i-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info)
& ilaggr,nlaggr,op_prol,prec%ag_data,info)
if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,&

12
mlprec/impl/mld_zcprecset.F90
Unescape Escape View File

@ -135,7 +135,7 @@ subroutine mld_zcprecseti(p,what,val,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1
@ -152,10 +152,10 @@ subroutine mld_zcprecseti(p,what,val,info,ilev,ilmax,pos,idx)
select case(psb_toupper(what))
case ('MIN_COARSE_SIZE')
p%min_coarse_size = max(val,-1)
p%ag_data%min_coarse_size = max(val,-1)
return
case('MAX_LEVS')
p%max_levs = max(val,1)
p%ag_data%max_levs = max(val,1)
return
case ('OUTER_SWEEPS')
p%outer_sweeps = max(val,1)
@ -536,7 +536,7 @@ subroutine mld_zcprecsetc(p,what,string,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1
@ -627,7 +627,7 @@ subroutine mld_zcprecsetr(p,what,val,info,ilev,ilmax,pos,idx)
select case(psb_toupper(what))
case ('MIN_CR_RATIO')
p%min_cr_ratio = max(done,val)
p%ag_data%min_cr_ratio = max(done,val)
return
end select
@ -649,7 +649,7 @@ subroutine mld_zcprecsetr(p,what,val,info,ilev,ilmax,pos,idx)
end if
else
ilev_ = 1
ilmax_ = ilev_
ilmax_ = nlev_
end if
if ((ilev_<1).or.(ilev_ > nlev_)) then
info = -1

4
mlprec/impl/mld_zprecinit.F90
Unescape Escape View File

@ -122,7 +122,7 @@ subroutine mld_zprecinit(ictxt,prec,ptype,info)
endif
endif
prec%ictxt = ictxt
prec%min_coarse_size = -1
prec%ag_data%min_coarse_size = -1
select case(psb_toupper(ptype(1:len_trim(ptype))))
case ('NOPREC','NONE')
@ -189,7 +189,7 @@ subroutine mld_zprecinit(ictxt,prec,ptype,info)
case ('ML')
nlev_ = prec%max_levs
nlev_ = prec%ag_data%max_levs
ilev_ = 1
allocate(prec%precv(nlev_),stat=info)

38
mlprec/mld_base_prec_type.F90
Unescape Escape View File

@ -120,6 +120,44 @@ module mld_base_prec_type
procedure, pass(pm) :: printout => d_ml_parms_printout
end type mld_dml_parms
type mld_saggr_data
!
! Aggregation data and defaults:
!
! 1. min_coarse_size = 0 Default target size will be computed as
! 40*(N_fine)**(1./3.)
! We are assuming that the coarse size fits in
! integer range of psb_ipk_, but this is
! not very restrictive
integer(psb_ipk_) :: min_coarse_size = izero
! 2. maximum number of levels. Defaults to 20
integer(psb_ipk_) :: max_levs = 20_psb_ipk_
! 3. min_cr_ratio = 1.5
real(psb_spk_) :: min_cr_ratio = 1.5_psb_spk_
real(psb_spk_) :: op_complexity = szero
real(psb_spk_) :: avg_cr = szero
end type mld_saggr_data
type mld_daggr_data
!
! Aggregation data and defaults:
!
!
! 1. min_coarse_size = 0 Default target size will be computed as
! 40*(N_fine)**(1./3.)
! We are assuming that the coarse size fits in
! integer range of psb_ipk_, but this is
! not very restrictive
integer(psb_ipk_) :: min_coarse_size = izero
! 2. maximum number of levels. Defaults to 20
integer(psb_ipk_) :: max_levs = 20_psb_ipk_
! 3. min_cr_ratio = 1.5
real(psb_dpk_) :: min_cr_ratio = 1.5_psb_dpk_
real(psb_dpk_) :: op_complexity = dzero
real(psb_dpk_) :: avg_cr = dzero
end type mld_daggr_data
!
! Entries in iprcparm

44
mlprec/mld_c_base_aggregator_mod.f90
Unescape Escape View File

@ -41,11 +41,11 @@
!
module mld_c_base_aggregator_mod
use mld_base_prec_type, only : mld_sml_parms
use mld_base_prec_type, only : mld_sml_parms, mld_saggr_data
use psb_base_mod, only : psb_cspmat_type, psb_lcspmat_type, psb_c_vect_type, &
& psb_c_base_vect_type, psb_clinmap_type, psb_spk_, &
& psb_ipk_, psb_epk_, psb_lpk_, psb_desc_type, psb_i_base_vect_type, &
& psb_erractionsave, psb_error_handler, psb_success_
& psb_erractionsave, psb_error_handler, psb_success_, psb_toupper
!
!
!
@ -79,7 +79,8 @@ module mld_c_base_aggregator_mod
!! cseti, csetr, csetc - Set internal parameters, if any
!
type mld_c_base_aggregator_type
! Do we want to purge explicit zeros when aggregating?
logical :: do_clean_zeros
contains
procedure, pass(ag) :: bld_tprol => mld_c_base_aggregator_build_tprol
procedure, pass(ag) :: mat_asb => mld_c_base_aggregator_mat_asb
@ -96,6 +97,19 @@ module mld_c_base_aggregator_mod
generic, public :: set => cseti, csetr, csetc
end type mld_c_base_aggregator_type
abstract interface
subroutine mld_c_soc_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
import :: psb_cspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_
implicit none
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_spk_), intent(in) :: theta
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
integer(psb_ipk_), intent(out) :: info
end subroutine mld_c_soc_map_bld
end interface
contains
@ -137,7 +151,16 @@ contains
character(len=*), intent(in) :: val
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: idx
! Do nothing
! Set clean zeros, or do nothing.
select case (psb_toupper(trim(what)))
case('AGGR_CLEAN_ZEROS')
select case (psb_toupper(trim(val)))
case('TRUE','T')
ag%do_clean_zeros = .true.
case('FALSE','F')
ag%do_clean_zeros = .false.
end select
end select
info = 0
end subroutine mld_c_base_aggregator_csetc
@ -181,8 +204,8 @@ contains
subroutine mld_c_base_aggregator_default(ag)
implicit none
class(mld_c_base_aggregator_type), intent(inout) :: ag
! Here we need do nothing
! Only one default setting
ag%do_clean_zeros = .true.
return
end subroutine mld_c_base_aggregator_default
@ -228,9 +251,12 @@ contains
!! will contribute, in global numbering.
!! Many aggregation produce a binary tentative prolongator, but some
!! do not, hence we also need the OP_PROL output.
!! AG_DATA is passed here just in case some of the
!! aggregators need it internally, most of them will ignore.
!!
!! \param ag The input aggregator object
!! \param parms The auxiliary parameters object
!! \param ag_data Auxiliary global aggregation info
!! \param a The local matrix part
!! \param desc_a The descriptor
!! \param ilaggr Output aggregation map
@ -239,11 +265,13 @@ contains
!! \param info Return code
!!
!
subroutine mld_c_base_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_c_base_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
implicit none
class(mld_c_base_aggregator_type), target, intent(inout) :: ag
type(mld_sml_parms), intent(inout) :: parms
type(mld_sml_parms), intent(inout) :: parms
type(mld_saggr_data), intent(in) :: ag_data
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)

38
mlprec/mld_c_dec_aggregator_mod.f90
Unescape Escape View File

@ -48,7 +48,7 @@ module mld_c_dec_aggregator_mod
!! \extends mld_c_base_aggregator_mod::mld_c_base_aggregator_type
!!
!! type, extends(mld_c_base_aggregator_type) :: mld_c_dec_aggregator_type
!! procedure(mld_c_map_bld), nopass, pointer :: map_bld => null()
!! procedure(mld_c_soc_map_bld), nopass, pointer :: soc_map_bld => null()
!! end type
!!
!! This is the simplest aggregation method: starting from the
@ -71,12 +71,12 @@ module mld_c_dec_aggregator_mod
!! PSBLAS-based parallel two-level Schwarz preconditioners, Appl. Num. Math.
!! 57 (2007), 1181-1196.
!!
!! The map_bld method is used inside the implementation of build_tprol
!! The soc_map_bld method is used inside the implementation of build_tprol
!!
!
!
type, extends(mld_c_base_aggregator_type) :: mld_c_dec_aggregator_type
procedure(mld_c_map_bld), nopass, pointer :: map_bld => null()
procedure(mld_c_soc_map_bld), nopass, pointer :: soc_map_bld => null()
contains
procedure, pass(ag) :: bld_tprol => mld_c_dec_aggregator_build_tprol
@ -88,28 +88,17 @@ module mld_c_dec_aggregator_mod
end type mld_c_dec_aggregator_type
abstract interface
subroutine mld_c_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
import :: psb_cspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_
implicit none
integer(psb_ipk_), intent(in) :: iorder
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_spk_), intent(in) :: theta
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
integer(psb_ipk_), intent(out) :: info
end subroutine mld_c_map_bld
end interface
procedure(mld_c_map_bld) :: mld_c_soc1_map_bld, mld_c_soc2_map_bld
procedure(mld_c_soc_map_bld) :: mld_c_soc1_map_bld, mld_c_soc2_map_bld
interface
subroutine mld_c_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_c_dec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: mld_c_dec_aggregator_type, psb_desc_type, psb_cspmat_type, psb_spk_, &
& psb_ipk_, psb_lpk_, psb_lcspmat_type, mld_sml_parms
& psb_ipk_, psb_lpk_, psb_lcspmat_type, mld_sml_parms, mld_saggr_data
implicit none
class(mld_c_dec_aggregator_type), target, intent(inout) :: ag
type(mld_sml_parms), intent(inout) :: parms
type(mld_saggr_data), intent(in) :: ag_data
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -147,14 +136,14 @@ contains
select case(parms%aggr_type)
case (mld_noalg_)
ag%map_bld => null()
ag%soc_map_bld => null()
case (mld_soc1_)
ag%map_bld => mld_c_soc1_map_bld
ag%soc_map_bld => mld_c_soc1_map_bld
case (mld_soc2_)
ag%map_bld => mld_c_soc2_map_bld
ag%soc_map_bld => mld_c_soc2_map_bld
case default
write(0,*) 'Unknown aggregation type, defaulting to SOC1'
ag%map_bld => mld_c_soc1_map_bld
ag%soc_map_bld => mld_c_soc1_map_bld
end select
return
@ -165,7 +154,8 @@ contains
implicit none
class(mld_c_dec_aggregator_type), intent(inout) :: ag
ag%map_bld => mld_c_soc1_map_bld
call ag%mld_c_base_aggregator_type%default()
ag%soc_map_bld => mld_c_soc1_map_bld
return
end subroutine mld_c_dec_aggregator_default

26
mlprec/mld_c_inner_mod.f90
Unescape Escape View File

@ -97,32 +97,6 @@ module mld_c_inner_mod
end subroutine mld_cmlprec_aply_vect
end interface mld_mlprec_aply
interface mld_aggrmap_bld
subroutine mld_c_lev_aggrmap_bld(p,a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: psb_cspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lcspmat_type
import :: mld_c_onelev_type
implicit none
type(mld_c_onelev_type), intent(inout), target :: p
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), intent(out) :: ilaggr(:),nlaggr(:)
type(psb_lcspmat_type), intent(out) :: op_prol
integer(psb_ipk_), intent(out) :: info
end subroutine mld_c_lev_aggrmap_bld
subroutine mld_caggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: psb_cspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lcspmat_type
implicit none
integer(psb_ipk_), intent(in) :: iorder
integer(psb_ipk_), intent(in) :: aggr_type
real(psb_spk_), intent(in) :: theta
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
type(psb_lcspmat_type), intent(out) :: op_prol
integer(psb_ipk_), intent(out) :: info
end subroutine mld_caggrmap_bld
end interface mld_aggrmap_bld
interface mld_map_to_tprol
subroutine mld_c_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod, only : psb_cspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lcspmat_type

10
mlprec/mld_c_onelev_mod.f90
Unescape Escape View File

@ -381,7 +381,7 @@ interface
interface
subroutine mld_c_base_onelev_dump(lv,level,info,prefix,head,ac,rp,smoother,&
& solver,global_num)
& solver,tprol,global_num)
import :: psb_cspmat_type, psb_c_vect_type, psb_c_base_vect_type, &
& psb_clinmap_type, psb_spk_, mld_c_onelev_type, &
& psb_ipk_, psb_epk_, psb_desc_type
@ -390,7 +390,7 @@ interface
integer(psb_ipk_), intent(in) :: level
integer(psb_ipk_), intent(out) :: info
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: ac, rp, smoother, solver, global_num
logical, optional, intent(in) :: ac, rp, smoother, solver, tprol, global_num
end subroutine mld_c_base_onelev_dump
end interface
@ -484,16 +484,18 @@ contains
end subroutine c_base_onelev_default
subroutine c_base_onelev_bld_tprol(lv,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine c_base_onelev_bld_tprol(lv,a,desc_a,&
& ilaggr,nlaggr,op_prol,ag_data,info)
implicit none
class(mld_c_onelev_type), intent(inout), target :: lv
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
type(psb_lcspmat_type), intent(out) :: op_prol
type(mld_saggr_data), intent(in) :: ag_data
integer(psb_ipk_), intent(out) :: info
call lv%aggr%bld_tprol(lv%parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
call lv%aggr%bld_tprol(lv%parms,ag_data,a,desc_a,ilaggr,nlaggr,op_prol,info)
end subroutine c_base_onelev_bld_tprol

45
mlprec/mld_c_prec_type.f90
Unescape Escape View File

@ -85,18 +85,8 @@ module mld_c_prec_type
integer, parameter, private :: wv_size_=4
type, extends(psb_cprec_type) :: mld_cprec_type
! integer(psb_ipk_) :: ictxt ! Now it's in the PSBLAS prec.
!
! Aggregation defaults:
!
! 1. min_coarse_size = 0 Default target size will be computed as 40*(N_fine)**(1./3.)
integer(psb_ipk_) :: min_coarse_size = izero
! 2. maximum number of levels. Defaults to 20
integer(psb_ipk_) :: max_levs = 20_psb_ipk_
! 3. min_cr_ratio = 1.5
real(psb_spk_) :: min_cr_ratio = 1.5_psb_spk_
real(psb_spk_) :: op_complexity = szero
real(psb_spk_) :: avg_cr = szero
! integer(psb_ipk_) :: ictxt ! Now it's in the PSBLAS prec.
type(mld_saggr_data) :: ag_data
!
! Number of outer sweeps. Sometimes 2 V-cycles may be better than 1 W-cycle.
!
@ -445,7 +435,7 @@ contains
class(mld_cprec_type), intent(in) :: prec
complex(psb_spk_) :: val
val = prec%op_complexity
val = prec%ag_data%op_complexity
end function mld_c_get_compl
@ -480,7 +470,7 @@ contains
call psb_sum(ictxt,num)
call psb_sum(ictxt,den)
end if
prec%op_complexity = num/den
prec%ag_data%op_complexity = num/den
end subroutine mld_c_cmp_compl
!
@ -492,7 +482,7 @@ contains
class(mld_cprec_type), intent(in) :: prec
complex(psb_spk_) :: val
val = prec%avg_cr
val = prec%ag_data%avg_cr
end function mld_c_get_avg_cr
@ -517,7 +507,7 @@ contains
avgcr = avgcr / (nl-1)
end if
call psb_sum(ictxt,avgcr)
prec%avg_cr = avgcr/np
prec%ag_data%avg_cr = avgcr/np
end subroutine mld_c_cmp_avg_cr
!
@ -722,14 +712,15 @@ contains
end subroutine mld_c_apply1v
subroutine mld_c_dump(prec,info,istart,iend,prefix,head,ac,rp,smoother,solver,global_num)
subroutine mld_c_dump(prec,info,istart,iend,prefix,head,ac,rp,smoother,solver,tprol,&
& global_num)
implicit none
class(mld_cprec_type), intent(in) :: prec
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: istart, iend
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: smoother, solver,ac, rp, global_num
logical, optional, intent(in) :: smoother, solver,ac, rp, tprol, global_num
integer(psb_ipk_) :: i, j, il1, iln, lname, lev
integer(psb_ipk_) :: icontxt,iam, np
character(len=80) :: prefix_
@ -750,7 +741,8 @@ contains
do lev=il1, iln
call prec%precv(lev)%dump(lev,info,prefix=prefix,head=head,&
& ac=ac,smoother=smoother,solver=solver,rp=rp,global_num=global_num)
& ac=ac,smoother=smoother,solver=solver,rp=rp,tprol=tprol, &
& global_num=global_num)
end do
end subroutine mld_c_dump
@ -801,13 +793,9 @@ contains
info = psb_success_
select type(pout => precout)
class is (mld_cprec_type)
pout%ictxt = prec%ictxt
pout%max_levs = prec%max_levs
pout%min_coarse_size = prec%min_coarse_size
pout%min_cr_ratio = prec%min_cr_ratio
pout%outer_sweeps = prec%outer_sweeps
pout%op_complexity = prec%op_complexity
pout%avg_cr = prec%avg_cr
pout%ictxt = prec%ictxt
pout%ag_data = prec%ag_data
pout%outer_sweeps = prec%outer_sweeps
if (allocated(prec%precv)) then
ln = size(prec%precv)
allocate(pout%precv(ln),stat=info)
@ -853,7 +841,10 @@ contains
!!$ return
endif
end if
b%ictxt = prec%ictxt
b%ag_data = prec%ag_data
b%outer_sweeps = prec%outer_sweeps
call move_alloc(prec%precv,b%precv)
! Fix the pointers except on level 1.
do i=2, size(b%precv)

6
mlprec/mld_c_symdec_aggregator_mod.f90
Unescape Escape View File

@ -62,12 +62,14 @@ module mld_c_symdec_aggregator_mod
interface
subroutine mld_c_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_c_symdec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: mld_c_symdec_aggregator_type, psb_desc_type, psb_cspmat_type, psb_spk_, &
& psb_ipk_, psb_lpk_, psb_lcspmat_type, mld_sml_parms
& psb_ipk_, psb_lpk_, psb_lcspmat_type, mld_sml_parms, mld_saggr_data
implicit none
class(mld_c_symdec_aggregator_type), target, intent(inout) :: ag
type(mld_sml_parms), intent(inout) :: parms
type(mld_saggr_data), intent(in) :: ag_data
type(psb_cspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)

44
mlprec/mld_d_base_aggregator_mod.f90
Unescape Escape View File

@ -41,11 +41,11 @@
!
module mld_d_base_aggregator_mod
use mld_base_prec_type, only : mld_dml_parms
use mld_base_prec_type, only : mld_dml_parms, mld_daggr_data
use psb_base_mod, only : psb_dspmat_type, psb_ldspmat_type, psb_d_vect_type, &
& psb_d_base_vect_type, psb_dlinmap_type, psb_dpk_, &
& psb_ipk_, psb_epk_, psb_lpk_, psb_desc_type, psb_i_base_vect_type, &
& psb_erractionsave, psb_error_handler, psb_success_
& psb_erractionsave, psb_error_handler, psb_success_, psb_toupper
!
!
!
@ -79,7 +79,8 @@ module mld_d_base_aggregator_mod
!! cseti, csetr, csetc - Set internal parameters, if any
!
type mld_d_base_aggregator_type
! Do we want to purge explicit zeros when aggregating?
logical :: do_clean_zeros
contains
procedure, pass(ag) :: bld_tprol => mld_d_base_aggregator_build_tprol
procedure, pass(ag) :: mat_asb => mld_d_base_aggregator_mat_asb
@ -96,6 +97,19 @@ module mld_d_base_aggregator_mod
generic, public :: set => cseti, csetr, csetc
end type mld_d_base_aggregator_type
abstract interface
subroutine mld_d_soc_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
import :: psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_
implicit none
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_dpk_), intent(in) :: theta
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
integer(psb_ipk_), intent(out) :: info
end subroutine mld_d_soc_map_bld
end interface
contains
@ -137,7 +151,16 @@ contains
character(len=*), intent(in) :: val
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: idx
! Do nothing
! Set clean zeros, or do nothing.
select case (psb_toupper(trim(what)))
case('AGGR_CLEAN_ZEROS')
select case (psb_toupper(trim(val)))
case('TRUE','T')
ag%do_clean_zeros = .true.
case('FALSE','F')
ag%do_clean_zeros = .false.
end select
end select
info = 0
end subroutine mld_d_base_aggregator_csetc
@ -181,8 +204,8 @@ contains
subroutine mld_d_base_aggregator_default(ag)
implicit none
class(mld_d_base_aggregator_type), intent(inout) :: ag
! Here we need do nothing
! Only one default setting
ag%do_clean_zeros = .true.
return
end subroutine mld_d_base_aggregator_default
@ -228,9 +251,12 @@ contains
!! will contribute, in global numbering.
!! Many aggregation produce a binary tentative prolongator, but some
!! do not, hence we also need the OP_PROL output.
!! AG_DATA is passed here just in case some of the
!! aggregators need it internally, most of them will ignore.
!!
!! \param ag The input aggregator object
!! \param parms The auxiliary parameters object
!! \param ag_data Auxiliary global aggregation info
!! \param a The local matrix part
!! \param desc_a The descriptor
!! \param ilaggr Output aggregation map
@ -239,11 +265,13 @@ contains
!! \param info Return code
!!
!
subroutine mld_d_base_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_d_base_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
implicit none
class(mld_d_base_aggregator_type), target, intent(inout) :: ag
type(mld_dml_parms), intent(inout) :: parms
type(mld_dml_parms), intent(inout) :: parms
type(mld_daggr_data), intent(in) :: ag_data
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)

38
mlprec/mld_d_dec_aggregator_mod.f90
Unescape Escape View File

@ -48,7 +48,7 @@ module mld_d_dec_aggregator_mod
!! \extends mld_d_base_aggregator_mod::mld_d_base_aggregator_type
!!
!! type, extends(mld_d_base_aggregator_type) :: mld_d_dec_aggregator_type
!! procedure(mld_d_map_bld), nopass, pointer :: map_bld => null()
!! procedure(mld_d_soc_map_bld), nopass, pointer :: soc_map_bld => null()
!! end type
!!
!! This is the simplest aggregation method: starting from the
@ -71,12 +71,12 @@ module mld_d_dec_aggregator_mod
!! PSBLAS-based parallel two-level Schwarz preconditioners, Appl. Num. Math.
!! 57 (2007), 1181-1196.
!!
!! The map_bld method is used inside the implementation of build_tprol
!! The soc_map_bld method is used inside the implementation of build_tprol
!!
!
!
type, extends(mld_d_base_aggregator_type) :: mld_d_dec_aggregator_type
procedure(mld_d_map_bld), nopass, pointer :: map_bld => null()
procedure(mld_d_soc_map_bld), nopass, pointer :: soc_map_bld => null()
contains
procedure, pass(ag) :: bld_tprol => mld_d_dec_aggregator_build_tprol
@ -88,28 +88,17 @@ module mld_d_dec_aggregator_mod
end type mld_d_dec_aggregator_type
abstract interface
subroutine mld_d_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
import :: psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_
implicit none
integer(psb_ipk_), intent(in) :: iorder
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_dpk_), intent(in) :: theta
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
integer(psb_ipk_), intent(out) :: info
end subroutine mld_d_map_bld
end interface
procedure(mld_d_map_bld) :: mld_d_soc1_map_bld, mld_d_soc2_map_bld
procedure(mld_d_soc_map_bld) :: mld_d_soc1_map_bld, mld_d_soc2_map_bld
interface
subroutine mld_d_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_d_dec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: mld_d_dec_aggregator_type, psb_desc_type, psb_dspmat_type, psb_dpk_, &
& psb_ipk_, psb_lpk_, psb_ldspmat_type, mld_dml_parms
& psb_ipk_, psb_lpk_, psb_ldspmat_type, mld_dml_parms, mld_daggr_data
implicit none
class(mld_d_dec_aggregator_type), target, intent(inout) :: ag
type(mld_dml_parms), intent(inout) :: parms
type(mld_daggr_data), intent(in) :: ag_data
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -147,14 +136,14 @@ contains
select case(parms%aggr_type)
case (mld_noalg_)
ag%map_bld => null()
ag%soc_map_bld => null()
case (mld_soc1_)
ag%map_bld => mld_d_soc1_map_bld
ag%soc_map_bld => mld_d_soc1_map_bld
case (mld_soc2_)
ag%map_bld => mld_d_soc2_map_bld
ag%soc_map_bld => mld_d_soc2_map_bld
case default
write(0,*) 'Unknown aggregation type, defaulting to SOC1'
ag%map_bld => mld_d_soc1_map_bld
ag%soc_map_bld => mld_d_soc1_map_bld
end select
return
@ -165,7 +154,8 @@ contains
implicit none
class(mld_d_dec_aggregator_type), intent(inout) :: ag
ag%map_bld => mld_d_soc1_map_bld
call ag%mld_d_base_aggregator_type%default()
ag%soc_map_bld => mld_d_soc1_map_bld
return
end subroutine mld_d_dec_aggregator_default

26
mlprec/mld_d_inner_mod.f90
Unescape Escape View File

@ -97,32 +97,6 @@ module mld_d_inner_mod
end subroutine mld_dmlprec_aply_vect
end interface mld_mlprec_aply
interface mld_aggrmap_bld
subroutine mld_d_lev_aggrmap_bld(p,a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_ldspmat_type
import :: mld_d_onelev_type
implicit none
type(mld_d_onelev_type), intent(inout), target :: p
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), intent(out) :: ilaggr(:),nlaggr(:)
type(psb_ldspmat_type), intent(out) :: op_prol
integer(psb_ipk_), intent(out) :: info
end subroutine mld_d_lev_aggrmap_bld
subroutine mld_daggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_ldspmat_type
implicit none
integer(psb_ipk_), intent(in) :: iorder
integer(psb_ipk_), intent(in) :: aggr_type
real(psb_dpk_), intent(in) :: theta
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
type(psb_ldspmat_type), intent(out) :: op_prol
integer(psb_ipk_), intent(out) :: info
end subroutine mld_daggrmap_bld
end interface mld_aggrmap_bld
interface mld_map_to_tprol
subroutine mld_d_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod, only : psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_ldspmat_type

12
mlprec/mld_d_onelev_mod.f90
Unescape Escape View File

@ -381,7 +381,7 @@ interface
interface
subroutine mld_d_base_onelev_dump(lv,level,info,prefix,head,ac,rp,smoother,&
& solver,global_num)
& solver,tprol,global_num)
import :: psb_dspmat_type, psb_d_vect_type, psb_d_base_vect_type, &
& psb_dlinmap_type, psb_dpk_, mld_d_onelev_type, &
& psb_ipk_, psb_epk_, psb_desc_type
@ -390,7 +390,7 @@ interface
integer(psb_ipk_), intent(in) :: level
integer(psb_ipk_), intent(out) :: info
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: ac, rp, smoother, solver, global_num
logical, optional, intent(in) :: ac, rp, smoother, solver, tprol, global_num
end subroutine mld_d_base_onelev_dump
end interface
@ -422,7 +422,7 @@ contains
val = val + lv%desc_ac%sizeof()
val = val + lv%ac%sizeof()
val = val + lv%tprol%sizeof()
val = val + lv%map%sizeof()
val = val + lv%map%sizeof()
if (allocated(lv%sm)) val = val + lv%sm%sizeof()
if (allocated(lv%sm2a)) val = val + lv%sm2a%sizeof()
end function d_base_onelev_sizeof
@ -484,16 +484,18 @@ contains
end subroutine d_base_onelev_default
subroutine d_base_onelev_bld_tprol(lv,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine d_base_onelev_bld_tprol(lv,a,desc_a,&
& ilaggr,nlaggr,op_prol,ag_data,info)
implicit none
class(mld_d_onelev_type), intent(inout), target :: lv
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
type(psb_ldspmat_type), intent(out) :: op_prol
type(mld_daggr_data), intent(in) :: ag_data
integer(psb_ipk_), intent(out) :: info
call lv%aggr%bld_tprol(lv%parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
call lv%aggr%bld_tprol(lv%parms,ag_data,a,desc_a,ilaggr,nlaggr,op_prol,info)
end subroutine d_base_onelev_bld_tprol

45
mlprec/mld_d_prec_type.f90
Unescape Escape View File

@ -85,18 +85,8 @@ module mld_d_prec_type
integer, parameter, private :: wv_size_=4
type, extends(psb_dprec_type) :: mld_dprec_type
! integer(psb_ipk_) :: ictxt ! Now it's in the PSBLAS prec.
!
! Aggregation defaults:
!
! 1. min_coarse_size = 0 Default target size will be computed as 40*(N_fine)**(1./3.)
integer(psb_ipk_) :: min_coarse_size = izero
! 2. maximum number of levels. Defaults to 20
integer(psb_ipk_) :: max_levs = 20_psb_ipk_
! 3. min_cr_ratio = 1.5
real(psb_dpk_) :: min_cr_ratio = 1.5_psb_dpk_
real(psb_dpk_) :: op_complexity = dzero
real(psb_dpk_) :: avg_cr = dzero
! integer(psb_ipk_) :: ictxt ! Now it's in the PSBLAS prec.
type(mld_daggr_data) :: ag_data
!
! Number of outer sweeps. Sometimes 2 V-cycles may be better than 1 W-cycle.
!
@ -445,7 +435,7 @@ contains
class(mld_dprec_type), intent(in) :: prec
real(psb_dpk_) :: val
val = prec%op_complexity
val = prec%ag_data%op_complexity
end function mld_d_get_compl
@ -480,7 +470,7 @@ contains
call psb_sum(ictxt,num)
call psb_sum(ictxt,den)
end if
prec%op_complexity = num/den
prec%ag_data%op_complexity = num/den
end subroutine mld_d_cmp_compl
!
@ -492,7 +482,7 @@ contains
class(mld_dprec_type), intent(in) :: prec
real(psb_dpk_) :: val
val = prec%avg_cr
val = prec%ag_data%avg_cr
end function mld_d_get_avg_cr
@ -517,7 +507,7 @@ contains
avgcr = avgcr / (nl-1)
end if
call psb_sum(ictxt,avgcr)
prec%avg_cr = avgcr/np
prec%ag_data%avg_cr = avgcr/np
end subroutine mld_d_cmp_avg_cr
!
@ -722,14 +712,15 @@ contains
end subroutine mld_d_apply1v
subroutine mld_d_dump(prec,info,istart,iend,prefix,head,ac,rp,smoother,solver,global_num)
subroutine mld_d_dump(prec,info,istart,iend,prefix,head,ac,rp,smoother,solver,tprol,&
& global_num)
implicit none
class(mld_dprec_type), intent(in) :: prec
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: istart, iend
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: smoother, solver,ac, rp, global_num
logical, optional, intent(in) :: smoother, solver,ac, rp, tprol, global_num
integer(psb_ipk_) :: i, j, il1, iln, lname, lev
integer(psb_ipk_) :: icontxt,iam, np
character(len=80) :: prefix_
@ -750,7 +741,8 @@ contains
do lev=il1, iln
call prec%precv(lev)%dump(lev,info,prefix=prefix,head=head,&
& ac=ac,smoother=smoother,solver=solver,rp=rp,global_num=global_num)
& ac=ac,smoother=smoother,solver=solver,rp=rp,tprol=tprol, &
& global_num=global_num)
end do
end subroutine mld_d_dump
@ -801,13 +793,9 @@ contains
info = psb_success_
select type(pout => precout)
class is (mld_dprec_type)
pout%ictxt = prec%ictxt
pout%max_levs = prec%max_levs
pout%min_coarse_size = prec%min_coarse_size
pout%min_cr_ratio = prec%min_cr_ratio
pout%outer_sweeps = prec%outer_sweeps
pout%op_complexity = prec%op_complexity
pout%avg_cr = prec%avg_cr
pout%ictxt = prec%ictxt
pout%ag_data = prec%ag_data
pout%outer_sweeps = prec%outer_sweeps
if (allocated(prec%precv)) then
ln = size(prec%precv)
allocate(pout%precv(ln),stat=info)
@ -853,7 +841,10 @@ contains
!!$ return
endif
end if
b%ictxt = prec%ictxt
b%ag_data = prec%ag_data
b%outer_sweeps = prec%outer_sweeps
call move_alloc(prec%precv,b%precv)
! Fix the pointers except on level 1.
do i=2, size(b%precv)

6
mlprec/mld_d_symdec_aggregator_mod.f90
Unescape Escape View File

@ -62,12 +62,14 @@ module mld_d_symdec_aggregator_mod
interface
subroutine mld_d_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_d_symdec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: mld_d_symdec_aggregator_type, psb_desc_type, psb_dspmat_type, psb_dpk_, &
& psb_ipk_, psb_lpk_, psb_ldspmat_type, mld_dml_parms
& psb_ipk_, psb_lpk_, psb_ldspmat_type, mld_dml_parms, mld_daggr_data
implicit none
class(mld_d_symdec_aggregator_type), target, intent(inout) :: ag
type(mld_dml_parms), intent(inout) :: parms
type(mld_daggr_data), intent(in) :: ag_data
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)

44
mlprec/mld_s_base_aggregator_mod.f90
Unescape Escape View File

@ -41,11 +41,11 @@
!
module mld_s_base_aggregator_mod
use mld_base_prec_type, only : mld_sml_parms
use mld_base_prec_type, only : mld_sml_parms, mld_saggr_data
use psb_base_mod, only : psb_sspmat_type, psb_lsspmat_type, psb_s_vect_type, &
& psb_s_base_vect_type, psb_slinmap_type, psb_spk_, &
& psb_ipk_, psb_epk_, psb_lpk_, psb_desc_type, psb_i_base_vect_type, &
& psb_erractionsave, psb_error_handler, psb_success_
& psb_erractionsave, psb_error_handler, psb_success_, psb_toupper
!
!
!
@ -79,7 +79,8 @@ module mld_s_base_aggregator_mod
!! cseti, csetr, csetc - Set internal parameters, if any
!
type mld_s_base_aggregator_type
! Do we want to purge explicit zeros when aggregating?
logical :: do_clean_zeros
contains
procedure, pass(ag) :: bld_tprol => mld_s_base_aggregator_build_tprol
procedure, pass(ag) :: mat_asb => mld_s_base_aggregator_mat_asb
@ -96,6 +97,19 @@ module mld_s_base_aggregator_mod
generic, public :: set => cseti, csetr, csetc
end type mld_s_base_aggregator_type
abstract interface
subroutine mld_s_soc_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
import :: psb_sspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_
implicit none
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_spk_), intent(in) :: theta
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
integer(psb_ipk_), intent(out) :: info
end subroutine mld_s_soc_map_bld
end interface
contains
@ -137,7 +151,16 @@ contains
character(len=*), intent(in) :: val
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: idx
! Do nothing
! Set clean zeros, or do nothing.
select case (psb_toupper(trim(what)))
case('AGGR_CLEAN_ZEROS')
select case (psb_toupper(trim(val)))
case('TRUE','T')
ag%do_clean_zeros = .true.
case('FALSE','F')
ag%do_clean_zeros = .false.
end select
end select
info = 0
end subroutine mld_s_base_aggregator_csetc
@ -181,8 +204,8 @@ contains
subroutine mld_s_base_aggregator_default(ag)
implicit none
class(mld_s_base_aggregator_type), intent(inout) :: ag
! Here we need do nothing
! Only one default setting
ag%do_clean_zeros = .true.
return
end subroutine mld_s_base_aggregator_default
@ -228,9 +251,12 @@ contains
!! will contribute, in global numbering.
!! Many aggregation produce a binary tentative prolongator, but some
!! do not, hence we also need the OP_PROL output.
!! AG_DATA is passed here just in case some of the
!! aggregators need it internally, most of them will ignore.
!!
!! \param ag The input aggregator object
!! \param parms The auxiliary parameters object
!! \param ag_data Auxiliary global aggregation info
!! \param a The local matrix part
!! \param desc_a The descriptor
!! \param ilaggr Output aggregation map
@ -239,11 +265,13 @@ contains
!! \param info Return code
!!
!
subroutine mld_s_base_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_s_base_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
implicit none
class(mld_s_base_aggregator_type), target, intent(inout) :: ag
type(mld_sml_parms), intent(inout) :: parms
type(mld_sml_parms), intent(inout) :: parms
type(mld_saggr_data), intent(in) :: ag_data
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)

38
mlprec/mld_s_dec_aggregator_mod.f90
Unescape Escape View File

@ -48,7 +48,7 @@ module mld_s_dec_aggregator_mod
!! \extends mld_s_base_aggregator_mod::mld_s_base_aggregator_type
!!
!! type, extends(mld_s_base_aggregator_type) :: mld_s_dec_aggregator_type
!! procedure(mld_s_map_bld), nopass, pointer :: map_bld => null()
!! procedure(mld_s_soc_map_bld), nopass, pointer :: soc_map_bld => null()
!! end type
!!
!! This is the simplest aggregation method: starting from the
@ -71,12 +71,12 @@ module mld_s_dec_aggregator_mod
!! PSBLAS-based parallel two-level Schwarz preconditioners, Appl. Num. Math.
!! 57 (2007), 1181-1196.
!!
!! The map_bld method is used inside the implementation of build_tprol
!! The soc_map_bld method is used inside the implementation of build_tprol
!!
!
!
type, extends(mld_s_base_aggregator_type) :: mld_s_dec_aggregator_type
procedure(mld_s_map_bld), nopass, pointer :: map_bld => null()
procedure(mld_s_soc_map_bld), nopass, pointer :: soc_map_bld => null()
contains
procedure, pass(ag) :: bld_tprol => mld_s_dec_aggregator_build_tprol
@ -88,28 +88,17 @@ module mld_s_dec_aggregator_mod
end type mld_s_dec_aggregator_type
abstract interface
subroutine mld_s_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
import :: psb_sspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_
implicit none
integer(psb_ipk_), intent(in) :: iorder
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_spk_), intent(in) :: theta
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
integer(psb_ipk_), intent(out) :: info
end subroutine mld_s_map_bld
end interface
procedure(mld_s_map_bld) :: mld_s_soc1_map_bld, mld_s_soc2_map_bld
procedure(mld_s_soc_map_bld) :: mld_s_soc1_map_bld, mld_s_soc2_map_bld
interface
subroutine mld_s_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_s_dec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: mld_s_dec_aggregator_type, psb_desc_type, psb_sspmat_type, psb_spk_, &
& psb_ipk_, psb_lpk_, psb_lsspmat_type, mld_sml_parms
& psb_ipk_, psb_lpk_, psb_lsspmat_type, mld_sml_parms, mld_saggr_data
implicit none
class(mld_s_dec_aggregator_type), target, intent(inout) :: ag
type(mld_sml_parms), intent(inout) :: parms
type(mld_saggr_data), intent(in) :: ag_data
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -147,14 +136,14 @@ contains
select case(parms%aggr_type)
case (mld_noalg_)
ag%map_bld => null()
ag%soc_map_bld => null()
case (mld_soc1_)
ag%map_bld => mld_s_soc1_map_bld
ag%soc_map_bld => mld_s_soc1_map_bld
case (mld_soc2_)
ag%map_bld => mld_s_soc2_map_bld
ag%soc_map_bld => mld_s_soc2_map_bld
case default
write(0,*) 'Unknown aggregation type, defaulting to SOC1'
ag%map_bld => mld_s_soc1_map_bld
ag%soc_map_bld => mld_s_soc1_map_bld
end select
return
@ -165,7 +154,8 @@ contains
implicit none
class(mld_s_dec_aggregator_type), intent(inout) :: ag
ag%map_bld => mld_s_soc1_map_bld
call ag%mld_s_base_aggregator_type%default()
ag%soc_map_bld => mld_s_soc1_map_bld
return
end subroutine mld_s_dec_aggregator_default

26
mlprec/mld_s_inner_mod.f90
Unescape Escape View File

@ -97,32 +97,6 @@ module mld_s_inner_mod
end subroutine mld_smlprec_aply_vect
end interface mld_mlprec_aply
interface mld_aggrmap_bld
subroutine mld_s_lev_aggrmap_bld(p,a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: psb_sspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lsspmat_type
import :: mld_s_onelev_type
implicit none
type(mld_s_onelev_type), intent(inout), target :: p
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), intent(out) :: ilaggr(:),nlaggr(:)
type(psb_lsspmat_type), intent(out) :: op_prol
integer(psb_ipk_), intent(out) :: info
end subroutine mld_s_lev_aggrmap_bld
subroutine mld_saggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: psb_sspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lsspmat_type
implicit none
integer(psb_ipk_), intent(in) :: iorder
integer(psb_ipk_), intent(in) :: aggr_type
real(psb_spk_), intent(in) :: theta
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
type(psb_lsspmat_type), intent(out) :: op_prol
integer(psb_ipk_), intent(out) :: info
end subroutine mld_saggrmap_bld
end interface mld_aggrmap_bld
interface mld_map_to_tprol
subroutine mld_s_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod, only : psb_sspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lsspmat_type

10
mlprec/mld_s_onelev_mod.f90
Unescape Escape View File

@ -381,7 +381,7 @@ interface
interface
subroutine mld_s_base_onelev_dump(lv,level,info,prefix,head,ac,rp,smoother,&
& solver,global_num)
& solver,tprol,global_num)
import :: psb_sspmat_type, psb_s_vect_type, psb_s_base_vect_type, &
& psb_slinmap_type, psb_spk_, mld_s_onelev_type, &
& psb_ipk_, psb_epk_, psb_desc_type
@ -390,7 +390,7 @@ interface
integer(psb_ipk_), intent(in) :: level
integer(psb_ipk_), intent(out) :: info
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: ac, rp, smoother, solver, global_num
logical, optional, intent(in) :: ac, rp, smoother, solver, tprol, global_num
end subroutine mld_s_base_onelev_dump
end interface
@ -484,16 +484,18 @@ contains
end subroutine s_base_onelev_default
subroutine s_base_onelev_bld_tprol(lv,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine s_base_onelev_bld_tprol(lv,a,desc_a,&
& ilaggr,nlaggr,op_prol,ag_data,info)
implicit none
class(mld_s_onelev_type), intent(inout), target :: lv
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
type(psb_lsspmat_type), intent(out) :: op_prol
type(mld_saggr_data), intent(in) :: ag_data
integer(psb_ipk_), intent(out) :: info
call lv%aggr%bld_tprol(lv%parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
call lv%aggr%bld_tprol(lv%parms,ag_data,a,desc_a,ilaggr,nlaggr,op_prol,info)
end subroutine s_base_onelev_bld_tprol

45
mlprec/mld_s_prec_type.f90
Unescape Escape View File

@ -85,18 +85,8 @@ module mld_s_prec_type
integer, parameter, private :: wv_size_=4
type, extends(psb_sprec_type) :: mld_sprec_type
! integer(psb_ipk_) :: ictxt ! Now it's in the PSBLAS prec.
!
! Aggregation defaults:
!
! 1. min_coarse_size = 0 Default target size will be computed as 40*(N_fine)**(1./3.)
integer(psb_ipk_) :: min_coarse_size = izero
! 2. maximum number of levels. Defaults to 20
integer(psb_ipk_) :: max_levs = 20_psb_ipk_
! 3. min_cr_ratio = 1.5
real(psb_spk_) :: min_cr_ratio = 1.5_psb_spk_
real(psb_spk_) :: op_complexity = szero
real(psb_spk_) :: avg_cr = szero
! integer(psb_ipk_) :: ictxt ! Now it's in the PSBLAS prec.
type(mld_saggr_data) :: ag_data
!
! Number of outer sweeps. Sometimes 2 V-cycles may be better than 1 W-cycle.
!
@ -445,7 +435,7 @@ contains
class(mld_sprec_type), intent(in) :: prec
real(psb_spk_) :: val
val = prec%op_complexity
val = prec%ag_data%op_complexity
end function mld_s_get_compl
@ -480,7 +470,7 @@ contains
call psb_sum(ictxt,num)
call psb_sum(ictxt,den)
end if
prec%op_complexity = num/den
prec%ag_data%op_complexity = num/den
end subroutine mld_s_cmp_compl
!
@ -492,7 +482,7 @@ contains
class(mld_sprec_type), intent(in) :: prec
real(psb_spk_) :: val
val = prec%avg_cr
val = prec%ag_data%avg_cr
end function mld_s_get_avg_cr
@ -517,7 +507,7 @@ contains
avgcr = avgcr / (nl-1)
end if
call psb_sum(ictxt,avgcr)
prec%avg_cr = avgcr/np
prec%ag_data%avg_cr = avgcr/np
end subroutine mld_s_cmp_avg_cr
!
@ -722,14 +712,15 @@ contains
end subroutine mld_s_apply1v
subroutine mld_s_dump(prec,info,istart,iend,prefix,head,ac,rp,smoother,solver,global_num)
subroutine mld_s_dump(prec,info,istart,iend,prefix,head,ac,rp,smoother,solver,tprol,&
& global_num)
implicit none
class(mld_sprec_type), intent(in) :: prec
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: istart, iend
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: smoother, solver,ac, rp, global_num
logical, optional, intent(in) :: smoother, solver,ac, rp, tprol, global_num
integer(psb_ipk_) :: i, j, il1, iln, lname, lev
integer(psb_ipk_) :: icontxt,iam, np
character(len=80) :: prefix_
@ -750,7 +741,8 @@ contains
do lev=il1, iln
call prec%precv(lev)%dump(lev,info,prefix=prefix,head=head,&
& ac=ac,smoother=smoother,solver=solver,rp=rp,global_num=global_num)
& ac=ac,smoother=smoother,solver=solver,rp=rp,tprol=tprol, &
& global_num=global_num)
end do
end subroutine mld_s_dump
@ -801,13 +793,9 @@ contains
info = psb_success_
select type(pout => precout)
class is (mld_sprec_type)
pout%ictxt = prec%ictxt
pout%max_levs = prec%max_levs
pout%min_coarse_size = prec%min_coarse_size
pout%min_cr_ratio = prec%min_cr_ratio
pout%outer_sweeps = prec%outer_sweeps
pout%op_complexity = prec%op_complexity
pout%avg_cr = prec%avg_cr
pout%ictxt = prec%ictxt
pout%ag_data = prec%ag_data
pout%outer_sweeps = prec%outer_sweeps
if (allocated(prec%precv)) then
ln = size(prec%precv)
allocate(pout%precv(ln),stat=info)
@ -853,7 +841,10 @@ contains
!!$ return
endif
end if
b%ictxt = prec%ictxt
b%ag_data = prec%ag_data
b%outer_sweeps = prec%outer_sweeps
call move_alloc(prec%precv,b%precv)
! Fix the pointers except on level 1.
do i=2, size(b%precv)

6
mlprec/mld_s_symdec_aggregator_mod.f90
Unescape Escape View File

@ -62,12 +62,14 @@ module mld_s_symdec_aggregator_mod
interface
subroutine mld_s_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_s_symdec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: mld_s_symdec_aggregator_type, psb_desc_type, psb_sspmat_type, psb_spk_, &
& psb_ipk_, psb_lpk_, psb_lsspmat_type, mld_sml_parms
& psb_ipk_, psb_lpk_, psb_lsspmat_type, mld_sml_parms, mld_saggr_data
implicit none
class(mld_s_symdec_aggregator_type), target, intent(inout) :: ag
type(mld_sml_parms), intent(inout) :: parms
type(mld_saggr_data), intent(in) :: ag_data
type(psb_sspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)

44
mlprec/mld_z_base_aggregator_mod.f90
Unescape Escape View File

@ -41,11 +41,11 @@
!
module mld_z_base_aggregator_mod
use mld_base_prec_type, only : mld_dml_parms
use mld_base_prec_type, only : mld_dml_parms, mld_daggr_data
use psb_base_mod, only : psb_zspmat_type, psb_lzspmat_type, psb_z_vect_type, &
& psb_z_base_vect_type, psb_zlinmap_type, psb_dpk_, &
& psb_ipk_, psb_epk_, psb_lpk_, psb_desc_type, psb_i_base_vect_type, &
& psb_erractionsave, psb_error_handler, psb_success_
& psb_erractionsave, psb_error_handler, psb_success_, psb_toupper
!
!
!
@ -79,7 +79,8 @@ module mld_z_base_aggregator_mod
!! cseti, csetr, csetc - Set internal parameters, if any
!
type mld_z_base_aggregator_type
! Do we want to purge explicit zeros when aggregating?
logical :: do_clean_zeros
contains
procedure, pass(ag) :: bld_tprol => mld_z_base_aggregator_build_tprol
procedure, pass(ag) :: mat_asb => mld_z_base_aggregator_mat_asb
@ -96,6 +97,19 @@ module mld_z_base_aggregator_mod
generic, public :: set => cseti, csetr, csetc
end type mld_z_base_aggregator_type
abstract interface
subroutine mld_z_soc_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,info)
import :: psb_zspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_
implicit none
integer(psb_ipk_), intent(in) :: iorder
logical, intent(in) :: clean_zeros
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_dpk_), intent(in) :: theta
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
integer(psb_ipk_), intent(out) :: info
end subroutine mld_z_soc_map_bld
end interface
contains
@ -137,7 +151,16 @@ contains
character(len=*), intent(in) :: val
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: idx
! Do nothing
! Set clean zeros, or do nothing.
select case (psb_toupper(trim(what)))
case('AGGR_CLEAN_ZEROS')
select case (psb_toupper(trim(val)))
case('TRUE','T')
ag%do_clean_zeros = .true.
case('FALSE','F')
ag%do_clean_zeros = .false.
end select
end select
info = 0
end subroutine mld_z_base_aggregator_csetc
@ -181,8 +204,8 @@ contains
subroutine mld_z_base_aggregator_default(ag)
implicit none
class(mld_z_base_aggregator_type), intent(inout) :: ag
! Here we need do nothing
! Only one default setting
ag%do_clean_zeros = .true.
return
end subroutine mld_z_base_aggregator_default
@ -228,9 +251,12 @@ contains
!! will contribute, in global numbering.
!! Many aggregation produce a binary tentative prolongator, but some
!! do not, hence we also need the OP_PROL output.
!! AG_DATA is passed here just in case some of the
!! aggregators need it internally, most of them will ignore.
!!
!! \param ag The input aggregator object
!! \param parms The auxiliary parameters object
!! \param ag_data Auxiliary global aggregation info
!! \param a The local matrix part
!! \param desc_a The descriptor
!! \param ilaggr Output aggregation map
@ -239,11 +265,13 @@ contains
!! \param info Return code
!!
!
subroutine mld_z_base_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_z_base_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod
implicit none
class(mld_z_base_aggregator_type), target, intent(inout) :: ag
type(mld_dml_parms), intent(inout) :: parms
type(mld_dml_parms), intent(inout) :: parms
type(mld_daggr_data), intent(in) :: ag_data
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)

38
mlprec/mld_z_dec_aggregator_mod.f90
Unescape Escape View File

@ -48,7 +48,7 @@ module mld_z_dec_aggregator_mod
!! \extends mld_z_base_aggregator_mod::mld_z_base_aggregator_type
!!
!! type, extends(mld_z_base_aggregator_type) :: mld_z_dec_aggregator_type
!! procedure(mld_z_map_bld), nopass, pointer :: map_bld => null()
!! procedure(mld_z_soc_map_bld), nopass, pointer :: soc_map_bld => null()
!! end type
!!
!! This is the simplest aggregation method: starting from the
@ -71,12 +71,12 @@ module mld_z_dec_aggregator_mod
!! PSBLAS-based parallel two-level Schwarz preconditioners, Appl. Num. Math.
!! 57 (2007), 1181-1196.
!!
!! The map_bld method is used inside the implementation of build_tprol
!! The soc_map_bld method is used inside the implementation of build_tprol
!!
!
!
type, extends(mld_z_base_aggregator_type) :: mld_z_dec_aggregator_type
procedure(mld_z_map_bld), nopass, pointer :: map_bld => null()
procedure(mld_z_soc_map_bld), nopass, pointer :: soc_map_bld => null()
contains
procedure, pass(ag) :: bld_tprol => mld_z_dec_aggregator_build_tprol
@ -88,28 +88,17 @@ module mld_z_dec_aggregator_mod
end type mld_z_dec_aggregator_type
abstract interface
subroutine mld_z_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
import :: psb_zspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_
implicit none
integer(psb_ipk_), intent(in) :: iorder
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
real(psb_dpk_), intent(in) :: theta
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
integer(psb_ipk_), intent(out) :: info
end subroutine mld_z_map_bld
end interface
procedure(mld_z_map_bld) :: mld_z_soc1_map_bld, mld_z_soc2_map_bld
procedure(mld_z_soc_map_bld) :: mld_z_soc1_map_bld, mld_z_soc2_map_bld
interface
subroutine mld_z_dec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_z_dec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: mld_z_dec_aggregator_type, psb_desc_type, psb_zspmat_type, psb_dpk_, &
& psb_ipk_, psb_lpk_, psb_lzspmat_type, mld_dml_parms
& psb_ipk_, psb_lpk_, psb_lzspmat_type, mld_dml_parms, mld_daggr_data
implicit none
class(mld_z_dec_aggregator_type), target, intent(inout) :: ag
type(mld_dml_parms), intent(inout) :: parms
type(mld_daggr_data), intent(in) :: ag_data
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
@ -147,14 +136,14 @@ contains
select case(parms%aggr_type)
case (mld_noalg_)
ag%map_bld => null()
ag%soc_map_bld => null()
case (mld_soc1_)
ag%map_bld => mld_z_soc1_map_bld
ag%soc_map_bld => mld_z_soc1_map_bld
case (mld_soc2_)
ag%map_bld => mld_z_soc2_map_bld
ag%soc_map_bld => mld_z_soc2_map_bld
case default
write(0,*) 'Unknown aggregation type, defaulting to SOC1'
ag%map_bld => mld_z_soc1_map_bld
ag%soc_map_bld => mld_z_soc1_map_bld
end select
return
@ -165,7 +154,8 @@ contains
implicit none
class(mld_z_dec_aggregator_type), intent(inout) :: ag
ag%map_bld => mld_z_soc1_map_bld
call ag%mld_z_base_aggregator_type%default()
ag%soc_map_bld => mld_z_soc1_map_bld
return
end subroutine mld_z_dec_aggregator_default

26
mlprec/mld_z_inner_mod.f90
Unescape Escape View File

@ -97,32 +97,6 @@ module mld_z_inner_mod
end subroutine mld_zmlprec_aply_vect
end interface mld_mlprec_aply
interface mld_aggrmap_bld
subroutine mld_z_lev_aggrmap_bld(p,a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: psb_zspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_lzspmat_type
import :: mld_z_onelev_type
implicit none
type(mld_z_onelev_type), intent(inout), target :: p
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), intent(out) :: ilaggr(:),nlaggr(:)
type(psb_lzspmat_type), intent(out) :: op_prol
integer(psb_ipk_), intent(out) :: info
end subroutine mld_z_lev_aggrmap_bld
subroutine mld_zaggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: psb_zspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_lzspmat_type
implicit none
integer(psb_ipk_), intent(in) :: iorder
integer(psb_ipk_), intent(in) :: aggr_type
real(psb_dpk_), intent(in) :: theta
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
type(psb_lzspmat_type), intent(out) :: op_prol
integer(psb_ipk_), intent(out) :: info
end subroutine mld_zaggrmap_bld
end interface mld_aggrmap_bld
interface mld_map_to_tprol
subroutine mld_z_map_to_tprol(desc_a,ilaggr,nlaggr,op_prol,info)
use psb_base_mod, only : psb_zspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_lzspmat_type

10
mlprec/mld_z_onelev_mod.f90
Unescape Escape View File

@ -381,7 +381,7 @@ interface
interface
subroutine mld_z_base_onelev_dump(lv,level,info,prefix,head,ac,rp,smoother,&
& solver,global_num)
& solver,tprol,global_num)
import :: psb_zspmat_type, psb_z_vect_type, psb_z_base_vect_type, &
& psb_zlinmap_type, psb_dpk_, mld_z_onelev_type, &
& psb_ipk_, psb_epk_, psb_desc_type
@ -390,7 +390,7 @@ interface
integer(psb_ipk_), intent(in) :: level
integer(psb_ipk_), intent(out) :: info
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: ac, rp, smoother, solver, global_num
logical, optional, intent(in) :: ac, rp, smoother, solver, tprol, global_num
end subroutine mld_z_base_onelev_dump
end interface
@ -484,16 +484,18 @@ contains
end subroutine z_base_onelev_default
subroutine z_base_onelev_bld_tprol(lv,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine z_base_onelev_bld_tprol(lv,a,desc_a,&
& ilaggr,nlaggr,op_prol,ag_data,info)
implicit none
class(mld_z_onelev_type), intent(inout), target :: lv
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)
type(psb_lzspmat_type), intent(out) :: op_prol
type(mld_daggr_data), intent(in) :: ag_data
integer(psb_ipk_), intent(out) :: info
call lv%aggr%bld_tprol(lv%parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
call lv%aggr%bld_tprol(lv%parms,ag_data,a,desc_a,ilaggr,nlaggr,op_prol,info)
end subroutine z_base_onelev_bld_tprol

45
mlprec/mld_z_prec_type.f90
Unescape Escape View File

@ -85,18 +85,8 @@ module mld_z_prec_type
integer, parameter, private :: wv_size_=4
type, extends(psb_zprec_type) :: mld_zprec_type
! integer(psb_ipk_) :: ictxt ! Now it's in the PSBLAS prec.
!
! Aggregation defaults:
!
! 1. min_coarse_size = 0 Default target size will be computed as 40*(N_fine)**(1./3.)
integer(psb_ipk_) :: min_coarse_size = izero
! 2. maximum number of levels. Defaults to 20
integer(psb_ipk_) :: max_levs = 20_psb_ipk_
! 3. min_cr_ratio = 1.5
real(psb_dpk_) :: min_cr_ratio = 1.5_psb_dpk_
real(psb_dpk_) :: op_complexity = dzero
real(psb_dpk_) :: avg_cr = dzero
! integer(psb_ipk_) :: ictxt ! Now it's in the PSBLAS prec.
type(mld_daggr_data) :: ag_data
!
! Number of outer sweeps. Sometimes 2 V-cycles may be better than 1 W-cycle.
!
@ -445,7 +435,7 @@ contains
class(mld_zprec_type), intent(in) :: prec
complex(psb_dpk_) :: val
val = prec%op_complexity
val = prec%ag_data%op_complexity
end function mld_z_get_compl
@ -480,7 +470,7 @@ contains
call psb_sum(ictxt,num)
call psb_sum(ictxt,den)
end if
prec%op_complexity = num/den
prec%ag_data%op_complexity = num/den
end subroutine mld_z_cmp_compl
!
@ -492,7 +482,7 @@ contains
class(mld_zprec_type), intent(in) :: prec
complex(psb_dpk_) :: val
val = prec%avg_cr
val = prec%ag_data%avg_cr
end function mld_z_get_avg_cr
@ -517,7 +507,7 @@ contains
avgcr = avgcr / (nl-1)
end if
call psb_sum(ictxt,avgcr)
prec%avg_cr = avgcr/np
prec%ag_data%avg_cr = avgcr/np
end subroutine mld_z_cmp_avg_cr
!
@ -722,14 +712,15 @@ contains
end subroutine mld_z_apply1v
subroutine mld_z_dump(prec,info,istart,iend,prefix,head,ac,rp,smoother,solver,global_num)
subroutine mld_z_dump(prec,info,istart,iend,prefix,head,ac,rp,smoother,solver,tprol,&
& global_num)
implicit none
class(mld_zprec_type), intent(in) :: prec
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(in), optional :: istart, iend
character(len=*), intent(in), optional :: prefix, head
logical, optional, intent(in) :: smoother, solver,ac, rp, global_num
logical, optional, intent(in) :: smoother, solver,ac, rp, tprol, global_num
integer(psb_ipk_) :: i, j, il1, iln, lname, lev
integer(psb_ipk_) :: icontxt,iam, np
character(len=80) :: prefix_
@ -750,7 +741,8 @@ contains
do lev=il1, iln
call prec%precv(lev)%dump(lev,info,prefix=prefix,head=head,&
& ac=ac,smoother=smoother,solver=solver,rp=rp,global_num=global_num)
& ac=ac,smoother=smoother,solver=solver,rp=rp,tprol=tprol, &
& global_num=global_num)
end do
end subroutine mld_z_dump
@ -801,13 +793,9 @@ contains
info = psb_success_
select type(pout => precout)
class is (mld_zprec_type)
pout%ictxt = prec%ictxt
pout%max_levs = prec%max_levs
pout%min_coarse_size = prec%min_coarse_size
pout%min_cr_ratio = prec%min_cr_ratio
pout%outer_sweeps = prec%outer_sweeps
pout%op_complexity = prec%op_complexity
pout%avg_cr = prec%avg_cr
pout%ictxt = prec%ictxt
pout%ag_data = prec%ag_data
pout%outer_sweeps = prec%outer_sweeps
if (allocated(prec%precv)) then
ln = size(prec%precv)
allocate(pout%precv(ln),stat=info)
@ -853,7 +841,10 @@ contains
!!$ return
endif
end if
b%ictxt = prec%ictxt
b%ag_data = prec%ag_data
b%outer_sweeps = prec%outer_sweeps
call move_alloc(prec%precv,b%precv)
! Fix the pointers except on level 1.
do i=2, size(b%precv)

6
mlprec/mld_z_symdec_aggregator_mod.f90
Unescape Escape View File

@ -62,12 +62,14 @@ module mld_z_symdec_aggregator_mod
interface
subroutine mld_z_symdec_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr,op_prol,info)
subroutine mld_z_symdec_aggregator_build_tprol(ag,parms,ag_data,&
& a,desc_a,ilaggr,nlaggr,op_prol,info)
import :: mld_z_symdec_aggregator_type, psb_desc_type, psb_zspmat_type, psb_dpk_, &
& psb_ipk_, psb_lpk_, psb_lzspmat_type, mld_dml_parms
& psb_ipk_, psb_lpk_, psb_lzspmat_type, mld_dml_parms, mld_daggr_data
implicit none
class(mld_z_symdec_aggregator_type), target, intent(inout) :: ag
type(mld_dml_parms), intent(inout) :: parms
type(mld_daggr_data), intent(in) :: ag_data
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:)

Write Preview
Loading…
Cancel
Save