fdurastante/amg4psblas
1
0
Fork 0
mirror of https://github.com/sfilippone/amg4psblas.git synced 2026-10-06 14:44:55 +00:00
Code Issues Packages Projects Releases Wiki Activity

mld2p4-2:

Browse Source 
docs/src/userinterface.tex
 mlprec/impl/level/mld_c_base_onelev_cseti.f90
 mlprec/impl/level/mld_c_base_onelev_seti.f90
 mlprec/impl/level/mld_d_base_onelev_cseti.f90
 mlprec/impl/level/mld_d_base_onelev_seti.f90
 mlprec/impl/level/mld_s_base_onelev_cseti.f90
 mlprec/impl/level/mld_s_base_onelev_seti.f90
 mlprec/impl/level/mld_z_base_onelev_cseti.f90
 mlprec/impl/level/mld_z_base_onelev_seti.f90
 mlprec/impl/mld_c_dec_map_bld.f90
 mlprec/impl/mld_caggrmap_bld.f90
 mlprec/impl/mld_ccoarse_bld.f90
 mlprec/impl/mld_ccprecset.F90
 mlprec/impl/mld_cprecinit.F90
 mlprec/impl/mld_cprecset.F90
 mlprec/impl/mld_d_dec_map_bld.f90
 mlprec/impl/mld_daggrmap_bld.f90
 mlprec/impl/mld_dcoarse_bld.f90
 mlprec/impl/mld_dcprecset.F90
 mlprec/impl/mld_dprecinit.F90
 mlprec/impl/mld_dprecset.F90
 mlprec/impl/mld_s_dec_map_bld.f90
 mlprec/impl/mld_saggrmap_bld.f90
 mlprec/impl/mld_scoarse_bld.f90
 mlprec/impl/mld_scprecset.F90
 mlprec/impl/mld_sprecinit.F90
 mlprec/impl/mld_sprecset.F90
 mlprec/impl/mld_z_dec_map_bld.f90
 mlprec/impl/mld_zaggrmap_bld.f90
 mlprec/impl/mld_zcoarse_bld.f90
 mlprec/impl/mld_zcprecset.F90
 mlprec/impl/mld_zprecinit.F90
 mlprec/impl/mld_zprecset.F90
 mlprec/mld_base_prec_type.F90
 mlprec/mld_c_inner_mod.f90
 mlprec/mld_c_onelev_mod.f90
 mlprec/mld_d_base_solver_mod.f90
 mlprec/mld_d_inner_mod.f90
 mlprec/mld_d_onelev_mod.f90
 mlprec/mld_d_prec_type.f90
 mlprec/mld_d_umf_solver.F90
 mlprec/mld_s_inner_mod.f90
 mlprec/mld_s_onelev_mod.f90
 mlprec/mld_z_inner_mod.f90
 mlprec/mld_z_onelev_mod.f90
 tests/pdegen/ppde2d.f90
 tests/pdegen/ppde3d.f90
 tests/pdegen/runs/ppde.inp
 tests/pdegen/spde2d.f90
 tests/pdegen/spde3d.f90

