diff --git a/mlprec/impl/mld_c_hierarchy_bld.f90 b/mlprec/impl/mld_c_hierarchy_bld.f90 index 00219b23..305485b1 100644 --- a/mlprec/impl/mld_c_hierarchy_bld.f90 +++ b/mlprec/impl/mld_c_hierarchy_bld.f90 @@ -307,9 +307,6 @@ subroutine mld_c_hierarchy_bld(a,desc_a,p,info) ! Check for early termination of aggregation loop. ! iaggsize = sum(nlaggr) - if (iaggsize <= casize) then - newsz = i - end if sizeratio = iaggsize if (i==2) then @@ -318,6 +315,9 @@ subroutine mld_c_hierarchy_bld(a,desc_a,p,info) sizeratio = sum(p%precv(i-1)%map%naggr)/sizeratio end if p%precv(i)%szratio = sizeratio + if (iaggsize <= casize) then + newsz = i + end if if (i>2) then if (sizeratio < mnaggratio) then @@ -350,7 +350,16 @@ subroutine mld_c_hierarchy_bld(a,desc_a,p,info) if (newsz > 0) then if (info == 0) p%precv(newsz)%parms = coarseparms if (info == 0) call restore_smoothers(p%precv(newsz),coarse_sm,coarse_sm2,info) - + if (newsz < i) then + ! + ! We are going back and revisit a previous leve; + ! recover the aggregation. + ! + ilaggr = p%precv(newsz)%map%iaggr + nlaggr = p%precv(newsz)%map%naggr + call p%precv(newsz)%tprol%clone(op_prol,info) + end if + 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) diff --git a/mlprec/impl/mld_d_hierarchy_bld.f90 b/mlprec/impl/mld_d_hierarchy_bld.f90 index f2290180..e36ed39b 100644 --- a/mlprec/impl/mld_d_hierarchy_bld.f90 +++ b/mlprec/impl/mld_d_hierarchy_bld.f90 @@ -307,9 +307,6 @@ subroutine mld_d_hierarchy_bld(a,desc_a,p,info) ! Check for early termination of aggregation loop. ! iaggsize = sum(nlaggr) - if (iaggsize <= casize) then - newsz = i - end if sizeratio = iaggsize if (i==2) then @@ -318,6 +315,9 @@ subroutine mld_d_hierarchy_bld(a,desc_a,p,info) sizeratio = sum(p%precv(i-1)%map%naggr)/sizeratio end if p%precv(i)%szratio = sizeratio + if (iaggsize <= casize) then + newsz = i + end if if (i>2) then if (sizeratio < mnaggratio) then @@ -350,7 +350,16 @@ subroutine mld_d_hierarchy_bld(a,desc_a,p,info) if (newsz > 0) then if (info == 0) p%precv(newsz)%parms = coarseparms if (info == 0) call restore_smoothers(p%precv(newsz),coarse_sm,coarse_sm2,info) - + if (newsz < i) then + ! + ! We are going back and revisit a previous leve; + ! recover the aggregation. + ! + ilaggr = p%precv(newsz)%map%iaggr + nlaggr = p%precv(newsz)%map%naggr + call p%precv(newsz)%tprol%clone(op_prol,info) + end if + 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) diff --git a/mlprec/impl/mld_daggrmat_smth_asb.f90 b/mlprec/impl/mld_daggrmat_smth_asb.f90 index 92f99446..167645b2 100644 --- a/mlprec/impl/mld_daggrmat_smth_asb.f90 +++ b/mlprec/impl/mld_daggrmat_smth_asb.f90 @@ -280,8 +280,6 @@ subroutine mld_daggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_rest ! Doing it this way means to consider diag(Af_i) ! ! - write(0,*) 'Entering spspmm - filter',acsrf%get_ncols(),ptilde%get_nrows(),& - & desc_a%get_local_rows(),desc_a%get_local_cols() call psb_spspmm(acsrf,ptilde,acsr1,info) if(info /= psb_success_) then call psb_errpush(psb_err_from_subroutine_,name,a_err='spspmm 1') @@ -315,8 +313,6 @@ subroutine mld_daggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_rest ! Doing it this way means to consider diag(A_i) ! ! - write(0,*) 'Entering spspmm - nofilter',acsr3%get_ncols(),ptilde%get_nrows(),& - & desc_a%get_local_rows(),desc_a%get_local_cols() call psb_spspmm(acsr3,ptilde,acsr1,info) if(info /= psb_success_) then call psb_errpush(psb_err_from_subroutine_,name,a_err='spspmm 1') diff --git a/mlprec/impl/mld_s_hierarchy_bld.f90 b/mlprec/impl/mld_s_hierarchy_bld.f90 index b23b281d..86aaaa5a 100644 --- a/mlprec/impl/mld_s_hierarchy_bld.f90 +++ b/mlprec/impl/mld_s_hierarchy_bld.f90 @@ -307,9 +307,6 @@ subroutine mld_s_hierarchy_bld(a,desc_a,p,info) ! Check for early termination of aggregation loop. ! iaggsize = sum(nlaggr) - if (iaggsize <= casize) then - newsz = i - end if sizeratio = iaggsize if (i==2) then @@ -318,6 +315,9 @@ subroutine mld_s_hierarchy_bld(a,desc_a,p,info) sizeratio = sum(p%precv(i-1)%map%naggr)/sizeratio end if p%precv(i)%szratio = sizeratio + if (iaggsize <= casize) then + newsz = i + end if if (i>2) then if (sizeratio < mnaggratio) then @@ -350,7 +350,16 @@ subroutine mld_s_hierarchy_bld(a,desc_a,p,info) if (newsz > 0) then if (info == 0) p%precv(newsz)%parms = coarseparms if (info == 0) call restore_smoothers(p%precv(newsz),coarse_sm,coarse_sm2,info) - + if (newsz < i) then + ! + ! We are going back and revisit a previous leve; + ! recover the aggregation. + ! + ilaggr = p%precv(newsz)%map%iaggr + nlaggr = p%precv(newsz)%map%naggr + call p%precv(newsz)%tprol%clone(op_prol,info) + end if + 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) diff --git a/mlprec/impl/mld_z_hierarchy_bld.f90 b/mlprec/impl/mld_z_hierarchy_bld.f90 index c3ac7f86..27a32df4 100644 --- a/mlprec/impl/mld_z_hierarchy_bld.f90 +++ b/mlprec/impl/mld_z_hierarchy_bld.f90 @@ -307,9 +307,6 @@ subroutine mld_z_hierarchy_bld(a,desc_a,p,info) ! Check for early termination of aggregation loop. ! iaggsize = sum(nlaggr) - if (iaggsize <= casize) then - newsz = i - end if sizeratio = iaggsize if (i==2) then @@ -318,6 +315,9 @@ subroutine mld_z_hierarchy_bld(a,desc_a,p,info) sizeratio = sum(p%precv(i-1)%map%naggr)/sizeratio end if p%precv(i)%szratio = sizeratio + if (iaggsize <= casize) then + newsz = i + end if if (i>2) then if (sizeratio < mnaggratio) then @@ -350,7 +350,16 @@ subroutine mld_z_hierarchy_bld(a,desc_a,p,info) if (newsz > 0) then if (info == 0) p%precv(newsz)%parms = coarseparms if (info == 0) call restore_smoothers(p%precv(newsz),coarse_sm,coarse_sm2,info) - + if (newsz < i) then + ! + ! We are going back and revisit a previous leve; + ! recover the aggregation. + ! + ilaggr = p%precv(newsz)%map%iaggr + nlaggr = p%precv(newsz)%map%naggr + call p%precv(newsz)%tprol%clone(op_prol,info) + end if + 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) diff --git a/mlprec/mld_c_onelev_mod.f90 b/mlprec/mld_c_onelev_mod.f90 index 078657e1..1ea41f4a 100644 --- a/mlprec/mld_c_onelev_mod.f90 +++ b/mlprec/mld_c_onelev_mod.f90 @@ -129,6 +129,7 @@ module mld_c_onelev_mod type(psb_desc_type) :: desc_ac type(psb_cspmat_type), pointer :: base_a => null() type(psb_desc_type), pointer :: base_desc => null() + type(psb_cspmat_type) :: tprol type(psb_clinmap_type) :: map real(psb_spk_) :: szratio contains @@ -396,6 +397,7 @@ contains val = 0 val = val + lv%desc_ac%sizeof() val = val + lv%ac%sizeof() + val = val + lv%tprol%sizeof() val = val + lv%map%sizeof() if (allocated(lv%sm)) val = val + lv%sm%sizeof() if (allocated(lv%sm2a)) val = val + lv%sm2a%sizeof() @@ -488,6 +490,7 @@ contains end if if (info == psb_success_) call lv%parms%clone(lvout%parms,info) if (info == psb_success_) call lv%ac%clone(lvout%ac,info) + if (info == psb_success_) call lv%tprol%clone(lvout%tprol,info) if (info == psb_success_) call lv%desc_ac%clone(lvout%desc_ac,info) if (info == psb_success_) call lv%map%clone(lvout%map,info) lvout%base_a => lv%base_a @@ -517,7 +520,8 @@ contains b%sm2 =>b%sm end if - if (info == psb_success_) call psb_move_alloc(lv%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(lv%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(lv%tprol,b%tprol,info) if (info == psb_success_) call psb_move_alloc(lv%desc_ac,b%desc_ac,info) if (info == psb_success_) call psb_move_alloc(lv%map,b%map,info) b%base_a => lv%base_a diff --git a/mlprec/mld_d_onelev_mod.f90 b/mlprec/mld_d_onelev_mod.f90 index caceed28..805a9f84 100644 --- a/mlprec/mld_d_onelev_mod.f90 +++ b/mlprec/mld_d_onelev_mod.f90 @@ -129,6 +129,7 @@ module mld_d_onelev_mod type(psb_desc_type) :: desc_ac type(psb_dspmat_type), pointer :: base_a => null() type(psb_desc_type), pointer :: base_desc => null() + type(psb_dspmat_type) :: tprol type(psb_dlinmap_type) :: map real(psb_dpk_) :: szratio contains @@ -396,6 +397,7 @@ contains val = 0 val = val + lv%desc_ac%sizeof() val = val + lv%ac%sizeof() + val = val + lv%tprol%sizeof() val = val + lv%map%sizeof() if (allocated(lv%sm)) val = val + lv%sm%sizeof() if (allocated(lv%sm2a)) val = val + lv%sm2a%sizeof() @@ -488,6 +490,7 @@ contains end if if (info == psb_success_) call lv%parms%clone(lvout%parms,info) if (info == psb_success_) call lv%ac%clone(lvout%ac,info) + if (info == psb_success_) call lv%tprol%clone(lvout%tprol,info) if (info == psb_success_) call lv%desc_ac%clone(lvout%desc_ac,info) if (info == psb_success_) call lv%map%clone(lvout%map,info) lvout%base_a => lv%base_a @@ -517,7 +520,8 @@ contains b%sm2 =>b%sm end if - if (info == psb_success_) call psb_move_alloc(lv%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(lv%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(lv%tprol,b%tprol,info) if (info == psb_success_) call psb_move_alloc(lv%desc_ac,b%desc_ac,info) if (info == psb_success_) call psb_move_alloc(lv%map,b%map,info) b%base_a => lv%base_a diff --git a/mlprec/mld_s_onelev_mod.f90 b/mlprec/mld_s_onelev_mod.f90 index 1181e556..6acd74b4 100644 --- a/mlprec/mld_s_onelev_mod.f90 +++ b/mlprec/mld_s_onelev_mod.f90 @@ -129,6 +129,7 @@ module mld_s_onelev_mod type(psb_desc_type) :: desc_ac type(psb_sspmat_type), pointer :: base_a => null() type(psb_desc_type), pointer :: base_desc => null() + type(psb_sspmat_type) :: tprol type(psb_slinmap_type) :: map real(psb_spk_) :: szratio contains @@ -396,6 +397,7 @@ contains val = 0 val = val + lv%desc_ac%sizeof() val = val + lv%ac%sizeof() + val = val + lv%tprol%sizeof() val = val + lv%map%sizeof() if (allocated(lv%sm)) val = val + lv%sm%sizeof() if (allocated(lv%sm2a)) val = val + lv%sm2a%sizeof() @@ -488,6 +490,7 @@ contains end if if (info == psb_success_) call lv%parms%clone(lvout%parms,info) if (info == psb_success_) call lv%ac%clone(lvout%ac,info) + if (info == psb_success_) call lv%tprol%clone(lvout%tprol,info) if (info == psb_success_) call lv%desc_ac%clone(lvout%desc_ac,info) if (info == psb_success_) call lv%map%clone(lvout%map,info) lvout%base_a => lv%base_a @@ -517,7 +520,8 @@ contains b%sm2 =>b%sm end if - if (info == psb_success_) call psb_move_alloc(lv%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(lv%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(lv%tprol,b%tprol,info) if (info == psb_success_) call psb_move_alloc(lv%desc_ac,b%desc_ac,info) if (info == psb_success_) call psb_move_alloc(lv%map,b%map,info) b%base_a => lv%base_a diff --git a/mlprec/mld_z_onelev_mod.f90 b/mlprec/mld_z_onelev_mod.f90 index d50cd398..8be2597c 100644 --- a/mlprec/mld_z_onelev_mod.f90 +++ b/mlprec/mld_z_onelev_mod.f90 @@ -129,6 +129,7 @@ module mld_z_onelev_mod type(psb_desc_type) :: desc_ac type(psb_zspmat_type), pointer :: base_a => null() type(psb_desc_type), pointer :: base_desc => null() + type(psb_zspmat_type) :: tprol type(psb_zlinmap_type) :: map real(psb_dpk_) :: szratio contains @@ -396,6 +397,7 @@ contains val = 0 val = val + lv%desc_ac%sizeof() val = val + lv%ac%sizeof() + val = val + lv%tprol%sizeof() val = val + lv%map%sizeof() if (allocated(lv%sm)) val = val + lv%sm%sizeof() if (allocated(lv%sm2a)) val = val + lv%sm2a%sizeof() @@ -488,6 +490,7 @@ contains end if if (info == psb_success_) call lv%parms%clone(lvout%parms,info) if (info == psb_success_) call lv%ac%clone(lvout%ac,info) + if (info == psb_success_) call lv%tprol%clone(lvout%tprol,info) if (info == psb_success_) call lv%desc_ac%clone(lvout%desc_ac,info) if (info == psb_success_) call lv%map%clone(lvout%map,info) lvout%base_a => lv%base_a @@ -517,7 +520,8 @@ contains b%sm2 =>b%sm end if - if (info == psb_success_) call psb_move_alloc(lv%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(lv%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(lv%tprol,b%tprol,info) if (info == psb_success_) call psb_move_alloc(lv%desc_ac,b%desc_ac,info) if (info == psb_success_) call psb_move_alloc(lv%map,b%map,info) b%base_a => lv%base_a