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

mld2p4-extaggr:

Browse Source 
mlprec/impl/mld_c_hierarchy_bld.f90
 mlprec/impl/mld_ccprecset.F90
 mlprec/impl/mld_cprecset.F90
 mlprec/impl/mld_d_hierarchy_bld.f90
 mlprec/impl/mld_dcprecset.F90
 mlprec/impl/mld_dprecset.F90
 mlprec/impl/mld_s_hierarchy_bld.f90
 mlprec/impl/mld_scprecset.F90
 mlprec/impl/mld_sprecset.F90
 mlprec/impl/mld_z_hierarchy_bld.f90
 mlprec/impl/mld_zcprecset.F90
 mlprec/impl/mld_zprecset.F90

Added back SCALE in precset.
Fixed hierarchy buildup.
stopcriterion
Salvatore Filippone 8 years ago
parent 84562c67f4
commit 3e0040218a
12 changed files with 388 additions and 156 deletions
Show Stats Download Patch File Download Diff File

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

@ -82,7 +82,7 @@ subroutine mld_c_hierarchy_bld(a,desc_a,p,info)
integer(psb_ipk_) :: ictxt, me,np integer(psb_ipk_) :: ictxt, me,np
integer(psb_ipk_) :: err,i,k, err_act, iszv, newsz, casize, nplevs, mxplevs, iaggsize integer(psb_ipk_) :: err,i,k, err_act, iszv, newsz, casize, nplevs, mxplevs, iaggsize
real(psb_spk_) :: mnaggratio, sizeratio real(psb_spk_) :: mnaggratio, sizeratio
class(mld_c_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm class(mld_c_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm, base_sm2, med_sm2, coarse_sm2
type(mld_sml_parms) :: baseparms, medparms, coarseparms type(mld_sml_parms) :: baseparms, medparms, coarseparms
integer(psb_ipk_), allocatable :: ilaggr(:), nlaggr(:) integer(psb_ipk_), allocatable :: ilaggr(:), nlaggr(:)
type(psb_cspmat_type) :: op_prol type(psb_cspmat_type) :: op_prol
@ -215,50 +215,40 @@ subroutine mld_c_hierarchy_bld(a,desc_a,p,info)
baseparms = p%precv(1)%parms baseparms = p%precv(1)%parms
medparms = p%precv(2)%parms medparms = p%precv(2)%parms
allocate(coarse_sm, source=p%precv(iszv)%sm,stat=info) call save_smoothers(p%precv(iszv),coarse_sm,coarse_sm2,info)
if (info == psb_success_) & if (info == 0) call save_smoothers(p%precv(2),med_sm,med_sm2,info)
& allocate(med_sm, source=p%precv(2)%sm,stat=info) if (info == 0) call save_smoothers(p%precv(1),base_sm,base_sm2,info)
if (info == psb_success_) &
& allocate(base_sm, source=p%precv(1)%sm,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
write(0,*) 'Error in saving smoothers',info write(0,*) 'Error in saving smoothers',info
call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.')
goto 9999 goto 9999
end if end if
! !
! First set desired number of levels ! First set desired number of levels
! !
if (iszv /= nplevs) then if (iszv /= nplevs) then
allocate(tprecv(nplevs),stat=info) allocate(tprecv(nplevs),stat=info)
if (info /= psb_success_) then ! First all existing levels
call psb_errpush(psb_err_from_subroutine_,name,& if (info == 0) tprecv(1)%parms = baseparms
& a_err='prec reallocation') if (info == 0) call restore_smoothers(tprecv(1),p%precv(1)%sm,p%precv(1)%sm2a,info)
goto 9999 do i=2, min(iszv,nplevs) - 1
endif if (info == 0) tprecv(i)%parms = medparms
tprecv(1)%parms = baseparms if (info == 0) call restore_smoothers(tprecv(i),p%precv(i)%sm,p%precv(i)%sm2a,info)
allocate(tprecv(1)%sm,source=base_sm,stat=info) end do
if (info /= psb_success_) then ! Further intermediates, if any
call psb_errpush(psb_err_from_subroutine_,name,& do i=iszv-1, nplevs - 1
& a_err='prec reallocation') if (info == 0) tprecv(i)%parms = medparms
goto 9999 if (info == 0) call restore_smoothers(tprecv(i),med_sm,med_sm2,info)
endif
do i=2,nplevs-1
tprecv(i)%parms = medparms
allocate(tprecv(i)%sm,source=med_sm,stat=info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,&
& a_err='prec reallocation')
goto 9999
endif
end do end do
tprecv(nplevs)%parms = coarseparms ! Then coarse
allocate(tprecv(nplevs)%sm,source=coarse_sm,stat=info) if (info == 0) tprecv(nplevs)%parms = coarseparms
if (info == 0) call restore_smoothers(tprecv(nplevs),coarse_sm,coarse_sm2,info)
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,& call psb_errpush(psb_err_from_subroutine_,name,&
& a_err='prec reallocation') & a_err='prec reallocation')
goto 9999 goto 9999
endif endif
do i=1,iszv do i=1,iszv
call p%precv(i)%free(info) call p%precv(i)%free(info)
end do end do
@ -330,7 +320,7 @@ subroutine mld_c_hierarchy_bld(a,desc_a,p,info)
! !
! We are not gaining ! We are not gaining
! !
newsz = newsz-1 newsz = i-1
end if end if
end if end if
@ -350,20 +340,25 @@ subroutine mld_c_hierarchy_bld(a,desc_a,p,info)
end if end if
end if end if
call psb_bcast(ictxt,newsz) call psb_bcast(ictxt,newsz)
if (newsz > 0) & if (newsz > 0) then
& call p%precv(i)%parms%get_coarse(p%precv(iszv)%parms) if (info == 0) p%precv(newsz)%parms = coarseparms
if (info == 0) call restore_smoothers(p%precv(newsz),coarse_sm,coarse_sm2,info)
if (info == psb_success_) call mld_lev_mat_asb(p%precv(newsz),&
& p%precv(newsz-1)%base_a,p%precv(newsz-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info)
exit array_build_loop
else
if (info == psb_success_) call mld_lev_mat_asb(p%precv(i),& if (info == psb_success_) call mld_lev_mat_asb(p%precv(i),&
& p%precv(i-1)%base_a,p%precv(i-1)%base_desc,& & p%precv(i-1)%base_a,p%precv(i-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info) & ilaggr,nlaggr,op_prol,info)
end if
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,& call psb_errpush(psb_err_internal_error_,name,&
& a_err='Map build') & a_err='Map build')
goto 9999 goto 9999
endif endif
if (newsz > 0) exit array_build_loop
end do array_build_loop end do array_build_loop
if (newsz > 0) then if (newsz > 0) then
@ -399,7 +394,6 @@ subroutine mld_c_hierarchy_bld(a,desc_a,p,info)
end do end do
end if end if
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,& call psb_errpush(psb_err_internal_error_,name,&
& a_err='Internal hierarchy build' ) & a_err='Internal hierarchy build' )
@ -419,4 +413,58 @@ subroutine mld_c_hierarchy_bld(a,desc_a,p,info)
return return
contains
subroutine save_smoothers(level,save1, save2,info)
type(mld_c_onelev_type), intent(in) :: level
class(mld_c_base_smoother_type), allocatable , intent(inout) :: save1, save2
integer(psb_ipk_), intent(out) :: info
info = 0
if (allocated(save1)) then
call save1%free(info)
if (info == 0) deallocate(save1,stat=info)
if (info /= 0) return
end if
if (allocated(save2)) then
call save2%free(info)
if (info == 0) deallocate(save2,stat=info)
if (info /= 0) return
end if
allocate(save1, source=level%sm,stat=info)
if ((info == 0).and.allocated(level%sm2a)) allocate(save2, source=level%sm2a,stat=info)
return
end subroutine save_smoothers
subroutine restore_smoothers(level,save1, save2,info)
type(mld_c_onelev_type), intent(inout), target :: level
class(mld_c_base_smoother_type), allocatable, intent(in) :: save1, save2
integer(psb_ipk_), intent(out) :: info
info = 0
if (allocated(level%sm)) then
if (info == 0) call level%sm%free(info)
if (info == 0) deallocate(level%sm,stat=info)
end if
if (allocated(save1)) then
if (info == 0) allocate(level%sm,source=save1,stat=info)
end if
if (info /= 0) return
if (allocated(level%sm2a)) then
if (info == 0) call level%sm2a%free(info)
if (info == 0) deallocate(level%sm2a,stat=info)
end if
if (allocated(save2)) then
if (info == 0) allocate(level%sm2a,source=save2,stat=info)
if (info == 0) level%sm2 => level%sm2a
else
if (allocated(level%sm)) level%sm2 => level%sm
end if
return
end subroutine restore_smoothers
end subroutine mld_c_hierarchy_bld end subroutine mld_c_hierarchy_bld

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

@ -519,6 +519,11 @@ subroutine mld_ccprecsetr(p,what,val,info,ilev,pos)
ilev_=nlev_ ilev_=nlev_
call p%precv(ilev_)%set('SUB_ILUTHRS',val,info,pos=pos) call p%precv(ilev_)%set('SUB_ILUTHRS',val,info,pos=pos)
case('AGGR_SCALE')
do ilev_ = 2, nlev_
call p%precv(ilev_)%set('AGGR_SCALE',val,info,pos=pos)
end do
case('AGGR_THRESH') case('AGGR_THRESH')
thr = val thr = val
do ilev_ = 2, nlev_ do ilev_ = 2, nlev_

5
mlprec/impl/mld_cprecset.F90
Unescape Escape View File

@ -621,6 +621,11 @@ subroutine mld_cprecsetr(p,what,val,info,ilev,pos)
ilev_=nlev_ ilev_=nlev_
call p%precv(ilev_)%set(mld_sub_iluthrs_,val,info,pos=pos) call p%precv(ilev_)%set(mld_sub_iluthrs_,val,info,pos=pos)
case(mld_aggr_scale_)
do ilev_ = 2, nlev_
call p%precv(ilev_)%set(mld_aggr_scale_,val,info,pos=pos)
end do
case(mld_aggr_thresh_) case(mld_aggr_thresh_)
thr = val thr = val
do ilev_ = 2, nlev_ do ilev_ = 2, nlev_

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

@ -82,7 +82,7 @@ subroutine mld_d_hierarchy_bld(a,desc_a,p,info)
integer(psb_ipk_) :: ictxt, me,np integer(psb_ipk_) :: ictxt, me,np
integer(psb_ipk_) :: err,i,k, err_act, iszv, newsz, casize, nplevs, mxplevs, iaggsize integer(psb_ipk_) :: err,i,k, err_act, iszv, newsz, casize, nplevs, mxplevs, iaggsize
real(psb_dpk_) :: mnaggratio, sizeratio real(psb_dpk_) :: mnaggratio, sizeratio
class(mld_d_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm class(mld_d_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm, base_sm2, med_sm2, coarse_sm2
type(mld_dml_parms) :: baseparms, medparms, coarseparms type(mld_dml_parms) :: baseparms, medparms, coarseparms
integer(psb_ipk_), allocatable :: ilaggr(:), nlaggr(:) integer(psb_ipk_), allocatable :: ilaggr(:), nlaggr(:)
type(psb_dspmat_type) :: op_prol type(psb_dspmat_type) :: op_prol
@ -215,50 +215,40 @@ subroutine mld_d_hierarchy_bld(a,desc_a,p,info)
baseparms = p%precv(1)%parms baseparms = p%precv(1)%parms
medparms = p%precv(2)%parms medparms = p%precv(2)%parms
allocate(coarse_sm, source=p%precv(iszv)%sm,stat=info) call save_smoothers(p%precv(iszv),coarse_sm,coarse_sm2,info)
if (info == psb_success_) & if (info == 0) call save_smoothers(p%precv(2),med_sm,med_sm2,info)
& allocate(med_sm, source=p%precv(2)%sm,stat=info) if (info == 0) call save_smoothers(p%precv(1),base_sm,base_sm2,info)
if (info == psb_success_) &
& allocate(base_sm, source=p%precv(1)%sm,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
write(0,*) 'Error in saving smoothers',info write(0,*) 'Error in saving smoothers',info
call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.')
goto 9999 goto 9999
end if end if
! !
! First set desired number of levels ! First set desired number of levels
! !
if (iszv /= nplevs) then if (iszv /= nplevs) then
allocate(tprecv(nplevs),stat=info) allocate(tprecv(nplevs),stat=info)
if (info /= psb_success_) then ! First all existing levels
call psb_errpush(psb_err_from_subroutine_,name,& if (info == 0) tprecv(1)%parms = baseparms
& a_err='prec reallocation') if (info == 0) call restore_smoothers(tprecv(1),p%precv(1)%sm,p%precv(1)%sm2a,info)
goto 9999 do i=2, min(iszv,nplevs) - 1
endif if (info == 0) tprecv(i)%parms = medparms
tprecv(1)%parms = baseparms if (info == 0) call restore_smoothers(tprecv(i),p%precv(i)%sm,p%precv(i)%sm2a,info)
allocate(tprecv(1)%sm,source=base_sm,stat=info) end do
if (info /= psb_success_) then ! Further intermediates, if any
call psb_errpush(psb_err_from_subroutine_,name,& do i=iszv-1, nplevs - 1
& a_err='prec reallocation') if (info == 0) tprecv(i)%parms = medparms
goto 9999 if (info == 0) call restore_smoothers(tprecv(i),med_sm,med_sm2,info)
endif
do i=2,nplevs-1
tprecv(i)%parms = medparms
allocate(tprecv(i)%sm,source=med_sm,stat=info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,&
& a_err='prec reallocation')
goto 9999
endif
end do end do
tprecv(nplevs)%parms = coarseparms ! Then coarse
allocate(tprecv(nplevs)%sm,source=coarse_sm,stat=info) if (info == 0) tprecv(nplevs)%parms = coarseparms
if (info == 0) call restore_smoothers(tprecv(nplevs),coarse_sm,coarse_sm2,info)
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,& call psb_errpush(psb_err_from_subroutine_,name,&
& a_err='prec reallocation') & a_err='prec reallocation')
goto 9999 goto 9999
endif endif
do i=1,iszv do i=1,iszv
call p%precv(i)%free(info) call p%precv(i)%free(info)
end do end do
@ -330,7 +320,7 @@ subroutine mld_d_hierarchy_bld(a,desc_a,p,info)
! !
! We are not gaining ! We are not gaining
! !
newsz = newsz-1 newsz = i-1
end if end if
end if end if
@ -350,20 +340,25 @@ subroutine mld_d_hierarchy_bld(a,desc_a,p,info)
end if end if
end if end if
call psb_bcast(ictxt,newsz) call psb_bcast(ictxt,newsz)
if (newsz > 0) & if (newsz > 0) then
& call p%precv(i)%parms%get_coarse(p%precv(iszv)%parms) if (info == 0) p%precv(newsz)%parms = coarseparms
if (info == 0) call restore_smoothers(p%precv(newsz),coarse_sm,coarse_sm2,info)
if (info == psb_success_) call mld_lev_mat_asb(p%precv(newsz),&
& p%precv(newsz-1)%base_a,p%precv(newsz-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info)
exit array_build_loop
else
if (info == psb_success_) call mld_lev_mat_asb(p%precv(i),& if (info == psb_success_) call mld_lev_mat_asb(p%precv(i),&
& p%precv(i-1)%base_a,p%precv(i-1)%base_desc,& & p%precv(i-1)%base_a,p%precv(i-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info) & ilaggr,nlaggr,op_prol,info)
end if
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,& call psb_errpush(psb_err_internal_error_,name,&
& a_err='Map build') & a_err='Map build')
goto 9999 goto 9999
endif endif
if (newsz > 0) exit array_build_loop
end do array_build_loop end do array_build_loop
if (newsz > 0) then if (newsz > 0) then
@ -399,7 +394,6 @@ subroutine mld_d_hierarchy_bld(a,desc_a,p,info)
end do end do
end if end if
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,& call psb_errpush(psb_err_internal_error_,name,&
& a_err='Internal hierarchy build' ) & a_err='Internal hierarchy build' )
@ -419,4 +413,58 @@ subroutine mld_d_hierarchy_bld(a,desc_a,p,info)
return return
contains
subroutine save_smoothers(level,save1, save2,info)
type(mld_d_onelev_type), intent(in) :: level
class(mld_d_base_smoother_type), allocatable , intent(inout) :: save1, save2
integer(psb_ipk_), intent(out) :: info
info = 0
if (allocated(save1)) then
call save1%free(info)
if (info == 0) deallocate(save1,stat=info)
if (info /= 0) return
end if
if (allocated(save2)) then
call save2%free(info)
if (info == 0) deallocate(save2,stat=info)
if (info /= 0) return
end if
allocate(save1, source=level%sm,stat=info)
if ((info == 0).and.allocated(level%sm2a)) allocate(save2, source=level%sm2a,stat=info)
return
end subroutine save_smoothers
subroutine restore_smoothers(level,save1, save2,info)
type(mld_d_onelev_type), intent(inout), target :: level
class(mld_d_base_smoother_type), allocatable, intent(in) :: save1, save2
integer(psb_ipk_), intent(out) :: info
info = 0
if (allocated(level%sm)) then
if (info == 0) call level%sm%free(info)
if (info == 0) deallocate(level%sm,stat=info)
end if
if (allocated(save1)) then
if (info == 0) allocate(level%sm,source=save1,stat=info)
end if
if (info /= 0) return
if (allocated(level%sm2a)) then
if (info == 0) call level%sm2a%free(info)
if (info == 0) deallocate(level%sm2a,stat=info)
end if
if (allocated(save2)) then
if (info == 0) allocate(level%sm2a,source=save2,stat=info)
if (info == 0) level%sm2 => level%sm2a
else
if (allocated(level%sm)) level%sm2 => level%sm
end if
return
end subroutine restore_smoothers
end subroutine mld_d_hierarchy_bld end subroutine mld_d_hierarchy_bld

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

@ -529,6 +529,11 @@ subroutine mld_dcprecsetr(p,what,val,info,ilev,pos)
ilev_=nlev_ ilev_=nlev_
call p%precv(ilev_)%set('SUB_ILUTHRS',val,info,pos=pos) call p%precv(ilev_)%set('SUB_ILUTHRS',val,info,pos=pos)
case('AGGR_SCALE')
do ilev_ = 2, nlev_
call p%precv(ilev_)%set('AGGR_SCALE',val,info,pos=pos)
end do
case('AGGR_THRESH') case('AGGR_THRESH')
thr = val thr = val
do ilev_ = 2, nlev_ do ilev_ = 2, nlev_

5
mlprec/impl/mld_dprecset.F90
Unescape Escape View File

@ -631,6 +631,11 @@ subroutine mld_dprecsetr(p,what,val,info,ilev,pos)
ilev_=nlev_ ilev_=nlev_
call p%precv(ilev_)%set(mld_sub_iluthrs_,val,info,pos=pos) call p%precv(ilev_)%set(mld_sub_iluthrs_,val,info,pos=pos)
case(mld_aggr_scale_)
do ilev_ = 2, nlev_
call p%precv(ilev_)%set(mld_aggr_scale_,val,info,pos=pos)
end do
case(mld_aggr_thresh_) case(mld_aggr_thresh_)
thr = val thr = val
do ilev_ = 2, nlev_ do ilev_ = 2, nlev_

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

@ -82,7 +82,7 @@ subroutine mld_s_hierarchy_bld(a,desc_a,p,info)
integer(psb_ipk_) :: ictxt, me,np integer(psb_ipk_) :: ictxt, me,np
integer(psb_ipk_) :: err,i,k, err_act, iszv, newsz, casize, nplevs, mxplevs, iaggsize integer(psb_ipk_) :: err,i,k, err_act, iszv, newsz, casize, nplevs, mxplevs, iaggsize
real(psb_spk_) :: mnaggratio, sizeratio real(psb_spk_) :: mnaggratio, sizeratio
class(mld_s_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm class(mld_s_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm, base_sm2, med_sm2, coarse_sm2
type(mld_sml_parms) :: baseparms, medparms, coarseparms type(mld_sml_parms) :: baseparms, medparms, coarseparms
integer(psb_ipk_), allocatable :: ilaggr(:), nlaggr(:) integer(psb_ipk_), allocatable :: ilaggr(:), nlaggr(:)
type(psb_sspmat_type) :: op_prol type(psb_sspmat_type) :: op_prol
@ -215,50 +215,40 @@ subroutine mld_s_hierarchy_bld(a,desc_a,p,info)
baseparms = p%precv(1)%parms baseparms = p%precv(1)%parms
medparms = p%precv(2)%parms medparms = p%precv(2)%parms
allocate(coarse_sm, source=p%precv(iszv)%sm,stat=info) call save_smoothers(p%precv(iszv),coarse_sm,coarse_sm2,info)
if (info == psb_success_) & if (info == 0) call save_smoothers(p%precv(2),med_sm,med_sm2,info)
& allocate(med_sm, source=p%precv(2)%sm,stat=info) if (info == 0) call save_smoothers(p%precv(1),base_sm,base_sm2,info)
if (info == psb_success_) &
& allocate(base_sm, source=p%precv(1)%sm,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
write(0,*) 'Error in saving smoothers',info write(0,*) 'Error in saving smoothers',info
call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.')
goto 9999 goto 9999
end if end if
! !
! First set desired number of levels ! First set desired number of levels
! !
if (iszv /= nplevs) then if (iszv /= nplevs) then
allocate(tprecv(nplevs),stat=info) allocate(tprecv(nplevs),stat=info)
if (info /= psb_success_) then ! First all existing levels
call psb_errpush(psb_err_from_subroutine_,name,& if (info == 0) tprecv(1)%parms = baseparms
& a_err='prec reallocation') if (info == 0) call restore_smoothers(tprecv(1),p%precv(1)%sm,p%precv(1)%sm2a,info)
goto 9999 do i=2, min(iszv,nplevs) - 1
endif if (info == 0) tprecv(i)%parms = medparms
tprecv(1)%parms = baseparms if (info == 0) call restore_smoothers(tprecv(i),p%precv(i)%sm,p%precv(i)%sm2a,info)
allocate(tprecv(1)%sm,source=base_sm,stat=info) end do
if (info /= psb_success_) then ! Further intermediates, if any
call psb_errpush(psb_err_from_subroutine_,name,& do i=iszv-1, nplevs - 1
& a_err='prec reallocation') if (info == 0) tprecv(i)%parms = medparms
goto 9999 if (info == 0) call restore_smoothers(tprecv(i),med_sm,med_sm2,info)
endif
do i=2,nplevs-1
tprecv(i)%parms = medparms
allocate(tprecv(i)%sm,source=med_sm,stat=info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,&
& a_err='prec reallocation')
goto 9999
endif
end do end do
tprecv(nplevs)%parms = coarseparms ! Then coarse
allocate(tprecv(nplevs)%sm,source=coarse_sm,stat=info) if (info == 0) tprecv(nplevs)%parms = coarseparms
if (info == 0) call restore_smoothers(tprecv(nplevs),coarse_sm,coarse_sm2,info)
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,& call psb_errpush(psb_err_from_subroutine_,name,&
& a_err='prec reallocation') & a_err='prec reallocation')
goto 9999 goto 9999
endif endif
do i=1,iszv do i=1,iszv
call p%precv(i)%free(info) call p%precv(i)%free(info)
end do end do
@ -330,7 +320,7 @@ subroutine mld_s_hierarchy_bld(a,desc_a,p,info)
! !
! We are not gaining ! We are not gaining
! !
newsz = newsz-1 newsz = i-1
end if end if
end if end if
@ -350,20 +340,25 @@ subroutine mld_s_hierarchy_bld(a,desc_a,p,info)
end if end if
end if end if
call psb_bcast(ictxt,newsz) call psb_bcast(ictxt,newsz)
if (newsz > 0) & if (newsz > 0) then
& call p%precv(i)%parms%get_coarse(p%precv(iszv)%parms) if (info == 0) p%precv(newsz)%parms = coarseparms
if (info == 0) call restore_smoothers(p%precv(newsz),coarse_sm,coarse_sm2,info)
if (info == psb_success_) call mld_lev_mat_asb(p%precv(newsz),&
& p%precv(newsz-1)%base_a,p%precv(newsz-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info)
exit array_build_loop
else
if (info == psb_success_) call mld_lev_mat_asb(p%precv(i),& if (info == psb_success_) call mld_lev_mat_asb(p%precv(i),&
& p%precv(i-1)%base_a,p%precv(i-1)%base_desc,& & p%precv(i-1)%base_a,p%precv(i-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info) & ilaggr,nlaggr,op_prol,info)
end if
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,& call psb_errpush(psb_err_internal_error_,name,&
& a_err='Map build') & a_err='Map build')
goto 9999 goto 9999
endif endif
if (newsz > 0) exit array_build_loop
end do array_build_loop end do array_build_loop
if (newsz > 0) then if (newsz > 0) then
@ -399,7 +394,6 @@ subroutine mld_s_hierarchy_bld(a,desc_a,p,info)
end do end do
end if end if
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,& call psb_errpush(psb_err_internal_error_,name,&
& a_err='Internal hierarchy build' ) & a_err='Internal hierarchy build' )
@ -419,4 +413,58 @@ subroutine mld_s_hierarchy_bld(a,desc_a,p,info)
return return
contains
subroutine save_smoothers(level,save1, save2,info)
type(mld_s_onelev_type), intent(in) :: level
class(mld_s_base_smoother_type), allocatable , intent(inout) :: save1, save2
integer(psb_ipk_), intent(out) :: info
info = 0
if (allocated(save1)) then
call save1%free(info)
if (info == 0) deallocate(save1,stat=info)
if (info /= 0) return
end if
if (allocated(save2)) then
call save2%free(info)
if (info == 0) deallocate(save2,stat=info)
if (info /= 0) return
end if
allocate(save1, source=level%sm,stat=info)
if ((info == 0).and.allocated(level%sm2a)) allocate(save2, source=level%sm2a,stat=info)
return
end subroutine save_smoothers
subroutine restore_smoothers(level,save1, save2,info)
type(mld_s_onelev_type), intent(inout), target :: level
class(mld_s_base_smoother_type), allocatable, intent(in) :: save1, save2
integer(psb_ipk_), intent(out) :: info
info = 0
if (allocated(level%sm)) then
if (info == 0) call level%sm%free(info)
if (info == 0) deallocate(level%sm,stat=info)
end if
if (allocated(save1)) then
if (info == 0) allocate(level%sm,source=save1,stat=info)
end if
if (info /= 0) return
if (allocated(level%sm2a)) then
if (info == 0) call level%sm2a%free(info)
if (info == 0) deallocate(level%sm2a,stat=info)
end if
if (allocated(save2)) then
if (info == 0) allocate(level%sm2a,source=save2,stat=info)
if (info == 0) level%sm2 => level%sm2a
else
if (allocated(level%sm)) level%sm2 => level%sm
end if
return
end subroutine restore_smoothers
end subroutine mld_s_hierarchy_bld end subroutine mld_s_hierarchy_bld

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

@ -519,6 +519,11 @@ subroutine mld_scprecsetr(p,what,val,info,ilev,pos)
ilev_=nlev_ ilev_=nlev_
call p%precv(ilev_)%set('SUB_ILUTHRS',val,info,pos=pos) call p%precv(ilev_)%set('SUB_ILUTHRS',val,info,pos=pos)
case('AGGR_SCALE')
do ilev_ = 2, nlev_
call p%precv(ilev_)%set('AGGR_SCALE',val,info,pos=pos)
end do
case('AGGR_THRESH') case('AGGR_THRESH')
thr = val thr = val
do ilev_ = 2, nlev_ do ilev_ = 2, nlev_

5
mlprec/impl/mld_sprecset.F90
Unescape Escape View File

@ -621,6 +621,11 @@ subroutine mld_sprecsetr(p,what,val,info,ilev,pos)
ilev_=nlev_ ilev_=nlev_
call p%precv(ilev_)%set(mld_sub_iluthrs_,val,info,pos=pos) call p%precv(ilev_)%set(mld_sub_iluthrs_,val,info,pos=pos)
case(mld_aggr_scale_)
do ilev_ = 2, nlev_
call p%precv(ilev_)%set(mld_aggr_scale_,val,info,pos=pos)
end do
case(mld_aggr_thresh_) case(mld_aggr_thresh_)
thr = val thr = val
do ilev_ = 2, nlev_ do ilev_ = 2, nlev_

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

@ -82,7 +82,7 @@ subroutine mld_z_hierarchy_bld(a,desc_a,p,info)
integer(psb_ipk_) :: ictxt, me,np integer(psb_ipk_) :: ictxt, me,np
integer(psb_ipk_) :: err,i,k, err_act, iszv, newsz, casize, nplevs, mxplevs, iaggsize integer(psb_ipk_) :: err,i,k, err_act, iszv, newsz, casize, nplevs, mxplevs, iaggsize
real(psb_dpk_) :: mnaggratio, sizeratio real(psb_dpk_) :: mnaggratio, sizeratio
class(mld_z_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm class(mld_z_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm, base_sm2, med_sm2, coarse_sm2
type(mld_dml_parms) :: baseparms, medparms, coarseparms type(mld_dml_parms) :: baseparms, medparms, coarseparms
integer(psb_ipk_), allocatable :: ilaggr(:), nlaggr(:) integer(psb_ipk_), allocatable :: ilaggr(:), nlaggr(:)
type(psb_zspmat_type) :: op_prol type(psb_zspmat_type) :: op_prol
@ -215,50 +215,40 @@ subroutine mld_z_hierarchy_bld(a,desc_a,p,info)
baseparms = p%precv(1)%parms baseparms = p%precv(1)%parms
medparms = p%precv(2)%parms medparms = p%precv(2)%parms
allocate(coarse_sm, source=p%precv(iszv)%sm,stat=info) call save_smoothers(p%precv(iszv),coarse_sm,coarse_sm2,info)
if (info == psb_success_) & if (info == 0) call save_smoothers(p%precv(2),med_sm,med_sm2,info)
& allocate(med_sm, source=p%precv(2)%sm,stat=info) if (info == 0) call save_smoothers(p%precv(1),base_sm,base_sm2,info)
if (info == psb_success_) &
& allocate(base_sm, source=p%precv(1)%sm,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
write(0,*) 'Error in saving smoothers',info write(0,*) 'Error in saving smoothers',info
call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.')
goto 9999 goto 9999
end if end if
! !
! First set desired number of levels ! First set desired number of levels
! !
if (iszv /= nplevs) then if (iszv /= nplevs) then
allocate(tprecv(nplevs),stat=info) allocate(tprecv(nplevs),stat=info)
if (info /= psb_success_) then ! First all existing levels
call psb_errpush(psb_err_from_subroutine_,name,& if (info == 0) tprecv(1)%parms = baseparms
& a_err='prec reallocation') if (info == 0) call restore_smoothers(tprecv(1),p%precv(1)%sm,p%precv(1)%sm2a,info)
goto 9999 do i=2, min(iszv,nplevs) - 1
endif if (info == 0) tprecv(i)%parms = medparms
tprecv(1)%parms = baseparms if (info == 0) call restore_smoothers(tprecv(i),p%precv(i)%sm,p%precv(i)%sm2a,info)
allocate(tprecv(1)%sm,source=base_sm,stat=info) end do
if (info /= psb_success_) then ! Further intermediates, if any
call psb_errpush(psb_err_from_subroutine_,name,& do i=iszv-1, nplevs - 1
& a_err='prec reallocation') if (info == 0) tprecv(i)%parms = medparms
goto 9999 if (info == 0) call restore_smoothers(tprecv(i),med_sm,med_sm2,info)
endif
do i=2,nplevs-1
tprecv(i)%parms = medparms
allocate(tprecv(i)%sm,source=med_sm,stat=info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,&
& a_err='prec reallocation')
goto 9999
endif
end do end do
tprecv(nplevs)%parms = coarseparms ! Then coarse
allocate(tprecv(nplevs)%sm,source=coarse_sm,stat=info) if (info == 0) tprecv(nplevs)%parms = coarseparms
if (info == 0) call restore_smoothers(tprecv(nplevs),coarse_sm,coarse_sm2,info)
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,& call psb_errpush(psb_err_from_subroutine_,name,&
& a_err='prec reallocation') & a_err='prec reallocation')
goto 9999 goto 9999
endif endif
do i=1,iszv do i=1,iszv
call p%precv(i)%free(info) call p%precv(i)%free(info)
end do end do
@ -330,7 +320,7 @@ subroutine mld_z_hierarchy_bld(a,desc_a,p,info)
! !
! We are not gaining ! We are not gaining
! !
newsz = newsz-1 newsz = i-1
end if end if
end if end if
@ -350,20 +340,25 @@ subroutine mld_z_hierarchy_bld(a,desc_a,p,info)
end if end if
end if end if
call psb_bcast(ictxt,newsz) call psb_bcast(ictxt,newsz)
if (newsz > 0) & if (newsz > 0) then
& call p%precv(i)%parms%get_coarse(p%precv(iszv)%parms) if (info == 0) p%precv(newsz)%parms = coarseparms
if (info == 0) call restore_smoothers(p%precv(newsz),coarse_sm,coarse_sm2,info)
if (info == psb_success_) call mld_lev_mat_asb(p%precv(newsz),&
& p%precv(newsz-1)%base_a,p%precv(newsz-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info)
exit array_build_loop
else
if (info == psb_success_) call mld_lev_mat_asb(p%precv(i),& if (info == psb_success_) call mld_lev_mat_asb(p%precv(i),&
& p%precv(i-1)%base_a,p%precv(i-1)%base_desc,& & p%precv(i-1)%base_a,p%precv(i-1)%base_desc,&
& ilaggr,nlaggr,op_prol,info) & ilaggr,nlaggr,op_prol,info)
end if
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,& call psb_errpush(psb_err_internal_error_,name,&
& a_err='Map build') & a_err='Map build')
goto 9999 goto 9999
endif endif
if (newsz > 0) exit array_build_loop
end do array_build_loop end do array_build_loop
if (newsz > 0) then if (newsz > 0) then
@ -399,7 +394,6 @@ subroutine mld_z_hierarchy_bld(a,desc_a,p,info)
end do end do
end if end if
if (info /= psb_success_) then if (info /= psb_success_) then
call psb_errpush(psb_err_internal_error_,name,& call psb_errpush(psb_err_internal_error_,name,&
& a_err='Internal hierarchy build' ) & a_err='Internal hierarchy build' )
@ -419,4 +413,58 @@ subroutine mld_z_hierarchy_bld(a,desc_a,p,info)
return return
contains
subroutine save_smoothers(level,save1, save2,info)
type(mld_z_onelev_type), intent(in) :: level
class(mld_z_base_smoother_type), allocatable , intent(inout) :: save1, save2
integer(psb_ipk_), intent(out) :: info
info = 0
if (allocated(save1)) then
call save1%free(info)
if (info == 0) deallocate(save1,stat=info)
if (info /= 0) return
end if
if (allocated(save2)) then
call save2%free(info)
if (info == 0) deallocate(save2,stat=info)
if (info /= 0) return
end if
allocate(save1, source=level%sm,stat=info)
if ((info == 0).and.allocated(level%sm2a)) allocate(save2, source=level%sm2a,stat=info)
return
end subroutine save_smoothers
subroutine restore_smoothers(level,save1, save2,info)
type(mld_z_onelev_type), intent(inout), target :: level
class(mld_z_base_smoother_type), allocatable, intent(in) :: save1, save2
integer(psb_ipk_), intent(out) :: info
info = 0
if (allocated(level%sm)) then
if (info == 0) call level%sm%free(info)
if (info == 0) deallocate(level%sm,stat=info)
end if
if (allocated(save1)) then
if (info == 0) allocate(level%sm,source=save1,stat=info)
end if
if (info /= 0) return
if (allocated(level%sm2a)) then
if (info == 0) call level%sm2a%free(info)
if (info == 0) deallocate(level%sm2a,stat=info)
end if
if (allocated(save2)) then
if (info == 0) allocate(level%sm2a,source=save2,stat=info)
if (info == 0) level%sm2 => level%sm2a
else
if (allocated(level%sm)) level%sm2 => level%sm
end if
return
end subroutine restore_smoothers
end subroutine mld_z_hierarchy_bld end subroutine mld_z_hierarchy_bld

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

@ -529,6 +529,11 @@ subroutine mld_zcprecsetr(p,what,val,info,ilev,pos)
ilev_=nlev_ ilev_=nlev_
call p%precv(ilev_)%set('SUB_ILUTHRS',val,info,pos=pos) call p%precv(ilev_)%set('SUB_ILUTHRS',val,info,pos=pos)
case('AGGR_SCALE')
do ilev_ = 2, nlev_
call p%precv(ilev_)%set('AGGR_SCALE',val,info,pos=pos)
end do
case('AGGR_THRESH') case('AGGR_THRESH')
thr = val thr = val
do ilev_ = 2, nlev_ do ilev_ = 2, nlev_

5
mlprec/impl/mld_zprecset.F90
Unescape Escape View File

@ -631,6 +631,11 @@ subroutine mld_zprecsetr(p,what,val,info,ilev,pos)
ilev_=nlev_ ilev_=nlev_
call p%precv(ilev_)%set(mld_sub_iluthrs_,val,info,pos=pos) call p%precv(ilev_)%set(mld_sub_iluthrs_,val,info,pos=pos)
case(mld_aggr_scale_)
do ilev_ = 2, nlev_
call p%precv(ilev_)%set(mld_aggr_scale_,val,info,pos=pos)
end do
case(mld_aggr_thresh_) case(mld_aggr_thresh_)
thr = val thr = val
do ilev_ = 2, nlev_ do ilev_ = 2, nlev_

Write Preview
Loading…
Cancel
Save