Added option to apply ordering while aggregating.
This commit is contained in:
Salvatore Filippone
2016-05-09 16:00:12 +00:00
parent c2ee7bd167
commit 875443efe7
50 changed files with 489 additions and 296 deletions
Download Patch File Download Diff File Expand all files Collapse all files
docs/src/userinterface.tex
+8
View File
@@ -305,6 +305,14 @@ according to their needs.
& \texttt{'DEC'}
& Aggregation algorithm. Currently, only the
decoupled aggregation is available. \\ \hline
\verb|mld_aggr_ord_| \break \verb|AGGR_ORD| & \verb|character(len=*)|
& \texttt{'NAT'}
& \texttt{'DEGREE'}
& Initial ordering of indices for aggregation
algorithm: natural ordering or sorted by
descending degree of the node in the
matrix graph. Since aggregation is
heuristics, results will be different. \\ \hline
\verb|mld_aggr_kind_| \break \verb|AGGR_KIND| & \verb|character(len=*)|
& \texttt{'SMOOTHED'} \hspace{2.5cm} \texttt{'NONSMOOTHED'}
& \texttt{'SMOOTHED'}
mlprec/impl/level/mld_c_base_onelev_cseti.f90
+3
View File
@@ -73,6 +73,9 @@ subroutine mld_c_base_onelev_cseti(lv,what,val,info)
case ('AGGR_ALG')
lv%parms%aggr_alg = val
case ('AGGR_ORD')
lv%parms%aggr_ord = val
case ('AGGR_KIND')
lv%parms%aggr_kind = val
mlprec/impl/level/mld_c_base_onelev_seti.f90
+3
View File
@@ -73,6 +73,9 @@ subroutine mld_c_base_onelev_seti(lv,what,val,info)
case (mld_aggr_alg_)
lv%parms%aggr_alg = val
case (mld_aggr_ord_)
lv%parms%aggr_ord = val
case (mld_aggr_kind_)
lv%parms%aggr_kind = val
mlprec/impl/level/mld_d_base_onelev_cseti.f90
+3
View File
@@ -73,6 +73,9 @@ subroutine mld_d_base_onelev_cseti(lv,what,val,info)
case ('AGGR_ALG')
lv%parms%aggr_alg = val
case ('AGGR_ORD')
lv%parms%aggr_ord = val
case ('AGGR_KIND')
lv%parms%aggr_kind = val
mlprec/impl/level/mld_d_base_onelev_seti.f90
+3
View File
@@ -73,6 +73,9 @@ subroutine mld_d_base_onelev_seti(lv,what,val,info)
case (mld_aggr_alg_)
lv%parms%aggr_alg = val
case (mld_aggr_ord_)
lv%parms%aggr_ord = val
case (mld_aggr_kind_)
lv%parms%aggr_kind = val
mlprec/impl/level/mld_s_base_onelev_cseti.f90
+3
View File
@@ -73,6 +73,9 @@ subroutine mld_s_base_onelev_cseti(lv,what,val,info)
case ('AGGR_ALG')
lv%parms%aggr_alg = val
case ('AGGR_ORD')
lv%parms%aggr_ord = val
case ('AGGR_KIND')
lv%parms%aggr_kind = val
mlprec/impl/level/mld_s_base_onelev_seti.f90
+3
View File
@@ -73,6 +73,9 @@ subroutine mld_s_base_onelev_seti(lv,what,val,info)
case (mld_aggr_alg_)
lv%parms%aggr_alg = val
case (mld_aggr_ord_)
lv%parms%aggr_ord = val
case (mld_aggr_kind_)
lv%parms%aggr_kind = val
mlprec/impl/level/mld_z_base_onelev_cseti.f90
+3
View File
@@ -73,6 +73,9 @@ subroutine mld_z_base_onelev_cseti(lv,what,val,info)
case ('AGGR_ALG')
lv%parms%aggr_alg = val
case ('AGGR_ORD')
lv%parms%aggr_ord = val
case ('AGGR_KIND')
lv%parms%aggr_kind = val
mlprec/impl/level/mld_z_base_onelev_seti.f90
+3
View File
@@ -73,6 +73,9 @@ subroutine mld_z_base_onelev_seti(lv,what,val,info)
case (mld_aggr_alg_)
lv%parms%aggr_alg = val
case (mld_aggr_ord_)
lv%parms%aggr_ord = val
case (mld_aggr_kind_)
lv%parms%aggr_kind = val
mlprec/impl/mld_c_dec_map_bld.f90
+33 -18
View File
@@ -37,7 +37,7 @@
!!$
!!$
subroutine mld_c_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_c_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_c_inner_mod, mld_protect_name => mld_c_dec_map_bld
@@ -45,6 +45,7 @@ subroutine mld_c_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
implicit none
! Arguments
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
@@ -52,9 +53,11 @@ subroutine mld_c_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
integer(psb_ipk_), intent(out) :: info
! Local variables
integer(psb_ipk_), allocatable :: ils(:), neigh(:), irow(:), icol(:)
integer(psb_ipk_), allocatable :: ils(:), neigh(:), irow(:), icol(:),&
& ideg(:), idxs(:)
complex(psb_spk_), allocatable :: val(:), diag(:)
integer(psb_ipk_) :: icnt,nlp,k,n,ia,isz,nr, naggr,i,j,m, nz, ilg
integer(psb_ipk_) :: icnt,nlp,k,n,ia,isz,nr, naggr,i,j,m, nz, ilg, ii
type(psb_c_csr_sparse_mat) :: acsr
real(psb_spk_) :: cpling, tcl
logical :: recovery
integer(psb_ipk_) :: debug_level, debug_unit,err_act
@@ -75,7 +78,7 @@ subroutine mld_c_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
ncol = desc_a%get_local_cols()
nr = a%get_nrows()
allocate(ilaggr(nr),neigh(nr),stat=info)
allocate(ilaggr(nr),neigh(nr),ideg(nr),idxs(nr),stat=info)
if(info /= psb_success_) then
info=psb_err_alloc_request_
call psb_errpush(info,name,i_err=(/2*nr,izero,izero,izero,izero/),&
@@ -90,11 +93,20 @@ subroutine mld_c_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
goto 9999
end if
do i=1, nr
ilaggr(i) = -(nr+1)
end do
if (iorder == mld_aggr_ord_nat_) then
do i=1, nr
ilaggr(i) = -(nr+1)
idxs(i) = i
end do
else
call a%cp_to(acsr)
do i=1, nr
ilaggr(i) = -(nr+1)
ideg(i) = acsr%irp(i+1) - acsr%irp(i)
end do
call acsr%free()
call psb_msort(ideg,ix=idxs,dir=psb_sort_down_)
end if
! Note: -(nr+1) Untouched as yet
! -i 1<=i<=nr Adjacent to aggregate i
! i 1<=i<=nr Belonging to aggregate i
@@ -106,7 +118,8 @@ subroutine mld_c_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
nlp = 0
do
icnt = 0
do i=1, nr
do ii=1, nr
i = idxs(ii)
if (ilaggr(i) == -(nr+1)) then
!
! 1. Untouched nodes are marked >0 together
@@ -124,14 +137,16 @@ subroutine mld_c_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
end if
do k=1, nz
j = icol(k)
ilg = ilaggr(j)
if ((ilg<0).and.(1<=j).and.(j<=nr).and.(i /= j)) then
if (abs(val(k)) > theta*sqrt(abs(diag(i)*diag(j)))) then
ilaggr(j) = naggr
else
ilaggr(j) = -naggr
endif
j = icol(k)
if ((1<=j).and.(j<=nr)) then
ilg = ilaggr(j)
if ((ilg<0).and.(i /= j)) then
if (abs(val(k)) > theta*sqrt(abs(diag(i)*diag(j)))) then
ilaggr(j) = naggr
else
ilaggr(j) = -naggr
endif
end if
end if
enddo
mlprec/impl/mld_caggrmap_bld.f90
+5 -4
View File
@@ -79,7 +79,7 @@
! info - integer, output.
! Error code.
!
subroutine mld_caggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
subroutine mld_caggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,info)
use psb_base_mod
use mld_c_inner_mod, mld_protect_name => mld_caggrmap_bld
@@ -87,6 +87,7 @@ subroutine mld_caggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
implicit none
! Arguments
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
@@ -118,7 +119,7 @@ subroutine mld_caggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
select case (aggr_type)
case (mld_dec_aggr_)
call mld_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
call mld_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
case (mld_sym_dec_aggr_)
nr = a%get_nrows()
@@ -126,14 +127,14 @@ subroutine mld_caggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
& rscale=.false.,cscale=.false.)
call atmp%set_nrows(nr)
call atmp%set_ncols(nr)
if (info == psb_success_) call atrans%transp(atmp)
if (info == psb_success_) call atmp%transp(atrans)
if (info == psb_success_) call atrans%cscnv(info,type='COO')
if (info == psb_success_) call psb_rwextd(nr,atmp,info,b=atrans,rowscale=.false.)
call atmp%set_nrows(nr)
call atmp%set_ncols(nr)
if (info == psb_success_) call atrans%free()
if (info == psb_success_) call atmp%cscnv(info,type='CSR')
if (info == psb_success_) call mld_dec_map_bld(theta,atmp,desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call mld_dec_map_bld(iorder,theta,atmp,desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call atmp%free()
case default
mlprec/impl/mld_ccoarse_bld.f90
+5 -2
View File
@@ -39,7 +39,7 @@
! File: mld_ccoarse_bld.f90
!
! Subroutine: mld_ccoarse_bld
! Version: complex
! Version: real
!
! This routine builds the matrix associated to the current level of the
! multilevel preconditioner from the matrix associated to the previous level,
@@ -95,6 +95,8 @@ subroutine mld_ccoarse_bld(a,desc_a,p,info)
& mld_mult_ml_,is_legal_ml_type)
call mld_check_def(p%parms%aggr_alg,'Aggregation',&
& mld_dec_aggr_,is_legal_ml_aggr_alg)
call mld_check_def(p%parms%aggr_ord,'Ordering',&
& mld_aggr_ord_nat_,is_legal_ml_aggr_ord)
call mld_check_def(p%parms%aggr_kind,'Smoother',&
& mld_smooth_prol_,is_legal_ml_aggr_kind)
call mld_check_def(p%parms%coarse_mat,'Coarse matrix',&
@@ -117,7 +119,7 @@ subroutine mld_ccoarse_bld(a,desc_a,p,info)
! aggregation algorithm. This also defines a tentative prolongator from
! the coarse to the fine level.
!
call mld_aggrmap_bld(p%parms%aggr_alg,p%parms%aggr_thresh,&
call mld_aggrmap_bld(p%parms%aggr_alg,p%parms%aggr_ord,p%parms%aggr_thresh,&
& a,desc_a,ilaggr,nlaggr,info)
if (info /= psb_success_) then
@@ -136,6 +138,7 @@ subroutine mld_ccoarse_bld(a,desc_a,p,info)
call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_aggrmat_asb')
goto 9999
end if
!
! Fix the base_a and base_desc pointers for handling of residuals.
! This is correct because this routine is only called at levels >=2.
mlprec/impl/mld_ccprecset.F90
+10 -9
View File
@@ -93,6 +93,7 @@ subroutine mld_ccprecseti(p,what,val,info,ilev)
use mld_c_mumps_solver
#endif
implicit none
! Arguments
@@ -148,9 +149,10 @@ subroutine mld_ccprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case('SUB_SOLVE')
call onelev_set_solver(p%precv(ilev_),val,info)
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_KIND',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG','AGGR_EIG',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_ORD',&
& 'AGGR_KIND','SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','SMOOTHER_SWEEPS_PRE',&
& 'SMOOTHER_SWEEPS_POST',&
& 'SUB_RESTR','SUB_PROL', &
& 'SUB_REN','SUB_OVR','SUB_FILLIN')
call p%precv(ilev_)%set(what,val,info)
@@ -166,9 +168,10 @@ subroutine mld_ccprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case('SUB_SOLVE')
call onelev_set_solver(p%precv(ilev_),val,info)
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_KIND',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG','AGGR_EIG',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_ORD',&
& 'AGGR_KIND','SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','SMOOTHER_SWEEPS_PRE',&
& 'SMOOTHER_SWEEPS_POST',&
& 'SUB_RESTR','SUB_PROL', &
& 'SUB_REN','SUB_OVR','SUB_FILLIN',&
& 'COARSE_MAT')
@@ -280,7 +283,7 @@ subroutine mld_ccprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
end do
case('ML_TYPE','AGGR_ALG','AGGR_KIND',&
case('ML_TYPE','AGGR_ALG','AGGR_ORD','AGGR_KIND',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','AGGR_FILTER')
@@ -302,7 +305,6 @@ subroutine mld_ccprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
#if defined(HAVE_SLU_)
call onelev_set_solver(p%precv(nlev_),mld_slu_,info)
call onelev_set_solver(p%precv(nlev_),mld_slu_,info)
#elif defined(HAVE_MUMPS_)
call onelev_set_solver(p%precv(nlev_),mld_mumps_,info)
#else
@@ -321,7 +323,6 @@ subroutine mld_ccprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
call onelev_set_solver(p%precv(nlev_),val,info)
call p%precv(nlev_)%set('COARSE_MAT',mld_distr_mat_,info)
case(mld_jac_)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
call onelev_set_solver(p%precv(nlev_),mld_diag_scale_,info)
mlprec/impl/mld_cprecinit.F90
+8 -12
View File
@@ -2,9 +2,9 @@
!!$
!!$ MLD2P4 version 2.0
!!$ MultiLevel Domain Decomposition Parallel Preconditioners Package
!!$ based on PSBLAS (Parallel Sparse BLAS version 3.0)
!!$ based on PSBLAS (Parallel Sparse BLAS version 3.3)
!!$
!!$ (C) Copyright 2008,2009,2010,2012,2013
!!$ (C) Copyright 2008, 2010, 2012, 2015
!!$
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari CNRS-IRIT, Toulouse
@@ -61,7 +61,7 @@
!
! 'ML' - Multilevel hybrid preconditioner (additive on the
! same level and multiplicative through the levels),
! with 2 levels and post-smoothing only. RAS with
! with 2 levels, pre and post-smoothing, RAS with
! overlap 1 and ILU(0) on the local blocks is
! applied as post-smoother at each level, but the
! coarsest one; four sweeps of the block-Jacobi solver,
@@ -97,9 +97,6 @@ subroutine mld_cprecinit(p,ptype,info,nlev)
use mld_c_id_solver
use mld_c_diag_solver
use mld_c_ilu_solver
#if defined(HAVE_UMF_) && 0
use mld_c_umf_solver
#endif
#if defined(HAVE_SLU_)
use mld_c_slu_solver
#endif
@@ -115,7 +112,7 @@ subroutine mld_cprecinit(p,ptype,info,nlev)
! Local variables
integer(psb_ipk_) :: nlev_, ilev_
real(psb_spk_) :: thr
real(psb_spk_) :: thr, scale
character(len=*), parameter :: name='mld_precinit'
info = psb_success_
@@ -191,9 +188,7 @@ subroutine mld_cprecinit(p,ptype,info,nlev)
ilev_ = nlev_
allocate(mld_c_jac_smoother_type :: p%precv(ilev_)%sm, stat=info)
if (info /= psb_success_) return
#if defined(HAVE_UMF_) && 0
allocate(mld_c_umf_solver_type :: p%precv(ilev_)%sm%sv, stat=info)
#elif defined(HAVE_SLU_)
#if defined(HAVE_SLU_)
allocate(mld_c_slu_solver_type :: p%precv(ilev_)%sm%sv, stat=info)
#else
allocate(mld_c_ilu_solver_type :: p%precv(ilev_)%sm%sv, stat=info)
@@ -205,10 +200,11 @@ subroutine mld_cprecinit(p,ptype,info,nlev)
call p%precv(ilev_)%set(mld_sub_prol_,psb_none_,info)
call p%precv(ilev_)%set(mld_sub_ovr_,izero,info)
thr = 0.16d0
thr = 0.05
scale = 1.0
do ilev_=1,nlev_
call p%precv(ilev_)%set(mld_aggr_thresh_,thr,info)
thr = thr/2
call p%precv(ilev_)%set(mld_aggr_scale_,scale,info)
end do
case default
mlprec/impl/mld_cprecset.F90
+26 -13
View File
@@ -92,6 +92,7 @@ subroutine mld_cprecseti(p,what,val,info,ilev)
#if defined(HAVE_MUMPS_)
use mld_c_mumps_solver
#endif
implicit none
! Arguments
@@ -148,8 +149,8 @@ subroutine mld_cprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case(mld_sub_solve_)
call onelev_set_solver(p%precv(ilev_),val,info)
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,&
& mld_aggr_kind_,mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_sub_restr_,mld_sub_prol_, &
& mld_sub_ren_,mld_sub_ovr_,mld_sub_fillin_)
@@ -166,8 +167,8 @@ subroutine mld_cprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case(mld_sub_solve_)
call onelev_set_solver(p%precv(ilev_),val,info)
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,&
& mld_aggr_kind_,mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_sub_restr_,mld_sub_prol_, &
& mld_sub_ren_,mld_sub_ovr_,mld_sub_fillin_,&
@@ -207,7 +208,7 @@ subroutine mld_cprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
call onelev_set_solver(p%precv(nlev_),val,info)
call p%precv(nlev_)%set(mld_coarse_mat_,mld_repl_mat_,info)
case(mld_sludist_, mld_mumps_)
case(mld_sludist_,mld_mumps_)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
call onelev_set_solver(p%precv(nlev_),val,info)
call p%precv(nlev_)%set(mld_coarse_mat_,mld_distr_mat_,info)
@@ -276,7 +277,7 @@ subroutine mld_cprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
end do
case(mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
case(mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,mld_aggr_kind_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,&
& mld_aggr_eig_,mld_aggr_filter_)
@@ -298,8 +299,8 @@ subroutine mld_cprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
#if defined(HAVE_SLU_)
call onelev_set_solver(p%precv(nlev_),mld_slu_,info)
#elif defined(HAVE_MUMPS_)
call onelev_set_solver(p%precv(nlev_),mld_mumps_,info)
#elif defined(HAVE_SLU_)
call onelev_set_solver(p%precv(nlev_),mld_slu_,info)
#else
call onelev_set_solver(p%precv(nlev_),mld_ilu_n_,info)
#endif
@@ -599,7 +600,6 @@ contains
end if
end if
#endif
case default
!
! Do nothing and hope for the best :)
@@ -723,13 +723,26 @@ subroutine mld_cprecsetsv(p,val,info,ilev)
do ilev_ = ilmin, ilmax
if (allocated(p%precv(ilev_)%sm)) then
if (allocated(p%precv(ilev_)%sm%sv)) &
& deallocate(p%precv(ilev_)%sm%sv)
if (allocated(p%precv(ilev_)%sm%sv)) then
if (.not.same_type_as(p%precv(ilev_)%sm%sv,val)) then
deallocate(p%precv(ilev_)%sm%sv,stat=info)
if (info /= 0) then
info = 3111
return
end if
end if
if (.not.allocated(p%precv(ilev_)%sm%sv)) then
#ifdef HAVE_MOLD
allocate(p%precv(ilev_)%sm%sv,mold=val)
allocate(p%precv(ilev_)%sm%sv,mold=val,stat=info)
#else
allocate(p%precv(ilev_)%sm%sv,source=val)
allocate(p%precv(ilev_)%sm%sv,source=val,stat=info)
#endif
if (info /= 0) then
info = 3111
return
end if
end if
end if
call p%precv(ilev_)%sm%sv%default()
else
info = 3111
mlprec/impl/mld_d_dec_map_bld.f90
+33 -18
View File
@@ -37,7 +37,7 @@
!!$
!!$
subroutine mld_d_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_d_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_d_inner_mod, mld_protect_name => mld_d_dec_map_bld
@@ -45,6 +45,7 @@ subroutine mld_d_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
implicit none
! Arguments
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
@@ -52,9 +53,11 @@ subroutine mld_d_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
integer(psb_ipk_), intent(out) :: info
! Local variables
integer(psb_ipk_), allocatable :: ils(:), neigh(:), irow(:), icol(:)
integer(psb_ipk_), allocatable :: ils(:), neigh(:), irow(:), icol(:),&
& ideg(:), idxs(:)
real(psb_dpk_), allocatable :: val(:), diag(:)
integer(psb_ipk_) :: icnt,nlp,k,n,ia,isz,nr, naggr,i,j,m, nz, ilg
integer(psb_ipk_) :: icnt,nlp,k,n,ia,isz,nr, naggr,i,j,m, nz, ilg, ii
type(psb_d_csr_sparse_mat) :: acsr
real(psb_dpk_) :: cpling, tcl
logical :: recovery
integer(psb_ipk_) :: debug_level, debug_unit,err_act
@@ -75,7 +78,7 @@ subroutine mld_d_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
ncol = desc_a%get_local_cols()
nr = a%get_nrows()
allocate(ilaggr(nr),neigh(nr),stat=info)
allocate(ilaggr(nr),neigh(nr),ideg(nr),idxs(nr),stat=info)
if(info /= psb_success_) then
info=psb_err_alloc_request_
call psb_errpush(info,name,i_err=(/2*nr,izero,izero,izero,izero/),&
@@ -90,11 +93,20 @@ subroutine mld_d_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
goto 9999
end if
do i=1, nr
ilaggr(i) = -(nr+1)
end do
if (iorder == mld_aggr_ord_nat_) then
do i=1, nr
ilaggr(i) = -(nr+1)
idxs(i) = i
end do
else
call a%cp_to(acsr)
do i=1, nr
ilaggr(i) = -(nr+1)
ideg(i) = acsr%irp(i+1) - acsr%irp(i)
end do
call acsr%free()
call psb_msort(ideg,ix=idxs,dir=psb_sort_down_)
end if
! Note: -(nr+1) Untouched as yet
! -i 1<=i<=nr Adjacent to aggregate i
! i 1<=i<=nr Belonging to aggregate i
@@ -106,7 +118,8 @@ subroutine mld_d_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
nlp = 0
do
icnt = 0
do i=1, nr
do ii=1, nr
i = idxs(ii)
if (ilaggr(i) == -(nr+1)) then
!
! 1. Untouched nodes are marked >0 together
@@ -124,14 +137,16 @@ subroutine mld_d_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
end if
do k=1, nz
j = icol(k)
ilg = ilaggr(j)
if ((ilg<0).and.(1<=j).and.(j<=nr).and.(i /= j)) then
if (abs(val(k)) > theta*sqrt(abs(diag(i)*diag(j)))) then
ilaggr(j) = naggr
else
ilaggr(j) = -naggr
endif
j = icol(k)
if ((1<=j).and.(j<=nr)) then
ilg = ilaggr(j)
if ((ilg<0).and.(i /= j)) then
if (abs(val(k)) > theta*sqrt(abs(diag(i)*diag(j)))) then
ilaggr(j) = naggr
else
ilaggr(j) = -naggr
endif
end if
end if
enddo
mlprec/impl/mld_daggrmap_bld.f90
+5 -4
View File
@@ -79,7 +79,7 @@
! info - integer, output.
! Error code.
!
subroutine mld_daggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
subroutine mld_daggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,info)
use psb_base_mod
use mld_d_inner_mod, mld_protect_name => mld_daggrmap_bld
@@ -87,6 +87,7 @@ subroutine mld_daggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
implicit none
! Arguments
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
@@ -118,7 +119,7 @@ subroutine mld_daggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
select case (aggr_type)
case (mld_dec_aggr_)
call mld_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
call mld_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
case (mld_sym_dec_aggr_)
nr = a%get_nrows()
@@ -126,14 +127,14 @@ subroutine mld_daggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
& rscale=.false.,cscale=.false.)
call atmp%set_nrows(nr)
call atmp%set_ncols(nr)
if (info == psb_success_) call atrans%transp(atmp)
if (info == psb_success_) call atmp%transp(atrans)
if (info == psb_success_) call atrans%cscnv(info,type='COO')
if (info == psb_success_) call psb_rwextd(nr,atmp,info,b=atrans,rowscale=.false.)
call atmp%set_nrows(nr)
call atmp%set_ncols(nr)
if (info == psb_success_) call atrans%free()
if (info == psb_success_) call atmp%cscnv(info,type='CSR')
if (info == psb_success_) call mld_dec_map_bld(theta,atmp,desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call mld_dec_map_bld(iorder,theta,atmp,desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call atmp%free()
case default
mlprec/impl/mld_dcoarse_bld.f90
+3 -1
View File
@@ -95,6 +95,8 @@ subroutine mld_dcoarse_bld(a,desc_a,p,info)
& mld_mult_ml_,is_legal_ml_type)
call mld_check_def(p%parms%aggr_alg,'Aggregation',&
& mld_dec_aggr_,is_legal_ml_aggr_alg)
call mld_check_def(p%parms%aggr_ord,'Ordering',&
& mld_aggr_ord_nat_,is_legal_ml_aggr_ord)
call mld_check_def(p%parms%aggr_kind,'Smoother',&
& mld_smooth_prol_,is_legal_ml_aggr_kind)
call mld_check_def(p%parms%coarse_mat,'Coarse matrix',&
@@ -117,7 +119,7 @@ subroutine mld_dcoarse_bld(a,desc_a,p,info)
! aggregation algorithm. This also defines a tentative prolongator from
! the coarse to the fine level.
!
call mld_aggrmap_bld(p%parms%aggr_alg,p%parms%aggr_thresh,&
call mld_aggrmap_bld(p%parms%aggr_alg,p%parms%aggr_ord,p%parms%aggr_thresh,&
& a,desc_a,ilaggr,nlaggr,info)
if (info /= psb_success_) then
mlprec/impl/mld_dcprecset.F90
+32 -28
View File
@@ -99,6 +99,7 @@ subroutine mld_dcprecseti(p,what,val,info,ilev)
use mld_d_mumps_solver
#endif
implicit none
! Arguments
@@ -113,6 +114,7 @@ subroutine mld_dcprecseti(p,what,val,info,ilev)
character(len=*), parameter :: name='mld_precseti'
info = psb_success_
if (.not.allocated(p%precv)) then
info = 3111
write(psb_err_unit,*) name,&
@@ -153,9 +155,10 @@ subroutine mld_dcprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case('SUB_SOLVE')
call onelev_set_solver(p%precv(ilev_),val,info)
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_KIND',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG','AGGR_EIG',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_ORD',&
& 'AGGR_KIND','SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','SMOOTHER_SWEEPS_PRE',&
& 'SMOOTHER_SWEEPS_POST',&
& 'SUB_RESTR','SUB_PROL', &
& 'SUB_REN','SUB_OVR','SUB_FILLIN')
call p%precv(ilev_)%set(what,val,info)
@@ -171,9 +174,10 @@ subroutine mld_dcprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case('SUB_SOLVE')
call onelev_set_solver(p%precv(ilev_),val,info)
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_KIND',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG','AGGR_EIG',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_ORD',&
& 'AGGR_KIND','SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','SMOOTHER_SWEEPS_PRE',&
& 'SMOOTHER_SWEEPS_POST',&
& 'SUB_RESTR','SUB_PROL', &
& 'SUB_REN','SUB_OVR','SUB_FILLIN',&
& 'COARSE_MAT')
@@ -287,7 +291,7 @@ subroutine mld_dcprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
end do
case('ML_TYPE','AGGR_ALG','AGGR_KIND',&
case('ML_TYPE','AGGR_ALG','AGGR_ORD','AGGR_KIND',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','AGGR_FILTER')
@@ -592,6 +596,27 @@ contains
info = -5
end if
#endif
#ifdef HAVE_MUMPS_
case (mld_mumps_)
if (allocated(level%sm%sv)) then
select type (sv => level%sm%sv)
class is (mld_d_mumps_solver_type)
! do nothing
class default
call level%sm%sv%free(info)
if (info == 0) deallocate(level%sm%sv)
if (info == 0) allocate(mld_d_mumps_solver_type ::&
& level%sm%sv, stat=info)
end select
else
allocate(mld_d_mumps_solver_type :: level%sm%sv, stat=info)
endif
if (allocated(level%sm)) then
if (allocated(level%sm%sv)) &
& call level%sm%sv%default()
end if
#endif
#ifdef HAVE_UMF_
case (mld_umf_)
if (allocated(level%sm)) then
@@ -642,27 +667,6 @@ contains
info = -5
end if
#endif
#ifdef HAVE_MUMPS_
case (mld_mumps_)
if (allocated(level%sm%sv)) then
select type (sv => level%sm%sv)
class is (mld_d_mumps_solver_type)
! do nothing
class default
call level%sm%sv%free(info)
if (info == 0) deallocate(level%sm%sv)
if (info == 0) allocate(mld_d_mumps_solver_type ::&
& level%sm%sv, stat=info)
end select
else
allocate(mld_d_mumps_solver_type :: level%sm%sv, stat=info)
endif
if (allocated(level%sm)) then
if (allocated(level%sm%sv)) &
& call level%sm%sv%default()
end if
#endif
case default
!
! Do nothing and hope for the best :)
mlprec/impl/mld_dprecinit.F90
+8 -6
View File
@@ -2,9 +2,9 @@
!!$
!!$ MLD2P4 version 2.0
!!$ MultiLevel Domain Decomposition Parallel Preconditioners Package
!!$ based on PSBLAS (Parallel Sparse BLAS version 3.0)
!!$ based on PSBLAS (Parallel Sparse BLAS version 3.3)
!!$
!!$ (C) Copyright 2008,2009,2010,2012,2013
!!$ (C) Copyright 2008, 2010, 2012, 2015
!!$
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari CNRS-IRIT, Toulouse
@@ -61,7 +61,7 @@
!
! 'ML' - Multilevel hybrid preconditioner (additive on the
! same level and multiplicative through the levels),
! with 2 levels and post-smoothing only. RAS with
! with 2 levels, pre and post-smoothing, RAS with
! overlap 1 and ILU(0) on the local blocks is
! applied as post-smoother at each level, but the
! coarsest one; four sweeps of the block-Jacobi solver,
@@ -104,6 +104,7 @@ subroutine mld_dprecinit(p,ptype,info,nlev)
use mld_d_slu_solver
#endif
implicit none
! Arguments
@@ -114,7 +115,7 @@ subroutine mld_dprecinit(p,ptype,info,nlev)
! Local variables
integer(psb_ipk_) :: nlev_, ilev_
real(psb_dpk_) :: thr
real(psb_dpk_) :: thr, scale
character(len=*), parameter :: name='mld_precinit'
info = psb_success_
@@ -204,10 +205,11 @@ subroutine mld_dprecinit(p,ptype,info,nlev)
call p%precv(ilev_)%set(mld_sub_prol_,psb_none_,info)
call p%precv(ilev_)%set(mld_sub_ovr_,izero,info)
thr = 0.16d0
thr = 0.05
scale = 1.0
do ilev_=1,nlev_
call p%precv(ilev_)%set(mld_aggr_thresh_,thr,info)
thr = thr/2
call p%precv(ilev_)%set(mld_aggr_scale_,scale,info)
end do
case default
mlprec/impl/mld_dprecset.F90
+6 -7
View File
@@ -155,8 +155,8 @@ subroutine mld_dprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case(mld_sub_solve_)
call onelev_set_solver(p%precv(ilev_),val,info)
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,&
& mld_aggr_kind_,mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_sub_restr_,mld_sub_prol_, &
& mld_sub_ren_,mld_sub_ovr_,mld_sub_fillin_)
@@ -173,8 +173,8 @@ subroutine mld_dprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case(mld_sub_solve_)
call onelev_set_solver(p%precv(ilev_),val,info)
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,&
& mld_aggr_kind_,mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_sub_restr_,mld_sub_prol_, &
& mld_sub_ren_,mld_sub_ovr_,mld_sub_fillin_,&
@@ -285,7 +285,7 @@ subroutine mld_dprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
end do
case(mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
case(mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,mld_aggr_kind_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,&
& mld_aggr_eig_,mld_aggr_filter_)
@@ -309,8 +309,6 @@ subroutine mld_dprecseti(p,what,val,info,ilev)
call onelev_set_solver(p%precv(nlev_),mld_umf_,info)
#elif defined(HAVE_SLU_)
call onelev_set_solver(p%precv(nlev_),mld_slu_,info)
#elif defined(HAVE_MUMPS_)
call onelev_set_solver(p%precv(nlev_),mld_mumps_,info)
#else
call onelev_set_solver(p%precv(nlev_),mld_ilu_n_,info)
#endif
@@ -1002,6 +1000,7 @@ subroutine mld_dprecsetr(p,what,val,info,ilev)
case(mld_coarse_iluthrs_)
ilev_=nlev_
call p%precv(ilev_)%set(mld_sub_iluthrs_,val,info)
case(mld_aggr_thresh_)
thr = val
do ilev_ = 2, nlev_
mlprec/impl/mld_s_dec_map_bld.f90
+33 -18
View File
@@ -37,7 +37,7 @@
!!$
!!$
subroutine mld_s_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_s_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_s_inner_mod, mld_protect_name => mld_s_dec_map_bld
@@ -45,6 +45,7 @@ subroutine mld_s_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
implicit none
! Arguments
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
@@ -52,9 +53,11 @@ subroutine mld_s_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
integer(psb_ipk_), intent(out) :: info
! Local variables
integer(psb_ipk_), allocatable :: ils(:), neigh(:), irow(:), icol(:)
integer(psb_ipk_), allocatable :: ils(:), neigh(:), irow(:), icol(:),&
& ideg(:), idxs(:)
real(psb_spk_), allocatable :: val(:), diag(:)
integer(psb_ipk_) :: icnt,nlp,k,n,ia,isz,nr, naggr,i,j,m, nz, ilg
integer(psb_ipk_) :: icnt,nlp,k,n,ia,isz,nr, naggr,i,j,m, nz, ilg, ii
type(psb_s_csr_sparse_mat) :: acsr
real(psb_spk_) :: cpling, tcl
logical :: recovery
integer(psb_ipk_) :: debug_level, debug_unit,err_act
@@ -75,7 +78,7 @@ subroutine mld_s_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
ncol = desc_a%get_local_cols()
nr = a%get_nrows()
allocate(ilaggr(nr),neigh(nr),stat=info)
allocate(ilaggr(nr),neigh(nr),ideg(nr),idxs(nr),stat=info)
if(info /= psb_success_) then
info=psb_err_alloc_request_
call psb_errpush(info,name,i_err=(/2*nr,izero,izero,izero,izero/),&
@@ -90,11 +93,20 @@ subroutine mld_s_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
goto 9999
end if
do i=1, nr
ilaggr(i) = -(nr+1)
end do
if (iorder == mld_aggr_ord_nat_) then
do i=1, nr
ilaggr(i) = -(nr+1)
idxs(i) = i
end do
else
call a%cp_to(acsr)
do i=1, nr
ilaggr(i) = -(nr+1)
ideg(i) = acsr%irp(i+1) - acsr%irp(i)
end do
call acsr%free()
call psb_msort(ideg,ix=idxs,dir=psb_sort_down_)
end if
! Note: -(nr+1) Untouched as yet
! -i 1<=i<=nr Adjacent to aggregate i
! i 1<=i<=nr Belonging to aggregate i
@@ -106,7 +118,8 @@ subroutine mld_s_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
nlp = 0
do
icnt = 0
do i=1, nr
do ii=1, nr
i = idxs(ii)
if (ilaggr(i) == -(nr+1)) then
!
! 1. Untouched nodes are marked >0 together
@@ -124,14 +137,16 @@ subroutine mld_s_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
end if
do k=1, nz
j = icol(k)
ilg = ilaggr(j)
if ((ilg<0).and.(1<=j).and.(j<=nr).and.(i /= j)) then
if (abs(val(k)) > theta*sqrt(abs(diag(i)*diag(j)))) then
ilaggr(j) = naggr
else
ilaggr(j) = -naggr
endif
j = icol(k)
if ((1<=j).and.(j<=nr)) then
ilg = ilaggr(j)
if ((ilg<0).and.(i /= j)) then
if (abs(val(k)) > theta*sqrt(abs(diag(i)*diag(j)))) then
ilaggr(j) = naggr
else
ilaggr(j) = -naggr
endif
end if
end if
enddo
mlprec/impl/mld_saggrmap_bld.f90
+5 -4
View File
@@ -79,7 +79,7 @@
! info - integer, output.
! Error code.
!
subroutine mld_saggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
subroutine mld_saggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,info)
use psb_base_mod
use mld_s_inner_mod, mld_protect_name => mld_saggrmap_bld
@@ -87,6 +87,7 @@ subroutine mld_saggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
implicit none
! Arguments
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
@@ -118,7 +119,7 @@ subroutine mld_saggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
select case (aggr_type)
case (mld_dec_aggr_)
call mld_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
call mld_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
case (mld_sym_dec_aggr_)
nr = a%get_nrows()
@@ -126,14 +127,14 @@ subroutine mld_saggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
& rscale=.false.,cscale=.false.)
call atmp%set_nrows(nr)
call atmp%set_ncols(nr)
if (info == psb_success_) call atrans%transp(atmp)
if (info == psb_success_) call atmp%transp(atrans)
if (info == psb_success_) call atrans%cscnv(info,type='COO')
if (info == psb_success_) call psb_rwextd(nr,atmp,info,b=atrans,rowscale=.false.)
call atmp%set_nrows(nr)
call atmp%set_ncols(nr)
if (info == psb_success_) call atrans%free()
if (info == psb_success_) call atmp%cscnv(info,type='CSR')
if (info == psb_success_) call mld_dec_map_bld(theta,atmp,desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call mld_dec_map_bld(iorder,theta,atmp,desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call atmp%free()
case default
mlprec/impl/mld_scoarse_bld.f90
+4 -1
View File
@@ -95,6 +95,8 @@ subroutine mld_scoarse_bld(a,desc_a,p,info)
& mld_mult_ml_,is_legal_ml_type)
call mld_check_def(p%parms%aggr_alg,'Aggregation',&
& mld_dec_aggr_,is_legal_ml_aggr_alg)
call mld_check_def(p%parms%aggr_ord,'Ordering',&
& mld_aggr_ord_nat_,is_legal_ml_aggr_ord)
call mld_check_def(p%parms%aggr_kind,'Smoother',&
& mld_smooth_prol_,is_legal_ml_aggr_kind)
call mld_check_def(p%parms%coarse_mat,'Coarse matrix',&
@@ -117,7 +119,7 @@ subroutine mld_scoarse_bld(a,desc_a,p,info)
! aggregation algorithm. This also defines a tentative prolongator from
! the coarse to the fine level.
!
call mld_aggrmap_bld(p%parms%aggr_alg,p%parms%aggr_thresh,&
call mld_aggrmap_bld(p%parms%aggr_alg,p%parms%aggr_ord,p%parms%aggr_thresh,&
& a,desc_a,ilaggr,nlaggr,info)
if (info /= psb_success_) then
@@ -136,6 +138,7 @@ subroutine mld_scoarse_bld(a,desc_a,p,info)
call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_aggrmat_asb')
goto 9999
end if
!
! Fix the base_a and base_desc pointers for handling of residuals.
! This is correct because this routine is only called at levels >=2.
mlprec/impl/mld_scprecset.F90
+10 -8
View File
@@ -93,6 +93,7 @@ subroutine mld_scprecseti(p,what,val,info,ilev)
use mld_s_mumps_solver
#endif
implicit none
! Arguments
@@ -148,9 +149,10 @@ subroutine mld_scprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case('SUB_SOLVE')
call onelev_set_solver(p%precv(ilev_),val,info)
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_KIND',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG','AGGR_EIG',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_ORD',&
& 'AGGR_KIND','SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','SMOOTHER_SWEEPS_PRE',&
& 'SMOOTHER_SWEEPS_POST',&
& 'SUB_RESTR','SUB_PROL', &
& 'SUB_REN','SUB_OVR','SUB_FILLIN')
call p%precv(ilev_)%set(what,val,info)
@@ -166,9 +168,10 @@ subroutine mld_scprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case('SUB_SOLVE')
call onelev_set_solver(p%precv(ilev_),val,info)
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_KIND',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG','AGGR_EIG',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_ORD',&
& 'AGGR_KIND','SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','SMOOTHER_SWEEPS_PRE',&
& 'SMOOTHER_SWEEPS_POST',&
& 'SUB_RESTR','SUB_PROL', &
& 'SUB_REN','SUB_OVR','SUB_FILLIN',&
& 'COARSE_MAT')
@@ -280,7 +283,7 @@ subroutine mld_scprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
end do
case('ML_TYPE','AGGR_ALG','AGGR_KIND',&
case('ML_TYPE','AGGR_ALG','AGGR_ORD','AGGR_KIND',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','AGGR_FILTER')
@@ -583,7 +586,6 @@ contains
info = -5
end if
#endif
#ifdef HAVE_MUMPS_
case (mld_mumps_)
if (allocated(level%sm%sv)) then
mlprec/impl/mld_sprecinit.F90
+8 -12
View File
@@ -2,9 +2,9 @@
!!$
!!$ MLD2P4 version 2.0
!!$ MultiLevel Domain Decomposition Parallel Preconditioners Package
!!$ based on PSBLAS (Parallel Sparse BLAS version 3.0)
!!$ based on PSBLAS (Parallel Sparse BLAS version 3.3)
!!$
!!$ (C) Copyright 2008,2009,2010,2012,2013
!!$ (C) Copyright 2008, 2010, 2012, 2015
!!$
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari CNRS-IRIT, Toulouse
@@ -61,7 +61,7 @@
!
! 'ML' - Multilevel hybrid preconditioner (additive on the
! same level and multiplicative through the levels),
! with 2 levels and post-smoothing only. RAS with
! with 2 levels, pre and post-smoothing, RAS with
! overlap 1 and ILU(0) on the local blocks is
! applied as post-smoother at each level, but the
! coarsest one; four sweeps of the block-Jacobi solver,
@@ -97,9 +97,6 @@ subroutine mld_sprecinit(p,ptype,info,nlev)
use mld_s_id_solver
use mld_s_diag_solver
use mld_s_ilu_solver
#if defined(HAVE_UMF_) && 0
use mld_s_umf_solver
#endif
#if defined(HAVE_SLU_)
use mld_s_slu_solver
#endif
@@ -115,7 +112,7 @@ subroutine mld_sprecinit(p,ptype,info,nlev)
! Local variables
integer(psb_ipk_) :: nlev_, ilev_
real(psb_spk_) :: thr
real(psb_spk_) :: thr, scale
character(len=*), parameter :: name='mld_precinit'
info = psb_success_
@@ -191,9 +188,7 @@ subroutine mld_sprecinit(p,ptype,info,nlev)
ilev_ = nlev_
allocate(mld_s_jac_smoother_type :: p%precv(ilev_)%sm, stat=info)
if (info /= psb_success_) return
#if defined(HAVE_UMF_) && 0
allocate(mld_s_umf_solver_type :: p%precv(ilev_)%sm%sv, stat=info)
#elif defined(HAVE_SLU_)
#if defined(HAVE_SLU_)
allocate(mld_s_slu_solver_type :: p%precv(ilev_)%sm%sv, stat=info)
#else
allocate(mld_s_ilu_solver_type :: p%precv(ilev_)%sm%sv, stat=info)
@@ -205,10 +200,11 @@ subroutine mld_sprecinit(p,ptype,info,nlev)
call p%precv(ilev_)%set(mld_sub_prol_,psb_none_,info)
call p%precv(ilev_)%set(mld_sub_ovr_,izero,info)
thr = 0.16d0
thr = 0.05
scale = 1.0
do ilev_=1,nlev_
call p%precv(ilev_)%set(mld_aggr_thresh_,thr,info)
thr = thr/2
call p%precv(ilev_)%set(mld_aggr_scale_,scale,info)
end do
case default
mlprec/impl/mld_sprecset.F90
+24 -14
View File
@@ -149,8 +149,8 @@ subroutine mld_sprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case(mld_sub_solve_)
call onelev_set_solver(p%precv(ilev_),val,info)
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,&
& mld_aggr_kind_,mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_sub_restr_,mld_sub_prol_, &
& mld_sub_ren_,mld_sub_ovr_,mld_sub_fillin_)
@@ -167,8 +167,8 @@ subroutine mld_sprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case(mld_sub_solve_)
call onelev_set_solver(p%precv(ilev_),val,info)
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,&
& mld_aggr_kind_,mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_sub_restr_,mld_sub_prol_, &
& mld_sub_ren_,mld_sub_ovr_,mld_sub_fillin_,&
@@ -277,7 +277,7 @@ subroutine mld_sprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
end do
case(mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
case(mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,mld_aggr_kind_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,&
& mld_aggr_eig_,mld_aggr_filter_)
@@ -299,8 +299,8 @@ subroutine mld_sprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
#if defined(HAVE_SLU_)
call onelev_set_solver(p%precv(nlev_),mld_slu_,info)
#elif defined(HAVE_MUMPS_)
call onelev_set_solver(p%precv(nlev_),mld_mumps_,info)
#elif defined(HAVE_SLU_)
call onelev_set_solver(p%precv(nlev_),mld_slu_,info)
#else
call onelev_set_solver(p%precv(nlev_),mld_ilu_n_,info)
#endif
@@ -579,7 +579,6 @@ contains
info = -5
end if
#endif
#ifdef HAVE_MUMPS_
case (mld_mumps_)
if (allocated(level%sm%sv)) then
@@ -601,8 +600,6 @@ contains
end if
end if
#endif
case default
!
! Do nothing and hope for the best :)
@@ -726,13 +723,26 @@ subroutine mld_sprecsetsv(p,val,info,ilev)
do ilev_ = ilmin, ilmax
if (allocated(p%precv(ilev_)%sm)) then
if (allocated(p%precv(ilev_)%sm%sv)) &
& deallocate(p%precv(ilev_)%sm%sv)
if (allocated(p%precv(ilev_)%sm%sv)) then
if (.not.same_type_as(p%precv(ilev_)%sm%sv,val)) then
deallocate(p%precv(ilev_)%sm%sv,stat=info)
if (info /= 0) then
info = 3111
return
end if
end if
if (.not.allocated(p%precv(ilev_)%sm%sv)) then
#ifdef HAVE_MOLD
allocate(p%precv(ilev_)%sm%sv,mold=val)
allocate(p%precv(ilev_)%sm%sv,mold=val,stat=info)
#else
allocate(p%precv(ilev_)%sm%sv,source=val)
allocate(p%precv(ilev_)%sm%sv,source=val,stat=info)
#endif
if (info /= 0) then
info = 3111
return
end if
end if
end if
call p%precv(ilev_)%sm%sv%default()
else
info = 3111
mlprec/impl/mld_z_dec_map_bld.f90
+33 -18
View File
@@ -37,7 +37,7 @@
!!$
!!$
subroutine mld_z_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_z_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use mld_z_inner_mod, mld_protect_name => mld_z_dec_map_bld
@@ -45,6 +45,7 @@ subroutine mld_z_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
implicit none
! Arguments
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
@@ -52,9 +53,11 @@ subroutine mld_z_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
integer(psb_ipk_), intent(out) :: info
! Local variables
integer(psb_ipk_), allocatable :: ils(:), neigh(:), irow(:), icol(:)
integer(psb_ipk_), allocatable :: ils(:), neigh(:), irow(:), icol(:),&
& ideg(:), idxs(:)
complex(psb_dpk_), allocatable :: val(:), diag(:)
integer(psb_ipk_) :: icnt,nlp,k,n,ia,isz,nr, naggr,i,j,m, nz, ilg
integer(psb_ipk_) :: icnt,nlp,k,n,ia,isz,nr, naggr,i,j,m, nz, ilg, ii
type(psb_z_csr_sparse_mat) :: acsr
real(psb_dpk_) :: cpling, tcl
logical :: recovery
integer(psb_ipk_) :: debug_level, debug_unit,err_act
@@ -75,7 +78,7 @@ subroutine mld_z_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
ncol = desc_a%get_local_cols()
nr = a%get_nrows()
allocate(ilaggr(nr),neigh(nr),stat=info)
allocate(ilaggr(nr),neigh(nr),ideg(nr),idxs(nr),stat=info)
if(info /= psb_success_) then
info=psb_err_alloc_request_
call psb_errpush(info,name,i_err=(/2*nr,izero,izero,izero,izero/),&
@@ -90,11 +93,20 @@ subroutine mld_z_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
goto 9999
end if
do i=1, nr
ilaggr(i) = -(nr+1)
end do
if (iorder == mld_aggr_ord_nat_) then
do i=1, nr
ilaggr(i) = -(nr+1)
idxs(i) = i
end do
else
call a%cp_to(acsr)
do i=1, nr
ilaggr(i) = -(nr+1)
ideg(i) = acsr%irp(i+1) - acsr%irp(i)
end do
call acsr%free()
call psb_msort(ideg,ix=idxs,dir=psb_sort_down_)
end if
! Note: -(nr+1) Untouched as yet
! -i 1<=i<=nr Adjacent to aggregate i
! i 1<=i<=nr Belonging to aggregate i
@@ -106,7 +118,8 @@ subroutine mld_z_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
nlp = 0
do
icnt = 0
do i=1, nr
do ii=1, nr
i = idxs(ii)
if (ilaggr(i) == -(nr+1)) then
!
! 1. Untouched nodes are marked >0 together
@@ -124,14 +137,16 @@ subroutine mld_z_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
end if
do k=1, nz
j = icol(k)
ilg = ilaggr(j)
if ((ilg<0).and.(1<=j).and.(j<=nr).and.(i /= j)) then
if (abs(val(k)) > theta*sqrt(abs(diag(i)*diag(j)))) then
ilaggr(j) = naggr
else
ilaggr(j) = -naggr
endif
j = icol(k)
if ((1<=j).and.(j<=nr)) then
ilg = ilaggr(j)
if ((ilg<0).and.(i /= j)) then
if (abs(val(k)) > theta*sqrt(abs(diag(i)*diag(j)))) then
ilaggr(j) = naggr
else
ilaggr(j) = -naggr
endif
end if
end if
enddo
mlprec/impl/mld_zaggrmap_bld.f90
+5 -4
View File
@@ -79,7 +79,7 @@
! info - integer, output.
! Error code.
!
subroutine mld_zaggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
subroutine mld_zaggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,info)
use psb_base_mod
use mld_z_inner_mod, mld_protect_name => mld_zaggrmap_bld
@@ -87,6 +87,7 @@ subroutine mld_zaggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
implicit none
! Arguments
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
@@ -118,7 +119,7 @@ subroutine mld_zaggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
select case (aggr_type)
case (mld_dec_aggr_)
call mld_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
call mld_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
case (mld_sym_dec_aggr_)
nr = a%get_nrows()
@@ -126,14 +127,14 @@ subroutine mld_zaggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
& rscale=.false.,cscale=.false.)
call atmp%set_nrows(nr)
call atmp%set_ncols(nr)
if (info == psb_success_) call atrans%transp(atmp)
if (info == psb_success_) call atmp%transp(atrans)
if (info == psb_success_) call atrans%cscnv(info,type='COO')
if (info == psb_success_) call psb_rwextd(nr,atmp,info,b=atrans,rowscale=.false.)
call atmp%set_nrows(nr)
call atmp%set_ncols(nr)
if (info == psb_success_) call atrans%free()
if (info == psb_success_) call atmp%cscnv(info,type='CSR')
if (info == psb_success_) call mld_dec_map_bld(theta,atmp,desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call mld_dec_map_bld(iorder,theta,atmp,desc_a,nlaggr,ilaggr,info)
if (info == psb_success_) call atmp%free()
case default
mlprec/impl/mld_zcoarse_bld.f90
+5 -3
View File
@@ -39,7 +39,7 @@
! File: mld_zcoarse_bld.f90
!
! Subroutine: mld_zcoarse_bld
! Version: complex
! Version: real
!
! This routine builds the matrix associated to the current level of the
! multilevel preconditioner from the matrix associated to the previous level,
@@ -95,6 +95,8 @@ subroutine mld_zcoarse_bld(a,desc_a,p,info)
& mld_mult_ml_,is_legal_ml_type)
call mld_check_def(p%parms%aggr_alg,'Aggregation',&
& mld_dec_aggr_,is_legal_ml_aggr_alg)
call mld_check_def(p%parms%aggr_ord,'Ordering',&
& mld_aggr_ord_nat_,is_legal_ml_aggr_ord)
call mld_check_def(p%parms%aggr_kind,'Smoother',&
& mld_smooth_prol_,is_legal_ml_aggr_kind)
call mld_check_def(p%parms%coarse_mat,'Coarse matrix',&
@@ -117,7 +119,7 @@ subroutine mld_zcoarse_bld(a,desc_a,p,info)
! aggregation algorithm. This also defines a tentative prolongator from
! the coarse to the fine level.
!
call mld_aggrmap_bld(p%parms%aggr_alg,p%parms%aggr_thresh,&
call mld_aggrmap_bld(p%parms%aggr_alg,p%parms%aggr_ord,p%parms%aggr_thresh,&
& a,desc_a,ilaggr,nlaggr,info)
if (info /= psb_success_) then
@@ -136,6 +138,7 @@ subroutine mld_zcoarse_bld(a,desc_a,p,info)
call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_aggrmat_asb')
goto 9999
end if
!
! Fix the base_a and base_desc pointers for handling of residuals.
! This is correct because this routine is only called at levels >=2.
@@ -147,7 +150,6 @@ subroutine mld_zcoarse_bld(a,desc_a,p,info)
return
9999 call psb_error_handler(err_act)
return
end subroutine mld_zcoarse_bld
mlprec/impl/mld_zcprecset.F90
+31 -30
View File
@@ -99,6 +99,7 @@ subroutine mld_zcprecseti(p,what,val,info,ilev)
use mld_z_mumps_solver
#endif
implicit none
! Arguments
@@ -154,9 +155,10 @@ subroutine mld_zcprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case('SUB_SOLVE')
call onelev_set_solver(p%precv(ilev_),val,info)
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_KIND',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG','AGGR_EIG',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_ORD',&
& 'AGGR_KIND','SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','SMOOTHER_SWEEPS_PRE',&
& 'SMOOTHER_SWEEPS_POST',&
& 'SUB_RESTR','SUB_PROL', &
& 'SUB_REN','SUB_OVR','SUB_FILLIN')
call p%precv(ilev_)%set(what,val,info)
@@ -172,9 +174,10 @@ subroutine mld_zcprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case('SUB_SOLVE')
call onelev_set_solver(p%precv(ilev_),val,info)
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_KIND',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG','AGGR_EIG',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
case('SMOOTHER_SWEEPS','ML_TYPE','AGGR_ALG','AGGR_ORD',&
& 'AGGR_KIND','SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','SMOOTHER_SWEEPS_PRE',&
& 'SMOOTHER_SWEEPS_POST',&
& 'SUB_RESTR','SUB_PROL', &
& 'SUB_REN','SUB_OVR','SUB_FILLIN',&
& 'COARSE_MAT')
@@ -288,7 +291,7 @@ subroutine mld_zcprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
end do
case('ML_TYPE','AGGR_ALG','AGGR_KIND',&
case('ML_TYPE','AGGR_ALG','AGGR_ORD','AGGR_KIND',&
& 'SMOOTHER_SWEEPS_PRE','SMOOTHER_SWEEPS_POST',&
& 'SMOOTHER_POS','AGGR_OMEGA_ALG',&
& 'AGGR_EIG','AGGR_FILTER')
@@ -330,7 +333,6 @@ subroutine mld_zcprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
call onelev_set_solver(p%precv(nlev_),val,info)
call p%precv(nlev_)%set('COARSE_MAT',mld_distr_mat_,info)
case(mld_jac_)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
call onelev_set_solver(p%precv(nlev_),mld_diag_scale_,info)
@@ -594,6 +596,27 @@ contains
info = -5
end if
#endif
#ifdef HAVE_MUMPS_
case (mld_mumps_)
if (allocated(level%sm%sv)) then
select type (sv => level%sm%sv)
class is (mld_z_mumps_solver_type)
! do nothing
class default
call level%sm%sv%free(info)
if (info == 0) deallocate(level%sm%sv)
if (info == 0) allocate(mld_z_mumps_solver_type ::&
& level%sm%sv, stat=info)
end select
else
allocate(mld_z_mumps_solver_type :: level%sm%sv, stat=info)
endif
if (allocated(level%sm)) then
if (allocated(level%sm%sv)) &
& call level%sm%sv%default()
end if
#endif
#ifdef HAVE_UMF_
case (mld_umf_)
if (allocated(level%sm)) then
@@ -644,28 +667,6 @@ contains
info = -5
end if
#endif
#ifdef HAVE_MUMPS_
case (mld_mumps_)
if (allocated(level%sm%sv)) then
select type (sv => level%sm%sv)
class is (mld_z_mumps_solver_type)
! do nothing
class default
call level%sm%sv%free(info)
if (info == 0) deallocate(level%sm%sv)
if (info == 0) allocate(mld_z_mumps_solver_type ::&
& level%sm%sv, stat=info)
end select
else
allocate(mld_z_mumps_solver_type :: level%sm%sv, stat=info)
endif
if (allocated(level%sm)) then
if (allocated(level%sm%sv)) &
& call level%sm%sv%default()
end if
#endif
case default
!
! Do nothing and hope for the best :)
mlprec/impl/mld_zprecinit.F90
+7 -6
View File
@@ -2,9 +2,9 @@
!!$
!!$ MLD2P4 version 2.0
!!$ MultiLevel Domain Decomposition Parallel Preconditioners Package
!!$ based on PSBLAS (Parallel Sparse BLAS version 3.0)
!!$ based on PSBLAS (Parallel Sparse BLAS version 3.3)
!!$
!!$ (C) Copyright 2008,2009,2010,2012,2013
!!$ (C) Copyright 2008, 2010, 2012, 2015
!!$
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari CNRS-IRIT, Toulouse
@@ -61,7 +61,7 @@
!
! 'ML' - Multilevel hybrid preconditioner (additive on the
! same level and multiplicative through the levels),
! with 2 levels and post-smoothing only. RAS with
! with 2 levels, pre and post-smoothing, RAS with
! overlap 1 and ILU(0) on the local blocks is
! applied as post-smoother at each level, but the
! coarsest one; four sweeps of the block-Jacobi solver,
@@ -115,7 +115,7 @@ subroutine mld_zprecinit(p,ptype,info,nlev)
! Local variables
integer(psb_ipk_) :: nlev_, ilev_
real(psb_dpk_) :: thr
real(psb_dpk_) :: thr, scale
character(len=*), parameter :: name='mld_precinit'
info = psb_success_
@@ -205,10 +205,11 @@ subroutine mld_zprecinit(p,ptype,info,nlev)
call p%precv(ilev_)%set(mld_sub_prol_,psb_none_,info)
call p%precv(ilev_)%set(mld_sub_ovr_,izero,info)
thr = 0.16d0
thr = 0.05
scale = 1.0
do ilev_=1,nlev_
call p%precv(ilev_)%set(mld_aggr_thresh_,thr,info)
thr = thr/2
call p%precv(ilev_)%set(mld_aggr_scale_,scale,info)
end do
case default
mlprec/impl/mld_zprecset.F90
+24 -12
View File
@@ -98,6 +98,7 @@ subroutine mld_zprecseti(p,what,val,info,ilev)
#if defined(HAVE_MUMPS_)
use mld_z_mumps_solver
#endif
implicit none
! Arguments
@@ -154,8 +155,8 @@ subroutine mld_zprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case(mld_sub_solve_)
call onelev_set_solver(p%precv(ilev_),val,info)
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,&
& mld_aggr_kind_,mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_sub_restr_,mld_sub_prol_, &
& mld_sub_ren_,mld_sub_ovr_,mld_sub_fillin_)
@@ -172,8 +173,8 @@ subroutine mld_zprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
case(mld_sub_solve_)
call onelev_set_solver(p%precv(ilev_),val,info)
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
case(mld_smoother_sweeps_,mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,&
& mld_aggr_kind_,mld_smoother_pos_,mld_aggr_omega_alg_,mld_aggr_eig_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_sub_restr_,mld_sub_prol_, &
& mld_sub_ren_,mld_sub_ovr_,mld_sub_fillin_,&
@@ -215,7 +216,7 @@ subroutine mld_zprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
call onelev_set_solver(p%precv(nlev_),val,info)
call p%precv(nlev_)%set(mld_coarse_mat_,mld_repl_mat_,info)
case(mld_sludist_, mld_mumps_)
case(mld_sludist_,mld_mumps_)
call onelev_set_smoother(p%precv(nlev_),mld_bjac_,info)
call onelev_set_solver(p%precv(nlev_),val,info)
call p%precv(nlev_)%set(mld_coarse_mat_,mld_distr_mat_,info)
@@ -284,7 +285,7 @@ subroutine mld_zprecseti(p,what,val,info,ilev)
call onelev_set_smoother(p%precv(ilev_),val,info)
end do
case(mld_ml_type_,mld_aggr_alg_,mld_aggr_kind_,&
case(mld_ml_type_,mld_aggr_alg_,mld_aggr_ord_,mld_aggr_kind_,&
& mld_smoother_sweeps_pre_,mld_smoother_sweeps_post_,&
& mld_smoother_pos_,mld_aggr_omega_alg_,&
& mld_aggr_eig_,mld_aggr_filter_)
@@ -308,8 +309,6 @@ subroutine mld_zprecseti(p,what,val,info,ilev)
call onelev_set_solver(p%precv(nlev_),mld_umf_,info)
#elif defined(HAVE_SLU_)
call onelev_set_solver(p%precv(nlev_),mld_slu_,info)
#elif defined(HAVE_MUMPS_)
call onelev_set_solver(p%precv(nlev_),mld_mumps_,info)
#else
call onelev_set_solver(p%precv(nlev_),mld_ilu_n_,info)
#endif
@@ -782,13 +781,26 @@ subroutine mld_zprecsetsv(p,val,info,ilev)
do ilev_ = ilmin, ilmax
if (allocated(p%precv(ilev_)%sm)) then
if (allocated(p%precv(ilev_)%sm%sv)) &
& deallocate(p%precv(ilev_)%sm%sv)
if (allocated(p%precv(ilev_)%sm%sv)) then
if (.not.same_type_as(p%precv(ilev_)%sm%sv,val)) then
deallocate(p%precv(ilev_)%sm%sv,stat=info)
if (info /= 0) then
info = 3111
return
end if
end if
if (.not.allocated(p%precv(ilev_)%sm%sv)) then
#ifdef HAVE_MOLD
allocate(p%precv(ilev_)%sm%sv,mold=val)
allocate(p%precv(ilev_)%sm%sv,mold=val,stat=info)
#else
allocate(p%precv(ilev_)%sm%sv,source=val)
allocate(p%precv(ilev_)%sm%sv,source=val,stat=info)
#endif
if (info /= 0) then
info = 3111
return
end if
end if
end if
call p%precv(ilev_)%sm%sv%default()
else
info = 3111
mlprec/mld_base_prec_type.F90
+38 -16
View File
@@ -98,7 +98,7 @@ module mld_base_prec_type
type mld_ml_parms
integer(psb_ipk_) :: sweeps, sweeps_pre, sweeps_post
integer(psb_ipk_) :: ml_type, smoother_pos
integer(psb_ipk_) :: aggr_alg, aggr_kind
integer(psb_ipk_) :: aggr_alg, aggr_ord, aggr_kind
integer(psb_ipk_) :: aggr_omega_alg, aggr_eig, aggr_filter
integer(psb_ipk_) :: coarse_mat, coarse_solve
logical :: clean_zeros=.true.
@@ -151,18 +151,19 @@ module mld_base_prec_type
integer(psb_ipk_), parameter :: mld_smoother_pos_ = 23
integer(psb_ipk_), parameter :: mld_aggr_kind_ = 24
integer(psb_ipk_), parameter :: mld_aggr_alg_ = 25
integer(psb_ipk_), parameter :: mld_aggr_omega_alg_ = 26
integer(psb_ipk_), parameter :: mld_aggr_eig_ = 27
integer(psb_ipk_), parameter :: mld_aggr_filter_ = 28
integer(psb_ipk_), parameter :: mld_coarse_mat_ = 29
integer(psb_ipk_), parameter :: mld_coarse_solve_ = 30
integer(psb_ipk_), parameter :: mld_coarse_sweeps_ = 31
integer(psb_ipk_), parameter :: mld_coarse_fillin_ = 32
integer(psb_ipk_), parameter :: mld_coarse_subsolve_ = 33
integer(psb_ipk_), parameter :: mld_smoother_sweeps_ = 34
integer(psb_ipk_), parameter :: mld_coarse_aggr_size_ = 35
integer(psb_ipk_), parameter :: mld_solver_sweeps_ = 36
integer(psb_ipk_), parameter :: mld_ifpsz_ = 37
integer(psb_ipk_), parameter :: mld_aggr_ord_ = 26
integer(psb_ipk_), parameter :: mld_aggr_omega_alg_ = 27
integer(psb_ipk_), parameter :: mld_aggr_eig_ = 28
integer(psb_ipk_), parameter :: mld_aggr_filter_ = 29
integer(psb_ipk_), parameter :: mld_coarse_mat_ = 30
integer(psb_ipk_), parameter :: mld_coarse_solve_ = 31
integer(psb_ipk_), parameter :: mld_coarse_sweeps_ = 32
integer(psb_ipk_), parameter :: mld_coarse_fillin_ = 33
integer(psb_ipk_), parameter :: mld_coarse_subsolve_ = 34
integer(psb_ipk_), parameter :: mld_smoother_sweeps_ = 35
integer(psb_ipk_), parameter :: mld_coarse_aggr_size_ = 36
integer(psb_ipk_), parameter :: mld_solver_sweeps_ = 37
integer(psb_ipk_), parameter :: mld_ifpsz_ = 38
!
! Legal values for entry: mld_smoother_type_
@@ -249,8 +250,13 @@ module mld_base_prec_type
integer(psb_ipk_), parameter :: mld_glb_aggr_=2
integer(psb_ipk_), parameter :: mld_new_dec_aggr_=3
integer(psb_ipk_), parameter :: mld_new_glb_aggr_=4
integer(psb_ipk_), parameter :: mld_max_aggr_alg_=mld_dec_aggr_
integer(psb_ipk_), parameter :: mld_max_aggr_alg_=mld_sym_dec_aggr_
!
! Legal values for entry: mld_aggr_ord_
!
integer(psb_ipk_), parameter :: mld_aggr_ord_nat_ = 0
integer(psb_ipk_), parameter :: mld_aggr_ord_desc_deg_ = 1
integer(psb_ipk_), parameter :: mld_max_aggr_ord_ = mld_aggr_ord_desc_deg_
!
! Legal values for entry: mld_aggr_omega_alg_
!
@@ -320,6 +326,8 @@ module mld_base_prec_type
character(len=18), parameter, private :: &
& aggr_names(0:4)=(/'local aggregation ','sym. local aggr. ',&
& 'global aggregation', 'new local aggr. ','new global aggr. '/)
character(len=18), parameter, private :: &
& ord_names(0:1)=(/'Natural ordering ','Desc. degree ord. '/)
character(len=6), parameter, private :: &
& restrict_names(0:4)=(/'none ','halo ',' ',' ',' '/)
character(len=12), parameter, private :: &
@@ -411,6 +419,10 @@ contains
val = mld_dec_aggr_
case('SYMDEC')
val = mld_sym_dec_aggr_
case('NAT','NATURAL')
val = mld_aggr_ord_nat_
case('DESC','RDEGREE','DEGREE')
val = mld_aggr_ord_desc_deg_
case('GLB')
val = mld_glb_aggr_
case('REPL')
@@ -467,7 +479,7 @@ contains
write(iout,*) 'Sweeps: ',pm%sweeps,pm%sweeps_pre,pm%sweeps_post
write(iout,*) 'ML : ',pm%ml_type,pm%smoother_pos
write(iout,*) 'AGGR : ',pm%aggr_alg,pm%aggr_kind
write(iout,*) 'AGGR : ',pm%aggr_alg,pm%aggr_kind, pm%aggr_ord
write(iout,*) ' : ',pm%aggr_omega_alg,pm%aggr_eig,pm%aggr_filter
write(iout,*) 'COARSE: ',pm%coarse_mat,pm%coarse_solve
end subroutine ml_parms_printout
@@ -533,6 +545,8 @@ contains
end if
write(iout,*) ' Aggregation: ', &
& aggr_names(pm%aggr_alg)
write(iout,*) ' with initial ordering: ',&
& ord_names(pm%aggr_ord)
write(iout,*) ' Aggregation type: ', &
& aggr_kinds(pm%aggr_kind)
if (pm%aggr_kind /= mld_no_smooth_) then
@@ -720,6 +734,14 @@ contains
is_legal_ml_aggr_alg = ((ip>=mld_dec_aggr_).and.(ip<=mld_max_aggr_alg_))
return
end function is_legal_ml_aggr_alg
function is_legal_ml_aggr_ord(ip)
implicit none
integer(psb_ipk_), intent(in) :: ip
logical :: is_legal_ml_aggr_ord
is_legal_ml_aggr_ord = ((mld_aggr_ord_nat_<=ip).and.(ip<=mld_max_aggr_ord_))
return
end function is_legal_ml_aggr_ord
function is_legal_ml_aggr_omega_alg(ip)
implicit none
integer(psb_ipk_), intent(in) :: ip
mlprec/mld_c_inner_mod.f90
+4 -2
View File
@@ -109,9 +109,10 @@ module mld_c_inner_mod
end interface mld_coarse_bld
interface mld_aggrmap_bld
subroutine mld_caggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
subroutine mld_caggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,info)
use psb_base_mod, only : psb_cspmat_type, psb_desc_type, psb_spk_, psb_ipk_
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
@@ -123,9 +124,10 @@ module mld_c_inner_mod
interface mld_dec_map_bld
subroutine mld_c_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_c_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod, only : psb_cspmat_type, psb_desc_type, psb_spk_, psb_ipk_
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
mlprec/mld_c_onelev_mod.f90
+1
View File
@@ -379,6 +379,7 @@ contains
lv%parms%sweeps_post = 1
lv%parms%ml_type = mld_mult_ml_
lv%parms%aggr_alg = mld_dec_aggr_
lv%parms%aggr_ord = mld_aggr_ord_nat_
lv%parms%aggr_kind = mld_smooth_prol_
lv%parms%coarse_mat = mld_distr_mat_
lv%parms%smoother_pos = mld_twoside_smooth_
mlprec/mld_d_base_solver_mod.f90
+1
View File
@@ -360,6 +360,7 @@ contains
integer(psb_long_int_k_) :: val
integer(psb_ipk_) :: i
val = 0
return
end function d_base_solver_sizeof
mlprec/mld_d_inner_mod.f90
+4 -2
View File
@@ -109,9 +109,10 @@ module mld_d_inner_mod
end interface mld_coarse_bld
interface mld_aggrmap_bld
subroutine mld_daggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
subroutine mld_daggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,info)
use psb_base_mod, only : psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_
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
@@ -123,9 +124,10 @@ module mld_d_inner_mod
interface mld_dec_map_bld
subroutine mld_d_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_d_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod, only : psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_
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
mlprec/mld_d_onelev_mod.f90
+2
View File
@@ -338,6 +338,7 @@ contains
class(mld_d_onelev_type), intent(in) :: lv
integer(psb_long_int_k_) :: val
integer(psb_ipk_) :: i
val = 0
val = val + lv%desc_ac%sizeof()
val = val + lv%ac%sizeof()
@@ -378,6 +379,7 @@ contains
lv%parms%sweeps_post = 1
lv%parms%ml_type = mld_mult_ml_
lv%parms%aggr_alg = mld_dec_aggr_
lv%parms%aggr_ord = mld_aggr_ord_nat_
lv%parms%aggr_kind = mld_smooth_prol_
lv%parms%coarse_mat = mld_distr_mat_
lv%parms%smoother_pos = mld_twoside_smooth_
mlprec/mld_d_prec_type.f90
+1 -1
View File
@@ -327,7 +327,6 @@ contains
class(mld_dprec_type), intent(in) :: prec
integer(psb_long_int_k_) :: val
integer(psb_ipk_) :: i
val = 0
val = val + psb_sizeof_int
if (allocated(prec%precv)) then
@@ -530,6 +529,7 @@ contains
! Local variables
integer(psb_ipk_) :: me,err_act,i
character(len=20) :: name
if(psb_get_errstatus().ne.0) return
info=psb_success_
name = 'mld_dprecfree'
mlprec/mld_d_umf_solver.F90
+1 -3
View File
@@ -133,6 +133,7 @@ contains
integer :: ictxt,np,me,i, err_act
character :: trans_
character(len=20) :: name='d_umf_solver_apply'
call psb_erractionsave(err_act)
info = psb_success_
@@ -216,7 +217,6 @@ contains
integer :: err_act
character(len=20) :: name='d_umf_solver_apply_vect'
call psb_erractionsave(err_act)
call psb_erractionsave(err_act)
info = psb_success_
@@ -258,8 +258,6 @@ contains
integer :: ictxt,np,me,i, err_act, debug_unit, debug_level
character(len=20) :: name='d_umf_solver_bld', ch_err
call psb_erractionsave(err_act)
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
mlprec/mld_s_inner_mod.f90
+4 -2
View File
@@ -109,9 +109,10 @@ module mld_s_inner_mod
end interface mld_coarse_bld
interface mld_aggrmap_bld
subroutine mld_saggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
subroutine mld_saggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,info)
use psb_base_mod, only : psb_sspmat_type, psb_desc_type, psb_spk_, psb_ipk_
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
@@ -123,9 +124,10 @@ module mld_s_inner_mod
interface mld_dec_map_bld
subroutine mld_s_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_s_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod, only : psb_sspmat_type, psb_desc_type, psb_spk_, psb_ipk_
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
mlprec/mld_s_onelev_mod.f90
+1
View File
@@ -379,6 +379,7 @@ contains
lv%parms%sweeps_post = 1
lv%parms%ml_type = mld_mult_ml_
lv%parms%aggr_alg = mld_dec_aggr_
lv%parms%aggr_ord = mld_aggr_ord_nat_
lv%parms%aggr_kind = mld_smooth_prol_
lv%parms%coarse_mat = mld_distr_mat_
lv%parms%smoother_pos = mld_twoside_smooth_
mlprec/mld_z_inner_mod.f90
+4 -2
View File
@@ -109,9 +109,10 @@ module mld_z_inner_mod
end interface mld_coarse_bld
interface mld_aggrmap_bld
subroutine mld_zaggrmap_bld(aggr_type,theta,a,desc_a,ilaggr,nlaggr,info)
subroutine mld_zaggrmap_bld(aggr_type,iorder,theta,a,desc_a,ilaggr,nlaggr,info)
use psb_base_mod, only : psb_zspmat_type, psb_desc_type, psb_dpk_, psb_ipk_
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
@@ -123,9 +124,10 @@ module mld_z_inner_mod
interface mld_dec_map_bld
subroutine mld_z_dec_map_bld(theta,a,desc_a,nlaggr,ilaggr,info)
subroutine mld_z_dec_map_bld(iorder,theta,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod, only : psb_zspmat_type, psb_desc_type, psb_dpk_, psb_ipk_
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
mlprec/mld_z_onelev_mod.f90
+1
View File
@@ -379,6 +379,7 @@ contains
lv%parms%sweeps_post = 1
lv%parms%ml_type = mld_mult_ml_
lv%parms%aggr_alg = mld_dec_aggr_
lv%parms%aggr_ord = mld_aggr_ord_nat_
lv%parms%aggr_kind = mld_smooth_prol_
lv%parms%coarse_mat = mld_distr_mat_
lv%parms%smoother_pos = mld_twoside_smooth_
tests/pdegen/ppde2d.f90
+4
View File
@@ -159,6 +159,7 @@ program ppde2d
integer(psb_ipk_) :: nlev ! Number of levels in multilevel prec.
character(len=16) :: aggrkind ! smoothed/raw aggregatin
character(len=16) :: aggr_alg ! local or global aggregation
character(len=16) :: aggr_ord ! Ordering for aggregation
character(len=16) :: mltype ! additive or multiplicative 2nd level prec
character(len=16) :: smthpos ! side: pre, post, both smoothing
integer(psb_ipk_) :: csize ! aggregation size at which to stop.
@@ -241,6 +242,7 @@ program ppde2d
call mld_precset(prec,'sub_iluthrs', prectype%thr1, info)
call mld_precset(prec,'aggr_kind', prectype%aggrkind,info)
call mld_precset(prec,'aggr_alg', prectype%aggr_alg,info)
call mld_precset(prec,'aggr_ord', prectype%aggr_ord,info)
call mld_precset(prec,'ml_type', prectype%mltype, info)
call mld_precset(prec,'smoother_pos', prectype%smthpos, info)
if (prectype%athres >= dzero) &
@@ -383,6 +385,7 @@ contains
call read_data(prectype%nlev,psb_inp_unit) ! Number of levels in multilevel prec.
call read_data(prectype%aggrkind,psb_inp_unit) ! smoothed/raw aggregatin
call read_data(prectype%aggr_alg,psb_inp_unit) ! local or global aggregation
call read_data(prectype%aggr_ord,psb_inp_unit) ! aggregation ordering
call read_data(prectype%mltype,psb_inp_unit) ! additive or multiplicative 2nd level prec
call read_data(prectype%smthpos,psb_inp_unit) ! side: pre, post, both smoothing
call read_data(prectype%cmat,psb_inp_unit) ! coarse mat
@@ -422,6 +425,7 @@ contains
call psb_bcast(ictxt,prectype%nlev) ! Number of levels in multilevel prec.
call psb_bcast(ictxt,prectype%aggrkind) ! smoothed/raw aggregatin
call psb_bcast(ictxt,prectype%aggr_alg) ! local or global aggregation
call psb_bcast(ictxt,prectype%aggr_ord) ! aggregation ordering
call psb_bcast(ictxt,prectype%mltype) ! additive or multiplicative 2nd level prec
call psb_bcast(ictxt,prectype%smthpos) ! side: pre, post, both smoothing
call psb_bcast(ictxt,prectype%cmat) ! coarse mat
tests/pdegen/ppde3d.f90
+4 -1
View File
@@ -171,6 +171,7 @@ program ppde3d
integer(psb_ipk_) :: nlev ! Number of levels in multilevel prec.
character(len=16) :: aggrkind ! smoothed/raw aggregatin
character(len=16) :: aggr_alg ! local or global aggregation
character(len=16) :: aggr_ord ! Ordering for aggregation
character(len=16) :: mltype ! additive or multiplicative 2nd level prec
character(len=16) :: smthpos ! side: pre, post, both smoothing
integer(psb_ipk_) :: csize ! aggregation size at which to stop.
@@ -253,6 +254,7 @@ program ppde3d
call mld_precset(prec,'sub_iluthrs', prectype%thr1, info)
call mld_precset(prec,'aggr_kind', prectype%aggrkind,info)
call mld_precset(prec,'aggr_alg', prectype%aggr_alg,info)
call mld_precset(prec,'aggr_ord', prectype%aggr_ord,info)
call mld_precset(prec,'ml_type', prectype%mltype, info)
call mld_precset(prec,'smoother_pos', prectype%smthpos, info)
if (prectype%athres >= dzero) &
@@ -396,6 +398,7 @@ contains
call read_data(prectype%nlev,psb_inp_unit) ! Number of levels in multilevel prec.
call read_data(prectype%aggrkind,psb_inp_unit) ! smoothed/raw aggregatin
call read_data(prectype%aggr_alg,psb_inp_unit) ! local or global aggregation
call read_data(prectype%aggr_ord,psb_inp_unit) ! aggregation ordering
call read_data(prectype%mltype,psb_inp_unit) ! additive or multiplicative 2nd level prec
call read_data(prectype%smthpos,psb_inp_unit) ! side: pre, post, both smoothing
call read_data(prectype%cmat,psb_inp_unit) ! coarse mat
@@ -435,6 +438,7 @@ contains
call psb_bcast(ictxt,prectype%nlev) ! Number of levels in multilevel prec.
call psb_bcast(ictxt,prectype%aggrkind) ! smoothed/raw aggregatin
call psb_bcast(ictxt,prectype%aggr_alg) ! local or global aggregation
call psb_bcast(ictxt,prectype%aggr_ord) ! aggregation ordering
call psb_bcast(ictxt,prectype%mltype) ! additive or multiplicative 2nd level prec
call psb_bcast(ictxt,prectype%smthpos) ! side: pre, post, both smoothing
call psb_bcast(ictxt,prectype%cmat) ! coarse mat
@@ -482,4 +486,3 @@ contains
end subroutine pr_usage
end program ppde3d
tests/pdegen/runs/ppde.inp
+5 -4
View File
@@ -3,15 +3,15 @@ CSR ! Storage format CSR COO JAD
0100 ! IDIM; domain size is idim**3
2 ! ISTOPC
0100 ! ITMAX
1 ! ITRACE
10 ! ITRACE
30 ! IRST (restart for RGMRES and BiCGSTABL)
1.d-6 ! EPS
3L-MUL-RAS-BJAC4-ILU ! Descriptive name for preconditioner (up to 40 chars)
ML ! Preconditioner NONE JACOBI BJAC AS ML
ML ! Preconditioner NONE JACOBI BJAC AS ML
1 ! Number of overlap layers for AS preconditioner at finest level
HALO ! Restriction operator NONE HALO
NONE ! Prolongation operator NONE SUM AVG
GS ! Subdomain solver DSCALE ILU MILU ILUT UMF SLU
GS ! Subdomain solver DSCALE ILU MILU ILUT UMF SLU
1 ! sweeps for GS
0 ! Level-set N for ILU(N), and P for ILUT
1.d-4 ! Threshold T for ILU(T,P)
@@ -20,11 +20,12 @@ BJAC ! Smoother type JACOBI BJAC AS; ignored for non-ML
2 ! Number of levels in a multilevel preconditioner
SMOOTHED ! Kind of aggregation: SMOOTHED, NONSMOOTHED, MINENERGY
DEC ! Type of aggregation DEC SYMDEC GLB
DEGREE ! Ordering of aggregation NATURAL DEGREE
MULT ! Type of multilevel correction: ADD MULT
TWOSIDE ! Side of correction PRE POST TWOSIDE (ignored for ADD)
REPL ! Coarse level: matrix distribution DIST REPL
BJAC ! Coarse level: solver JACOBI BJAC UMF SLU SLUDIST MUMPS
MUMPS ! Coarse level: subsolver DSCALE ILU UMF SLU SLUDIST MUMPS
MUMPS ! Coarse level: subsolver DSCALE ILU UMF SLU SLUDIST MUMPS
1 ! Coarse level: Level-set N for ILU(N)
1.d-4 ! Coarse level: Threshold T for ILU(T,P)
4 ! Coarse level: Number of Jacobi sweeps
tests/pdegen/spde2d.f90
+4
View File
@@ -159,6 +159,7 @@ program spde2d
integer(psb_ipk_) :: nlev ! Number of levels in multilevel prec.
character(len=16) :: aggrkind ! smoothed/raw aggregatin
character(len=16) :: aggr_alg ! local or global aggregation
character(len=16) :: aggr_ord ! Ordering for aggregation
character(len=16) :: mltype ! additive or multiplicative 2nd level prec
character(len=16) :: smthpos ! side: pre, post, both smoothing
integer(psb_ipk_) :: csize ! aggregation size at which to stop.
@@ -241,6 +242,7 @@ program spde2d
call mld_precset(prec,'sub_iluthrs', prectype%thr1, info)
call mld_precset(prec,'aggr_kind', prectype%aggrkind,info)
call mld_precset(prec,'aggr_alg', prectype%aggr_alg,info)
call mld_precset(prec,'aggr_ord', prectype%aggr_ord,info)
call mld_precset(prec,'ml_type', prectype%mltype, info)
call mld_precset(prec,'smoother_pos', prectype%smthpos, info)
if (prectype%athres >= dzero) &
@@ -384,6 +386,7 @@ contains
call read_data(prectype%nlev,psb_inp_unit) ! Number of levels in multilevel prec.
call read_data(prectype%aggrkind,psb_inp_unit) ! smoothed/raw aggregatin
call read_data(prectype%aggr_alg,psb_inp_unit) ! local or global aggregation
call read_data(prectype%aggr_ord,psb_inp_unit) ! aggregation ordering
call read_data(prectype%mltype,psb_inp_unit) ! additive or multiplicative 2nd level prec
call read_data(prectype%smthpos,psb_inp_unit) ! side: pre, post, both smoothing
call read_data(prectype%cmat,psb_inp_unit) ! coarse mat
@@ -423,6 +426,7 @@ contains
call psb_bcast(ictxt,prectype%nlev) ! Number of levels in multilevel prec.
call psb_bcast(ictxt,prectype%aggrkind) ! smoothed/raw aggregatin
call psb_bcast(ictxt,prectype%aggr_alg) ! local or global aggregation
call psb_bcast(ictxt,prectype%aggr_ord) ! aggregation ordering
call psb_bcast(ictxt,prectype%mltype) ! additive or multiplicative 2nd level prec
call psb_bcast(ictxt,prectype%smthpos) ! side: pre, post, both smoothing
call psb_bcast(ictxt,prectype%cmat) ! coarse mat
tests/pdegen/spde3d.f90
+4
View File
@@ -171,6 +171,7 @@ program spde3d
integer(psb_ipk_) :: nlev ! Number of levels in multilevel prec.
character(len=16) :: aggrkind ! smoothed/raw aggregatin
character(len=16) :: aggr_alg ! local or global aggregation
character(len=16) :: aggr_ord ! Ordering for aggregation
character(len=16) :: mltype ! additive or multiplicative 2nd level prec
character(len=16) :: smthpos ! side: pre, post, both smoothing
integer(psb_ipk_) :: csize ! aggregation size at which to stop.
@@ -254,6 +255,7 @@ program spde3d
call mld_precset(prec,'sub_iluthrs', prectype%thr1, info)
call mld_precset(prec,'aggr_kind', prectype%aggrkind,info)
call mld_precset(prec,'aggr_alg', prectype%aggr_alg,info)
call mld_precset(prec,'aggr_ord', prectype%aggr_ord,info)
call mld_precset(prec,'ml_type', prectype%mltype, info)
call mld_precset(prec,'smoother_pos', prectype%smthpos, info)
if (prectype%athres >= dzero) &
@@ -397,6 +399,7 @@ contains
call read_data(prectype%nlev,psb_inp_unit) ! Number of levels in multilevel prec.
call read_data(prectype%aggrkind,psb_inp_unit) ! smoothed/raw aggregatin
call read_data(prectype%aggr_alg,psb_inp_unit) ! local or global aggregation
call read_data(prectype%aggr_ord,psb_inp_unit) ! aggregation ordering
call read_data(prectype%mltype,psb_inp_unit) ! additive or multiplicative 2nd level prec
call read_data(prectype%smthpos,psb_inp_unit) ! side: pre, post, both smoothing
call read_data(prectype%cmat,psb_inp_unit) ! coarse mat
@@ -436,6 +439,7 @@ contains
call psb_bcast(ictxt,prectype%nlev) ! Number of levels in multilevel prec.
call psb_bcast(ictxt,prectype%aggrkind) ! smoothed/raw aggregatin
call psb_bcast(ictxt,prectype%aggr_alg) ! local or global aggregation
call psb_bcast(ictxt,prectype%aggr_ord) ! aggregation ordering
call psb_bcast(ictxt,prectype%mltype) ! additive or multiplicative 2nd level prec
call psb_bcast(ictxt,prectype%smthpos) ! side: pre, post, both smoothing
call psb_bcast(ictxt,prectype%cmat) ! coarse mat