From 96ddd0a2bb1ef41ffc6cb634d8422ec94f2155ec Mon Sep 17 00:00:00 2001 From: Salvatore Filippone Date: Mon, 30 Apr 2012 15:11:33 +0000 Subject: [PATCH] mld2p4-NewNL: mlprec/Makefile mlprec/impl/mld_c_onelev_impl.f90 mlprec/impl/mld_cmlprec_bld.f90 mlprec/impl/mld_cprecinit.F90 mlprec/impl/mld_cprecset.F90 mlprec/impl/mld_d_onelev_impl.f90 mlprec/impl/mld_dmlprec_bld.f90 mlprec/impl/mld_s_onelev_impl.f90 mlprec/impl/mld_smlprec_bld.f90 mlprec/impl/mld_sprecinit.F90 mlprec/impl/mld_sprecset.F90 mlprec/impl/mld_z_onelev_impl.f90 mlprec/impl/mld_zmlprec_bld.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_move_alloc_mod.f90 mlprec/mld_c_onelev_mod.f90 mlprec/mld_c_prec_mod.f90 mlprec/mld_c_prec_type.f90 mlprec/mld_d_inner_mod.f90 mlprec/mld_d_move_alloc_mod.f90 mlprec/mld_d_onelev_mod.f90 mlprec/mld_d_prec_mod.f90 mlprec/mld_d_prec_type.f90 mlprec/mld_s_inner_mod.f90 mlprec/mld_s_move_alloc_mod.f90 mlprec/mld_s_onelev_mod.f90 mlprec/mld_s_prec_mod.f90 mlprec/mld_s_prec_type.f90 mlprec/mld_z_inner_mod.f90 mlprec/mld_z_move_alloc_mod.f90 mlprec/mld_z_onelev_mod.f90 mlprec/mld_z_prec_mod.f90 mlprec/mld_z_prec_type.f90 tests/pdegen/ppde2d.f90 tests/pdegen/ppde3d.f90 tests/pdegen/runs/ppde.inp tests/pdegen/spde2d.f90 tests/pdegen/spde3d.f90 Working version of choice of levels with coarse space size. --- mlprec/Makefile | 32 ++-- mlprec/impl/mld_c_onelev_impl.f90 | 2 +- mlprec/impl/mld_cmlprec_bld.f90 | 267 ++++++++++++++++++++---------- mlprec/impl/mld_cprecinit.F90 | 1 + mlprec/impl/mld_cprecset.F90 | 5 + mlprec/impl/mld_d_onelev_impl.f90 | 2 +- mlprec/impl/mld_dmlprec_bld.f90 | 267 ++++++++++++++++++++---------- mlprec/impl/mld_s_onelev_impl.f90 | 2 +- mlprec/impl/mld_smlprec_bld.f90 | 267 ++++++++++++++++++++---------- mlprec/impl/mld_sprecinit.F90 | 1 + mlprec/impl/mld_sprecset.F90 | 5 + mlprec/impl/mld_z_onelev_impl.f90 | 2 +- mlprec/impl/mld_zmlprec_bld.f90 | 267 ++++++++++++++++++++---------- mlprec/impl/mld_zprecinit.F90 | 1 + mlprec/impl/mld_zprecset.F90 | 5 + mlprec/mld_base_prec_type.F90 | 37 +++++ mlprec/mld_c_inner_mod.f90 | 2 - mlprec/mld_c_move_alloc_mod.f90 | 1 + mlprec/mld_c_onelev_mod.f90 | 26 +++ mlprec/mld_c_prec_mod.f90 | 1 - mlprec/mld_c_prec_type.f90 | 34 +++- mlprec/mld_d_inner_mod.f90 | 2 - mlprec/mld_d_move_alloc_mod.f90 | 1 + mlprec/mld_d_onelev_mod.f90 | 26 +++ mlprec/mld_d_prec_mod.f90 | 1 - mlprec/mld_d_prec_type.f90 | 35 +++- mlprec/mld_s_inner_mod.f90 | 2 - mlprec/mld_s_move_alloc_mod.f90 | 1 + mlprec/mld_s_onelev_mod.f90 | 26 +++ mlprec/mld_s_prec_mod.f90 | 1 - mlprec/mld_s_prec_type.f90 | 34 +++- mlprec/mld_z_inner_mod.f90 | 2 - mlprec/mld_z_move_alloc_mod.f90 | 1 + mlprec/mld_z_onelev_mod.f90 | 26 +++ mlprec/mld_z_prec_mod.f90 | 1 - mlprec/mld_z_prec_type.f90 | 34 +++- tests/pdegen/ppde2d.f90 | 4 + tests/pdegen/ppde3d.f90 | 4 + tests/pdegen/runs/ppde.inp | 5 +- tests/pdegen/spde2d.f90 | 4 + tests/pdegen/spde3d.f90 | 4 + 41 files changed, 1046 insertions(+), 395 deletions(-) diff --git a/mlprec/Makefile b/mlprec/Makefile index 0650ba5d..e2cce546 100644 --- a/mlprec/Makefile +++ b/mlprec/Makefile @@ -7,22 +7,22 @@ HERE=. FINCLUDES=$(FMFLAG). $(FMFLAG)$(LIBDIR) $(FMFLAG)$(PSBINCDIR) $(FMFLAG)$(PSBLIBDIR) -DMODOBJS=mld_d_prec_type.o mld_d_prec_mod.o mld_d_move_alloc_mod.o mld_d_ilu_fact_mod.o \ +DMODOBJS=mld_d_prec_type.o mld_d_prec_mod.o mld_d_ilu_fact_mod.o \ mld_d_inner_mod.o mld_d_ilu_solver.o mld_d_diag_solver.o mld_d_jac_smoother.o mld_d_as_smoother.o \ mld_d_umf_solver.o mld_d_slu_solver.o mld_d_sludist_solver.o mld_d_id_solver.o\ mld_d_base_solver_mod.o mld_d_base_smoother_mod.o mld_d_onelev_mod.o -SMODOBJS=mld_s_prec_type.o mld_s_prec_mod.o mld_s_move_alloc_mod.o mld_s_ilu_fact_mod.o \ +SMODOBJS=mld_s_prec_type.o mld_s_prec_mod.o mld_s_ilu_fact_mod.o \ mld_s_inner_mod.o mld_s_ilu_solver.o mld_s_diag_solver.o mld_s_jac_smoother.o mld_s_as_smoother.o \ mld_s_slu_solver.o mld_s_sludist_solver.o mld_s_id_solver.o\ mld_s_base_solver_mod.o mld_s_base_smoother_mod.o mld_s_onelev_mod.o -ZMODOBJS=mld_z_prec_type.o mld_z_prec_mod.o mld_z_move_alloc_mod.o mld_z_ilu_fact_mod.o \ +ZMODOBJS=mld_z_prec_type.o mld_z_prec_mod.o mld_z_ilu_fact_mod.o \ mld_z_inner_mod.o mld_z_ilu_solver.o mld_z_diag_solver.o mld_z_jac_smoother.o mld_z_as_smoother.o \ mld_z_umf_solver.o mld_z_slu_solver.o mld_z_sludist_solver.o mld_z_id_solver.o\ mld_z_base_solver_mod.o mld_z_base_smoother_mod.o mld_z_onelev_mod.o -CMODOBJS=mld_c_prec_type.o mld_c_prec_mod.o mld_c_move_alloc_mod.o mld_c_ilu_fact_mod.o \ +CMODOBJS=mld_c_prec_type.o mld_c_prec_mod.o mld_c_ilu_fact_mod.o \ mld_c_inner_mod.o mld_c_ilu_solver.o mld_c_diag_solver.o mld_c_jac_smoother.o mld_c_as_smoother.o \ mld_c_slu_solver.o mld_c_sludist_solver.o mld_c_id_solver.o\ mld_c_base_solver_mod.o mld_c_base_smoother_mod.o mld_c_onelev_mod.o @@ -63,20 +63,15 @@ $(DINNEROBJS) $(DOUTEROBJS): $(DMODOBJS) $(CINNEROBJS) $(COUTEROBJS): $(CMODOBJS) $(ZINNEROBJS) $(ZOUTEROBJS): $(ZMODOBJS) -mld_s_inner_mod.o: mld_s_move_alloc_mod.o mld_s_prec_type.o -mld_d_inner_mod.o: mld_d_move_alloc_mod.o mld_d_prec_type.o -mld_c_inner_mod.o: mld_c_move_alloc_mod.o mld_c_prec_type.o -mld_z_inner_mod.o: mld_z_move_alloc_mod.o mld_z_prec_type.o +mld_s_inner_mod.o: mld_s_prec_type.o +mld_d_inner_mod.o: mld_d_prec_type.o +mld_c_inner_mod.o: mld_c_prec_type.o +mld_z_inner_mod.o: mld_z_prec_type.o -mld_s_move_alloc_mod.o: mld_s_prec_type.o -mld_d_move_alloc_mod.o: mld_d_prec_type.o -mld_c_move_alloc_mod.o: mld_c_prec_type.o -mld_z_move_alloc_mod.o: mld_z_prec_type.o - -mld_s_prec_mod.o: mld_s_move_alloc_mod.o -mld_d_prec_mod.o: mld_d_move_alloc_mod.o -mld_c_prec_mod.o: mld_c_move_alloc_mod.o -mld_z_prec_mod.o: mld_z_move_alloc_mod.o +mld_s_prec_mod.o: mld_s_prec_type.o +mld_d_prec_mod.o: mld_d_prec_type.o +mld_c_prec_mod.o: mld_c_prec_type.o +mld_z_prec_mod.o: mld_z_prec_type.o mld_s_prec_type.o: mld_s_onelev_mod.o @@ -97,9 +92,6 @@ mld_z_base_smoother_mod.o: mld_z_base_solver_mod.o mld_s_base_solver_mod.o mld_d_base_solver_mod.o mld_c_base_solver_mod.o mld_z_base_solver_mod.o: mld_base_prec_type.o - - - mld_d_id_solver.o mld_d_sludist_solver.o mld_d_slu_solver.o \ mld_d_umf_solver.o mld_d_diag_solver.o mld_d_ilu_solver.o: mld_d_base_solver_mod.o mld_d_prec_type.o diff --git a/mlprec/impl/mld_c_onelev_impl.f90 b/mlprec/impl/mld_c_onelev_impl.f90 index e5faa843..d25bd23f 100644 --- a/mlprec/impl/mld_c_onelev_impl.f90 +++ b/mlprec/impl/mld_c_onelev_impl.f90 @@ -154,7 +154,7 @@ subroutine mld_c_base_onelev_free(lv,info) & call lv%sm%free(info) call lv%ac%free() - if (psb_is_ok_desc(lv%desc_ac)) & + if (lv%desc_ac%is_ok()) & & call psb_cdfree(lv%desc_ac,info) call lv%map%free(info) diff --git a/mlprec/impl/mld_cmlprec_bld.f90 b/mlprec/impl/mld_cmlprec_bld.f90 index b3d6430a..a2fa21c9 100644 --- a/mlprec/impl/mld_cmlprec_bld.f90 +++ b/mlprec/impl/mld_cmlprec_bld.f90 @@ -93,11 +93,13 @@ subroutine mld_cmlprec_bld(a,desc_a,p,info,amold,vmold) ! Local Variables type(mld_cprec_type) :: t_prec - Integer :: err,i,k,ictxt, me,np, err_act, iszv, newsz + Integer :: err,i,k,ictxt, me,np, err_act, iszv, newsz, casize integer :: ipv(mld_ifpsz_), val integer :: int_err(5) character :: upd_ - type(mld_sml_parms) :: prm + class(mld_c_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm + type(mld_sml_parms) :: baseparms, medparms, coarseparms + type(mld_c_onelev_node), pointer :: head, tail, newnode, current integer :: debug_level, debug_unit character(len=20) :: name, ch_err @@ -145,12 +147,22 @@ subroutine mld_cmlprec_bld(a,desc_a,p,info,amold,vmold) ! Check to ensure all procs have the same ! newsz = -1 + casize = p%coarse_aggr_size iszv = size(p%precv) call psb_bcast(ictxt,iszv) - if (iszv /= size(p%precv)) then - info=psb_err_internal_error_ - call psb_errpush(info,name,a_err='Inconsistent size of precv') - goto 9999 + call psb_bcast(ictxt,casize) + if (casize > 0) then + if (casize /= p%coarse_aggr_size) then + info=psb_err_internal_error_ + call psb_errpush(info,name,a_err='Inconsistent coarse_aggr_size') + goto 9999 + end if + else + if (iszv /= size(p%precv)) then + info=psb_err_internal_error_ + call psb_errpush(info,name,a_err='Inconsistent size of precv') + goto 9999 + end if end if if (iszv <= 1) then @@ -162,7 +174,161 @@ subroutine mld_cmlprec_bld(a,desc_a,p,info,amold,vmold) endif - if (iszv > 1) then + + if (casize>0) then + ! + ! New strategy to build according to coarse size. + ! + coarseparms = p%precv(iszv)%parms + baseparms = p%precv(1)%parms + medparms = p%precv(2)%parms + + allocate(coarse_sm, source=p%precv(iszv)%sm,stat=info) + if (info == psb_success_) & + & allocate(med_sm, source=p%precv(2)%sm,stat=info) + if (info == psb_success_) & + & allocate(base_sm, source=p%precv(1)%sm,stat=info) + if (info /= psb_success_) then + write(0,*) 'Error in saving smoothers',info + call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') + goto 9999 + end if + ! + ! Replicated matrix should only ever happen at coarse level. + ! + call mld_check_def(baseparms%coarse_mat,'Coarse matrix',& + & mld_distr_mat_,is_distr_ml_coarse_mat) + ! + ! Now build a doubly linked list + ! + allocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='List start'); goto 9999 + end if + head => newnode + tail => newnode + newnode%item%base_a => a + newnode%item%base_desc => desc_a + newnode%item%parms = baseparms + newsz = 1 + current => head + list_build_loop: do + allocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='List start'); goto 9999 + end if + current%next => newnode + newnode%prev => current + newsz = newsz + 1 + newnode%item%parms = medparms + newnode%item%parms%aggr_thresh = current%item%parms%aggr_thresh/2 + call mld_coarse_bld(current%item%base_a, current%item%base_desc, & + & newnode%item,info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='build next level'); goto 9999 + end if + + if (newsz>2) then + if (all(current%item%map%naggr == newnode%item%map%naggr)) then + ! + ! We are not gaining anything + ! + newsz = newsz-1 + current%next => null() + call newnode%item%free(info) + if (info == psb_success_) deallocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Deallocate at list end'); goto 9999 + end if + end if + end if + + current => current%next + tail => current + if (sum(newnode%item%map%naggr) <= casize) then + ! + ! Target reached; but we may need to rebuild. + ! + exit list_build_loop + end if + end do list_build_loop + ! + ! At this point, we are at the list tail, + ! and it needs to be rebuilt in case the parms were + ! different. + ! + ! But the threshold has to be fixed before rebuliding + coarseparms%aggr_thresh = current%item%parms%aggr_thresh + current%item%parms = coarseparms + call mld_coarse_bld(current%prev%item%base_a,& + & current%prev%item%base_desc, & + & current%item,info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='build next level'); goto 9999 + end if + + ! + ! Ok, now allocate the output vector and fix items. + ! + do i=1,iszv + if (info == psb_success_) call p%precv(i)%free(info) + end do + if (info == psb_success_) deallocate(p%precv,stat=info) + if (info == psb_success_) allocate(p%precv(newsz),stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Reallocate precv'); goto 9999 + end if + newnode => head + do i=1, newsz + current => newnode + if (info == psb_success_) & + & call mld_move_alloc(current%item,p%precv(i),info) + if (info == psb_success_) then + if (i ==1) then + allocate(p%precv(i)%sm,source=base_sm,stat=info) + else if (i < newsz) then + allocate(p%precv(i)%sm,source=med_sm,stat=info) + else + allocate(p%precv(i)%sm,source=coarse_sm,stat=info) + end if + end if + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='list cpy'); goto 9999 + end if + if (i == 1) then + p%precv(i)%base_a => a + p%precv(i)%base_desc => desc_a + else + p%precv(i)%base_a => p%precv(i)%ac + p%precv(i)%base_desc => p%precv(i)%desc_ac + p%precv(i)%map%p_desc_X => p%precv(i-1)%base_desc + p%precv(i)%map%p_desc_Y => p%precv(i)%base_desc + end if + + newnode => current%next + deallocate(current) + end do + call base_sm%free(info) + if (info == psb_success_) call med_sm%free(info) + if (info == psb_success_) call coarse_sm%free(info) + if (info == psb_success_) deallocate(coarse_sm,med_sm,base_sm,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='final cleanup'); goto 9999 + end if + iszv = newsz + + else + ! + ! Default, oldstyle + ! ! ! Build the matrix and the transfer operators corresponding @@ -179,11 +345,6 @@ subroutine mld_cmlprec_bld(a,desc_a,p,info,amold,vmold) p%precv(1)%base_a => a p%precv(1)%base_desc => desc_a - if (info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') - goto 9999 - end if - do i=2, iszv ! @@ -201,11 +362,6 @@ subroutine mld_cmlprec_bld(a,desc_a,p,info,amold,vmold) ! call mld_check_def(p%precv(i)%parms%coarse_mat,'Coarse matrix',& & mld_distr_mat_,is_distr_ml_coarse_mat) - - else if (i == iszv) then - -!!$ call check_coarse_lev(p%precv(i)) - end if if (debug_level >= psb_debug_outer_) & @@ -277,9 +433,7 @@ subroutine mld_cmlprec_bld(a,desc_a,p,info,amold,vmold) p%precv(i)%map%p_desc_Y => p%precv(i)%base_desc end do - i = iszv - call check_coarse_lev(p%precv(i)) if (info == psb_success_) call mld_coarse_bld(p%precv(i-1)%base_a,& & p%precv(i-1)%base_desc, p%precv(i),info) if (info /= psb_success_) then @@ -289,6 +443,12 @@ subroutine mld_cmlprec_bld(a,desc_a,p,info,amold,vmold) end if end if + ! + ! The coarse space hierarchy has been build. + ! + ! Now do the preconditioner build. + ! + do i=1, iszv ! ! build the base preconditioner at level i @@ -316,10 +476,6 @@ subroutine mld_cmlprec_bld(a,desc_a,p,info,amold,vmold) goto 9999 end if - - ! - ! Test version for beginning of OO stuff. - ! call p%precv(i)%sm%build(p%precv(i)%base_a,p%precv(i)%base_desc,& & 'F',info,amold=amold,vmold=vmold) @@ -350,69 +506,4 @@ subroutine mld_cmlprec_bld(a,desc_a,p,info,amold,vmold) end if return -contains - - subroutine check_coarse_lev(prec) - type(mld_c_onelev_type) :: prec - - ! - ! At the coarsest level, check mld_coarse_solve_ - ! -!!$ val = prec%parms%coarse_solve -!!$ select case (val) -!!$ case(mld_jac_) -!!$ -!!$ if (prec%prec%iprcparm(mld_sub_solve_) /= mld_diag_scale_) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%prec%iprcparm(mld_sub_solve_) = mld_diag_scale_ -!!$ end if -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_jac_ -!!$ -!!$ case(mld_bjac_) -!!$ -!!$ if ((prec%prec%iprcparm(mld_sub_solve_) == mld_diag_scale_).or.& -!!$ & ( prec%prec%iprcparm(mld_smoother_type_) /= mld_bjac_)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$! !$#if defined(HAVE_UMF_) -!!$! !$ prec%prec%iprcparm(mld_sub_solve_) = mld_umf_ -!!$! !$#elif defined(HAVE_SLU_) -!!$! !$ prec%prec%iprcparm(mld_sub_solve_) = mld_slu_ -!!$! !$#else -!!$ prec%prec%iprcparm(mld_sub_solve_) = mld_ilu_n_ -!!$! !$#endif -!!$ end if -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ -!!$ case(mld_umf_, mld_slu_) -!!$ if ((prec%iprcparm(mld_coarse_mat_) /= mld_repl_mat_).or.& -!!$ & (prec%prec%iprcparm(mld_sub_solve_) /= val)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%iprcparm(mld_coarse_mat_) = mld_repl_mat_ -!!$ prec%prec%iprcparm(mld_sub_solve_) = val -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ end if -!!$ case(mld_sludist_) -!!$ if ((prec%iprcparm(mld_coarse_mat_) /= mld_distr_mat_).or.& -!!$ & (prec%prec%iprcparm(mld_sub_solve_) /= val)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%iprcparm(mld_coarse_mat_) = mld_distr_mat_ -!!$ prec%prec%iprcparm(mld_sub_solve_) = val -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ prec%prec%iprcparm(mld_smoother_sweeps_) = 1 -!!$ end if -!!$ end select - end subroutine check_coarse_lev - end subroutine mld_cmlprec_bld diff --git a/mlprec/impl/mld_cprecinit.F90 b/mlprec/impl/mld_cprecinit.F90 index 2890a83f..e486f6ba 100644 --- a/mlprec/impl/mld_cprecinit.F90 +++ b/mlprec/impl/mld_cprecinit.F90 @@ -125,6 +125,7 @@ subroutine mld_cprecinit(p,ptype,info,nlev) ! Do we want to do something? endif endif + p%coarse_aggr_size = -1 select case(psb_toupper(ptype(1:len_trim(ptype)))) case ('NOPREC','NONE') diff --git a/mlprec/impl/mld_cprecset.F90 b/mlprec/impl/mld_cprecset.F90 index 167f6bc4..3d11fad2 100644 --- a/mlprec/impl/mld_cprecset.F90 +++ b/mlprec/impl/mld_cprecset.F90 @@ -129,6 +129,11 @@ subroutine mld_cprecseti(p,what,val,info,ilev) return endif + if (what == mld_coarse_aggr_size_) then + p%coarse_aggr_size = max(val,-1) + return + end if + ! ! Set preconditioner parameters at level ilev. ! diff --git a/mlprec/impl/mld_d_onelev_impl.f90 b/mlprec/impl/mld_d_onelev_impl.f90 index 46d0ccab..c83ff1de 100644 --- a/mlprec/impl/mld_d_onelev_impl.f90 +++ b/mlprec/impl/mld_d_onelev_impl.f90 @@ -154,7 +154,7 @@ subroutine mld_d_base_onelev_free(lv,info) & call lv%sm%free(info) call lv%ac%free() - if (psb_is_ok_desc(lv%desc_ac)) & + if (lv%desc_ac%is_ok()) & & call psb_cdfree(lv%desc_ac,info) call lv%map%free(info) diff --git a/mlprec/impl/mld_dmlprec_bld.f90 b/mlprec/impl/mld_dmlprec_bld.f90 index 0051fe69..94481c66 100644 --- a/mlprec/impl/mld_dmlprec_bld.f90 +++ b/mlprec/impl/mld_dmlprec_bld.f90 @@ -93,11 +93,13 @@ subroutine mld_dmlprec_bld(a,desc_a,p,info,amold,vmold) ! Local Variables type(mld_dprec_type) :: t_prec - Integer :: err,i,k,ictxt, me,np, err_act, iszv, newsz + Integer :: err,i,k,ictxt, me,np, err_act, iszv, newsz, casize integer :: ipv(mld_ifpsz_), val integer :: int_err(5) character :: upd_ - type(mld_dml_parms) :: prm + class(mld_d_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm + type(mld_dml_parms) :: baseparms, medparms, coarseparms + type(mld_d_onelev_node), pointer :: head, tail, newnode, current integer :: debug_level, debug_unit character(len=20) :: name, ch_err @@ -145,12 +147,22 @@ subroutine mld_dmlprec_bld(a,desc_a,p,info,amold,vmold) ! Check to ensure all procs have the same ! newsz = -1 + casize = p%coarse_aggr_size iszv = size(p%precv) call psb_bcast(ictxt,iszv) - if (iszv /= size(p%precv)) then - info=psb_err_internal_error_ - call psb_errpush(info,name,a_err='Inconsistent size of precv') - goto 9999 + call psb_bcast(ictxt,casize) + if (casize > 0) then + if (casize /= p%coarse_aggr_size) then + info=psb_err_internal_error_ + call psb_errpush(info,name,a_err='Inconsistent coarse_aggr_size') + goto 9999 + end if + else + if (iszv /= size(p%precv)) then + info=psb_err_internal_error_ + call psb_errpush(info,name,a_err='Inconsistent size of precv') + goto 9999 + end if end if if (iszv <= 1) then @@ -162,7 +174,161 @@ subroutine mld_dmlprec_bld(a,desc_a,p,info,amold,vmold) endif - if (iszv > 1) then + + if (casize>0) then + ! + ! New strategy to build according to coarse size. + ! + coarseparms = p%precv(iszv)%parms + baseparms = p%precv(1)%parms + medparms = p%precv(2)%parms + + allocate(coarse_sm, source=p%precv(iszv)%sm,stat=info) + if (info == psb_success_) & + & allocate(med_sm, source=p%precv(2)%sm,stat=info) + if (info == psb_success_) & + & allocate(base_sm, source=p%precv(1)%sm,stat=info) + if (info /= psb_success_) then + write(0,*) 'Error in saving smoothers',info + call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') + goto 9999 + end if + ! + ! Replicated matrix should only ever happen at coarse level. + ! + call mld_check_def(baseparms%coarse_mat,'Coarse matrix',& + & mld_distr_mat_,is_distr_ml_coarse_mat) + ! + ! Now build a doubly linked list + ! + allocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='List start'); goto 9999 + end if + head => newnode + tail => newnode + newnode%item%base_a => a + newnode%item%base_desc => desc_a + newnode%item%parms = baseparms + newsz = 1 + current => head + list_build_loop: do + allocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='List start'); goto 9999 + end if + current%next => newnode + newnode%prev => current + newsz = newsz + 1 + newnode%item%parms = medparms + newnode%item%parms%aggr_thresh = current%item%parms%aggr_thresh/2 + call mld_coarse_bld(current%item%base_a, current%item%base_desc, & + & newnode%item,info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='build next level'); goto 9999 + end if + + if (newsz>2) then + if (all(current%item%map%naggr == newnode%item%map%naggr)) then + ! + ! We are not gaining anything + ! + newsz = newsz-1 + current%next => null() + call newnode%item%free(info) + if (info == psb_success_) deallocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Deallocate at list end'); goto 9999 + end if + end if + end if + + current => current%next + tail => current + if (sum(newnode%item%map%naggr) <= casize) then + ! + ! Target reached; but we may need to rebuild. + ! + exit list_build_loop + end if + end do list_build_loop + ! + ! At this point, we are at the list tail, + ! and it needs to be rebuilt in case the parms were + ! different. + ! + ! But the threshold has to be fixed before rebuliding + coarseparms%aggr_thresh = current%item%parms%aggr_thresh + current%item%parms = coarseparms + call mld_coarse_bld(current%prev%item%base_a,& + & current%prev%item%base_desc, & + & current%item,info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='build next level'); goto 9999 + end if + + ! + ! Ok, now allocate the output vector and fix items. + ! + do i=1,iszv + if (info == psb_success_) call p%precv(i)%free(info) + end do + if (info == psb_success_) deallocate(p%precv,stat=info) + if (info == psb_success_) allocate(p%precv(newsz),stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Reallocate precv'); goto 9999 + end if + newnode => head + do i=1, newsz + current => newnode + if (info == psb_success_) & + & call mld_move_alloc(current%item,p%precv(i),info) + if (info == psb_success_) then + if (i ==1) then + allocate(p%precv(i)%sm,source=base_sm,stat=info) + else if (i < newsz) then + allocate(p%precv(i)%sm,source=med_sm,stat=info) + else + allocate(p%precv(i)%sm,source=coarse_sm,stat=info) + end if + end if + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='list cpy'); goto 9999 + end if + if (i == 1) then + p%precv(i)%base_a => a + p%precv(i)%base_desc => desc_a + else + p%precv(i)%base_a => p%precv(i)%ac + p%precv(i)%base_desc => p%precv(i)%desc_ac + p%precv(i)%map%p_desc_X => p%precv(i-1)%base_desc + p%precv(i)%map%p_desc_Y => p%precv(i)%base_desc + end if + + newnode => current%next + deallocate(current) + end do + call base_sm%free(info) + if (info == psb_success_) call med_sm%free(info) + if (info == psb_success_) call coarse_sm%free(info) + if (info == psb_success_) deallocate(coarse_sm,med_sm,base_sm,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='final cleanup'); goto 9999 + end if + iszv = newsz + + else + ! + ! Default, oldstyle + ! ! ! Build the matrix and the transfer operators corresponding @@ -179,11 +345,6 @@ subroutine mld_dmlprec_bld(a,desc_a,p,info,amold,vmold) p%precv(1)%base_a => a p%precv(1)%base_desc => desc_a - if (info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') - goto 9999 - end if - do i=2, iszv ! @@ -201,11 +362,6 @@ subroutine mld_dmlprec_bld(a,desc_a,p,info,amold,vmold) ! call mld_check_def(p%precv(i)%parms%coarse_mat,'Coarse matrix',& & mld_distr_mat_,is_distr_ml_coarse_mat) - - else if (i == iszv) then - -!!$ call check_coarse_lev(p%precv(i)) - end if if (debug_level >= psb_debug_outer_) & @@ -277,9 +433,7 @@ subroutine mld_dmlprec_bld(a,desc_a,p,info,amold,vmold) p%precv(i)%map%p_desc_Y => p%precv(i)%base_desc end do - i = iszv - call check_coarse_lev(p%precv(i)) if (info == psb_success_) call mld_coarse_bld(p%precv(i-1)%base_a,& & p%precv(i-1)%base_desc, p%precv(i),info) if (info /= psb_success_) then @@ -289,6 +443,12 @@ subroutine mld_dmlprec_bld(a,desc_a,p,info,amold,vmold) end if end if + ! + ! The coarse space hierarchy has been build. + ! + ! Now do the preconditioner build. + ! + do i=1, iszv ! ! build the base preconditioner at level i @@ -316,10 +476,6 @@ subroutine mld_dmlprec_bld(a,desc_a,p,info,amold,vmold) goto 9999 end if - - ! - ! Test version for beginning of OO stuff. - ! call p%precv(i)%sm%build(p%precv(i)%base_a,p%precv(i)%base_desc,& & 'F',info,amold=amold,vmold=vmold) @@ -350,69 +506,4 @@ subroutine mld_dmlprec_bld(a,desc_a,p,info,amold,vmold) end if return -contains - - subroutine check_coarse_lev(prec) - type(mld_d_onelev_type) :: prec - - ! - ! At the coarsest level, check mld_coarse_solve_ - ! -!!$ val = prec%parms%coarse_solve -!!$ select case (val) -!!$ case(mld_jac_) -!!$ -!!$ if (prec%prec%iprcparm(mld_sub_solve_) /= mld_diag_scale_) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%prec%iprcparm(mld_sub_solve_) = mld_diag_scale_ -!!$ end if -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_jac_ -!!$ -!!$ case(mld_bjac_) -!!$ -!!$ if ((prec%prec%iprcparm(mld_sub_solve_) == mld_diag_scale_).or.& -!!$ & ( prec%prec%iprcparm(mld_smoother_type_) /= mld_bjac_)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$! !$#if defined(HAVE_UMF_) -!!$! !$ prec%prec%iprcparm(mld_sub_solve_) = mld_umf_ -!!$! !$#elif defined(HAVE_SLU_) -!!$! !$ prec%prec%iprcparm(mld_sub_solve_) = mld_slu_ -!!$! !$#else -!!$ prec%prec%iprcparm(mld_sub_solve_) = mld_ilu_n_ -!!$! !$#endif -!!$ end if -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ -!!$ case(mld_umf_, mld_slu_) -!!$ if ((prec%iprcparm(mld_coarse_mat_) /= mld_repl_mat_).or.& -!!$ & (prec%prec%iprcparm(mld_sub_solve_) /= val)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%iprcparm(mld_coarse_mat_) = mld_repl_mat_ -!!$ prec%prec%iprcparm(mld_sub_solve_) = val -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ end if -!!$ case(mld_sludist_) -!!$ if ((prec%iprcparm(mld_coarse_mat_) /= mld_distr_mat_).or.& -!!$ & (prec%prec%iprcparm(mld_sub_solve_) /= val)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%iprcparm(mld_coarse_mat_) = mld_distr_mat_ -!!$ prec%prec%iprcparm(mld_sub_solve_) = val -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ prec%prec%iprcparm(mld_smoother_sweeps_) = 1 -!!$ end if -!!$ end select - end subroutine check_coarse_lev - end subroutine mld_dmlprec_bld diff --git a/mlprec/impl/mld_s_onelev_impl.f90 b/mlprec/impl/mld_s_onelev_impl.f90 index 18cc7154..3a5d9d68 100644 --- a/mlprec/impl/mld_s_onelev_impl.f90 +++ b/mlprec/impl/mld_s_onelev_impl.f90 @@ -154,7 +154,7 @@ subroutine mld_s_base_onelev_free(lv,info) & call lv%sm%free(info) call lv%ac%free() - if (psb_is_ok_desc(lv%desc_ac)) & + if (lv%desc_ac%is_ok()) & & call psb_cdfree(lv%desc_ac,info) call lv%map%free(info) diff --git a/mlprec/impl/mld_smlprec_bld.f90 b/mlprec/impl/mld_smlprec_bld.f90 index d49576d6..de9e23fe 100644 --- a/mlprec/impl/mld_smlprec_bld.f90 +++ b/mlprec/impl/mld_smlprec_bld.f90 @@ -93,11 +93,13 @@ subroutine mld_smlprec_bld(a,desc_a,p,info,amold,vmold) ! Local Variables type(mld_sprec_type) :: t_prec - Integer :: err,i,k,ictxt, me,np, err_act, iszv, newsz + Integer :: err,i,k,ictxt, me,np, err_act, iszv, newsz, casize integer :: ipv(mld_ifpsz_), val integer :: int_err(5) character :: upd_ - type(mld_sml_parms) :: prm + class(mld_s_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm + type(mld_sml_parms) :: baseparms, medparms, coarseparms + type(mld_s_onelev_node), pointer :: head, tail, newnode, current integer :: debug_level, debug_unit character(len=20) :: name, ch_err @@ -145,12 +147,22 @@ subroutine mld_smlprec_bld(a,desc_a,p,info,amold,vmold) ! Check to ensure all procs have the same ! newsz = -1 + casize = p%coarse_aggr_size iszv = size(p%precv) call psb_bcast(ictxt,iszv) - if (iszv /= size(p%precv)) then - info=psb_err_internal_error_ - call psb_errpush(info,name,a_err='Inconsistent size of precv') - goto 9999 + call psb_bcast(ictxt,casize) + if (casize > 0) then + if (casize /= p%coarse_aggr_size) then + info=psb_err_internal_error_ + call psb_errpush(info,name,a_err='Inconsistent coarse_aggr_size') + goto 9999 + end if + else + if (iszv /= size(p%precv)) then + info=psb_err_internal_error_ + call psb_errpush(info,name,a_err='Inconsistent size of precv') + goto 9999 + end if end if if (iszv <= 1) then @@ -162,7 +174,161 @@ subroutine mld_smlprec_bld(a,desc_a,p,info,amold,vmold) endif - if (iszv > 1) then + + if (casize>0) then + ! + ! New strategy to build according to coarse size. + ! + coarseparms = p%precv(iszv)%parms + baseparms = p%precv(1)%parms + medparms = p%precv(2)%parms + + allocate(coarse_sm, source=p%precv(iszv)%sm,stat=info) + if (info == psb_success_) & + & allocate(med_sm, source=p%precv(2)%sm,stat=info) + if (info == psb_success_) & + & allocate(base_sm, source=p%precv(1)%sm,stat=info) + if (info /= psb_success_) then + write(0,*) 'Error in saving smoothers',info + call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') + goto 9999 + end if + ! + ! Replicated matrix should only ever happen at coarse level. + ! + call mld_check_def(baseparms%coarse_mat,'Coarse matrix',& + & mld_distr_mat_,is_distr_ml_coarse_mat) + ! + ! Now build a doubly linked list + ! + allocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='List start'); goto 9999 + end if + head => newnode + tail => newnode + newnode%item%base_a => a + newnode%item%base_desc => desc_a + newnode%item%parms = baseparms + newsz = 1 + current => head + list_build_loop: do + allocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='List start'); goto 9999 + end if + current%next => newnode + newnode%prev => current + newsz = newsz + 1 + newnode%item%parms = medparms + newnode%item%parms%aggr_thresh = current%item%parms%aggr_thresh/2 + call mld_coarse_bld(current%item%base_a, current%item%base_desc, & + & newnode%item,info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='build next level'); goto 9999 + end if + + if (newsz>2) then + if (all(current%item%map%naggr == newnode%item%map%naggr)) then + ! + ! We are not gaining anything + ! + newsz = newsz-1 + current%next => null() + call newnode%item%free(info) + if (info == psb_success_) deallocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Deallocate at list end'); goto 9999 + end if + end if + end if + + current => current%next + tail => current + if (sum(newnode%item%map%naggr) <= casize) then + ! + ! Target reached; but we may need to rebuild. + ! + exit list_build_loop + end if + end do list_build_loop + ! + ! At this point, we are at the list tail, + ! and it needs to be rebuilt in case the parms were + ! different. + ! + ! But the threshold has to be fixed before rebuliding + coarseparms%aggr_thresh = current%item%parms%aggr_thresh + current%item%parms = coarseparms + call mld_coarse_bld(current%prev%item%base_a,& + & current%prev%item%base_desc, & + & current%item,info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='build next level'); goto 9999 + end if + + ! + ! Ok, now allocate the output vector and fix items. + ! + do i=1,iszv + if (info == psb_success_) call p%precv(i)%free(info) + end do + if (info == psb_success_) deallocate(p%precv,stat=info) + if (info == psb_success_) allocate(p%precv(newsz),stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Reallocate precv'); goto 9999 + end if + newnode => head + do i=1, newsz + current => newnode + if (info == psb_success_) & + & call mld_move_alloc(current%item,p%precv(i),info) + if (info == psb_success_) then + if (i ==1) then + allocate(p%precv(i)%sm,source=base_sm,stat=info) + else if (i < newsz) then + allocate(p%precv(i)%sm,source=med_sm,stat=info) + else + allocate(p%precv(i)%sm,source=coarse_sm,stat=info) + end if + end if + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='list cpy'); goto 9999 + end if + if (i == 1) then + p%precv(i)%base_a => a + p%precv(i)%base_desc => desc_a + else + p%precv(i)%base_a => p%precv(i)%ac + p%precv(i)%base_desc => p%precv(i)%desc_ac + p%precv(i)%map%p_desc_X => p%precv(i-1)%base_desc + p%precv(i)%map%p_desc_Y => p%precv(i)%base_desc + end if + + newnode => current%next + deallocate(current) + end do + call base_sm%free(info) + if (info == psb_success_) call med_sm%free(info) + if (info == psb_success_) call coarse_sm%free(info) + if (info == psb_success_) deallocate(coarse_sm,med_sm,base_sm,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='final cleanup'); goto 9999 + end if + iszv = newsz + + else + ! + ! Default, oldstyle + ! ! ! Build the matrix and the transfer operators corresponding @@ -179,11 +345,6 @@ subroutine mld_smlprec_bld(a,desc_a,p,info,amold,vmold) p%precv(1)%base_a => a p%precv(1)%base_desc => desc_a - if (info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') - goto 9999 - end if - do i=2, iszv ! @@ -201,11 +362,6 @@ subroutine mld_smlprec_bld(a,desc_a,p,info,amold,vmold) ! call mld_check_def(p%precv(i)%parms%coarse_mat,'Coarse matrix',& & mld_distr_mat_,is_distr_ml_coarse_mat) - - else if (i == iszv) then - -!!$ call check_coarse_lev(p%precv(i)) - end if if (debug_level >= psb_debug_outer_) & @@ -277,9 +433,7 @@ subroutine mld_smlprec_bld(a,desc_a,p,info,amold,vmold) p%precv(i)%map%p_desc_Y => p%precv(i)%base_desc end do - i = iszv - call check_coarse_lev(p%precv(i)) if (info == psb_success_) call mld_coarse_bld(p%precv(i-1)%base_a,& & p%precv(i-1)%base_desc, p%precv(i),info) if (info /= psb_success_) then @@ -289,6 +443,12 @@ subroutine mld_smlprec_bld(a,desc_a,p,info,amold,vmold) end if end if + ! + ! The coarse space hierarchy has been build. + ! + ! Now do the preconditioner build. + ! + do i=1, iszv ! ! build the base preconditioner at level i @@ -316,10 +476,6 @@ subroutine mld_smlprec_bld(a,desc_a,p,info,amold,vmold) goto 9999 end if - - ! - ! Test version for beginning of OO stuff. - ! call p%precv(i)%sm%build(p%precv(i)%base_a,p%precv(i)%base_desc,& & 'F',info,amold=amold,vmold=vmold) @@ -350,69 +506,4 @@ subroutine mld_smlprec_bld(a,desc_a,p,info,amold,vmold) end if return -contains - - subroutine check_coarse_lev(prec) - type(mld_s_onelev_type) :: prec - - ! - ! At the coarsest level, check mld_coarse_solve_ - ! -!!$ val = prec%parms%coarse_solve -!!$ select case (val) -!!$ case(mld_jac_) -!!$ -!!$ if (prec%prec%iprcparm(mld_sub_solve_) /= mld_diag_scale_) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%prec%iprcparm(mld_sub_solve_) = mld_diag_scale_ -!!$ end if -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_jac_ -!!$ -!!$ case(mld_bjac_) -!!$ -!!$ if ((prec%prec%iprcparm(mld_sub_solve_) == mld_diag_scale_).or.& -!!$ & ( prec%prec%iprcparm(mld_smoother_type_) /= mld_bjac_)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$! !$#if defined(HAVE_UMF_) -!!$! !$ prec%prec%iprcparm(mld_sub_solve_) = mld_umf_ -!!$! !$#elif defined(HAVE_SLU_) -!!$! !$ prec%prec%iprcparm(mld_sub_solve_) = mld_slu_ -!!$! !$#else -!!$ prec%prec%iprcparm(mld_sub_solve_) = mld_ilu_n_ -!!$! !$#endif -!!$ end if -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ -!!$ case(mld_umf_, mld_slu_) -!!$ if ((prec%iprcparm(mld_coarse_mat_) /= mld_repl_mat_).or.& -!!$ & (prec%prec%iprcparm(mld_sub_solve_) /= val)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%iprcparm(mld_coarse_mat_) = mld_repl_mat_ -!!$ prec%prec%iprcparm(mld_sub_solve_) = val -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ end if -!!$ case(mld_sludist_) -!!$ if ((prec%iprcparm(mld_coarse_mat_) /= mld_distr_mat_).or.& -!!$ & (prec%prec%iprcparm(mld_sub_solve_) /= val)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%iprcparm(mld_coarse_mat_) = mld_distr_mat_ -!!$ prec%prec%iprcparm(mld_sub_solve_) = val -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ prec%prec%iprcparm(mld_smoother_sweeps_) = 1 -!!$ end if -!!$ end select - end subroutine check_coarse_lev - end subroutine mld_smlprec_bld diff --git a/mlprec/impl/mld_sprecinit.F90 b/mlprec/impl/mld_sprecinit.F90 index 12d10d9e..d8b66f83 100644 --- a/mlprec/impl/mld_sprecinit.F90 +++ b/mlprec/impl/mld_sprecinit.F90 @@ -125,6 +125,7 @@ subroutine mld_sprecinit(p,ptype,info,nlev) ! Do we want to do something? endif endif + p%coarse_aggr_size = -1 select case(psb_toupper(ptype(1:len_trim(ptype)))) case ('NOPREC','NONE') diff --git a/mlprec/impl/mld_sprecset.F90 b/mlprec/impl/mld_sprecset.F90 index 0c8f3fa6..895e5b89 100644 --- a/mlprec/impl/mld_sprecset.F90 +++ b/mlprec/impl/mld_sprecset.F90 @@ -129,6 +129,11 @@ subroutine mld_sprecseti(p,what,val,info,ilev) return endif + if (what == mld_coarse_aggr_size_) then + p%coarse_aggr_size = max(val,-1) + return + end if + ! ! Set preconditioner parameters at level ilev. ! diff --git a/mlprec/impl/mld_z_onelev_impl.f90 b/mlprec/impl/mld_z_onelev_impl.f90 index 1bdc04de..68f20bab 100644 --- a/mlprec/impl/mld_z_onelev_impl.f90 +++ b/mlprec/impl/mld_z_onelev_impl.f90 @@ -154,7 +154,7 @@ subroutine mld_z_base_onelev_free(lv,info) & call lv%sm%free(info) call lv%ac%free() - if (psb_is_ok_desc(lv%desc_ac)) & + if (lv%desc_ac%is_ok()) & & call psb_cdfree(lv%desc_ac,info) call lv%map%free(info) diff --git a/mlprec/impl/mld_zmlprec_bld.f90 b/mlprec/impl/mld_zmlprec_bld.f90 index 8f1a949a..24341f9b 100644 --- a/mlprec/impl/mld_zmlprec_bld.f90 +++ b/mlprec/impl/mld_zmlprec_bld.f90 @@ -93,11 +93,13 @@ subroutine mld_zmlprec_bld(a,desc_a,p,info,amold,vmold) ! Local Variables type(mld_zprec_type) :: t_prec - Integer :: err,i,k,ictxt, me,np, err_act, iszv, newsz + Integer :: err,i,k,ictxt, me,np, err_act, iszv, newsz, casize integer :: ipv(mld_ifpsz_), val integer :: int_err(5) character :: upd_ - type(mld_dml_parms) :: prm + class(mld_z_base_smoother_type), allocatable :: coarse_sm, base_sm, med_sm + type(mld_dml_parms) :: baseparms, medparms, coarseparms + type(mld_z_onelev_node), pointer :: head, tail, newnode, current integer :: debug_level, debug_unit character(len=20) :: name, ch_err @@ -145,12 +147,22 @@ subroutine mld_zmlprec_bld(a,desc_a,p,info,amold,vmold) ! Check to ensure all procs have the same ! newsz = -1 + casize = p%coarse_aggr_size iszv = size(p%precv) call psb_bcast(ictxt,iszv) - if (iszv /= size(p%precv)) then - info=psb_err_internal_error_ - call psb_errpush(info,name,a_err='Inconsistent size of precv') - goto 9999 + call psb_bcast(ictxt,casize) + if (casize > 0) then + if (casize /= p%coarse_aggr_size) then + info=psb_err_internal_error_ + call psb_errpush(info,name,a_err='Inconsistent coarse_aggr_size') + goto 9999 + end if + else + if (iszv /= size(p%precv)) then + info=psb_err_internal_error_ + call psb_errpush(info,name,a_err='Inconsistent size of precv') + goto 9999 + end if end if if (iszv <= 1) then @@ -162,7 +174,161 @@ subroutine mld_zmlprec_bld(a,desc_a,p,info,amold,vmold) endif - if (iszv > 1) then + + if (casize>0) then + ! + ! New strategy to build according to coarse size. + ! + coarseparms = p%precv(iszv)%parms + baseparms = p%precv(1)%parms + medparms = p%precv(2)%parms + + allocate(coarse_sm, source=p%precv(iszv)%sm,stat=info) + if (info == psb_success_) & + & allocate(med_sm, source=p%precv(2)%sm,stat=info) + if (info == psb_success_) & + & allocate(base_sm, source=p%precv(1)%sm,stat=info) + if (info /= psb_success_) then + write(0,*) 'Error in saving smoothers',info + call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') + goto 9999 + end if + ! + ! Replicated matrix should only ever happen at coarse level. + ! + call mld_check_def(baseparms%coarse_mat,'Coarse matrix',& + & mld_distr_mat_,is_distr_ml_coarse_mat) + ! + ! Now build a doubly linked list + ! + allocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='List start'); goto 9999 + end if + head => newnode + tail => newnode + newnode%item%base_a => a + newnode%item%base_desc => desc_a + newnode%item%parms = baseparms + newsz = 1 + current => head + list_build_loop: do + allocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='List start'); goto 9999 + end if + current%next => newnode + newnode%prev => current + newsz = newsz + 1 + newnode%item%parms = medparms + newnode%item%parms%aggr_thresh = current%item%parms%aggr_thresh/2 + call mld_coarse_bld(current%item%base_a, current%item%base_desc, & + & newnode%item,info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='build next level'); goto 9999 + end if + + if (newsz>2) then + if (all(current%item%map%naggr == newnode%item%map%naggr)) then + ! + ! We are not gaining anything + ! + newsz = newsz-1 + current%next => null() + call newnode%item%free(info) + if (info == psb_success_) deallocate(newnode,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Deallocate at list end'); goto 9999 + end if + end if + end if + + current => current%next + tail => current + if (sum(newnode%item%map%naggr) <= casize) then + ! + ! Target reached; but we may need to rebuild. + ! + exit list_build_loop + end if + end do list_build_loop + ! + ! At this point, we are at the list tail, + ! and it needs to be rebuilt in case the parms were + ! different. + ! + ! But the threshold has to be fixed before rebuliding + coarseparms%aggr_thresh = current%item%parms%aggr_thresh + current%item%parms = coarseparms + call mld_coarse_bld(current%prev%item%base_a,& + & current%prev%item%base_desc, & + & current%item,info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='build next level'); goto 9999 + end if + + ! + ! Ok, now allocate the output vector and fix items. + ! + do i=1,iszv + if (info == psb_success_) call p%precv(i)%free(info) + end do + if (info == psb_success_) deallocate(p%precv,stat=info) + if (info == psb_success_) allocate(p%precv(newsz),stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Reallocate precv'); goto 9999 + end if + newnode => head + do i=1, newsz + current => newnode + if (info == psb_success_) & + & call mld_move_alloc(current%item,p%precv(i),info) + if (info == psb_success_) then + if (i ==1) then + allocate(p%precv(i)%sm,source=base_sm,stat=info) + else if (i < newsz) then + allocate(p%precv(i)%sm,source=med_sm,stat=info) + else + allocate(p%precv(i)%sm,source=coarse_sm,stat=info) + end if + end if + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='list cpy'); goto 9999 + end if + if (i == 1) then + p%precv(i)%base_a => a + p%precv(i)%base_desc => desc_a + else + p%precv(i)%base_a => p%precv(i)%ac + p%precv(i)%base_desc => p%precv(i)%desc_ac + p%precv(i)%map%p_desc_X => p%precv(i-1)%base_desc + p%precv(i)%map%p_desc_Y => p%precv(i)%base_desc + end if + + newnode => current%next + deallocate(current) + end do + call base_sm%free(info) + if (info == psb_success_) call med_sm%free(info) + if (info == psb_success_) call coarse_sm%free(info) + if (info == psb_success_) deallocate(coarse_sm,med_sm,base_sm,stat=info) + if (info /= psb_success_) then + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='final cleanup'); goto 9999 + end if + iszv = newsz + + else + ! + ! Default, oldstyle + ! ! ! Build the matrix and the transfer operators corresponding @@ -179,11 +345,6 @@ subroutine mld_zmlprec_bld(a,desc_a,p,info,amold,vmold) p%precv(1)%base_a => a p%precv(1)%base_desc => desc_a - if (info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,a_err='Base level precbuild.') - goto 9999 - end if - do i=2, iszv ! @@ -201,11 +362,6 @@ subroutine mld_zmlprec_bld(a,desc_a,p,info,amold,vmold) ! call mld_check_def(p%precv(i)%parms%coarse_mat,'Coarse matrix',& & mld_distr_mat_,is_distr_ml_coarse_mat) - - else if (i == iszv) then - -!!$ call check_coarse_lev(p%precv(i)) - end if if (debug_level >= psb_debug_outer_) & @@ -277,9 +433,7 @@ subroutine mld_zmlprec_bld(a,desc_a,p,info,amold,vmold) p%precv(i)%map%p_desc_Y => p%precv(i)%base_desc end do - i = iszv - call check_coarse_lev(p%precv(i)) if (info == psb_success_) call mld_coarse_bld(p%precv(i-1)%base_a,& & p%precv(i-1)%base_desc, p%precv(i),info) if (info /= psb_success_) then @@ -289,6 +443,12 @@ subroutine mld_zmlprec_bld(a,desc_a,p,info,amold,vmold) end if end if + ! + ! The coarse space hierarchy has been build. + ! + ! Now do the preconditioner build. + ! + do i=1, iszv ! ! build the base preconditioner at level i @@ -316,10 +476,6 @@ subroutine mld_zmlprec_bld(a,desc_a,p,info,amold,vmold) goto 9999 end if - - ! - ! Test version for beginning of OO stuff. - ! call p%precv(i)%sm%build(p%precv(i)%base_a,p%precv(i)%base_desc,& & 'F',info,amold=amold,vmold=vmold) @@ -350,69 +506,4 @@ subroutine mld_zmlprec_bld(a,desc_a,p,info,amold,vmold) end if return -contains - - subroutine check_coarse_lev(prec) - type(mld_z_onelev_type) :: prec - - ! - ! At the coarsest level, check mld_coarse_solve_ - ! -!!$ val = prec%parms%coarse_solve -!!$ select case (val) -!!$ case(mld_jac_) -!!$ -!!$ if (prec%prec%iprcparm(mld_sub_solve_) /= mld_diag_scale_) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%prec%iprcparm(mld_sub_solve_) = mld_diag_scale_ -!!$ end if -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_jac_ -!!$ -!!$ case(mld_bjac_) -!!$ -!!$ if ((prec%prec%iprcparm(mld_sub_solve_) == mld_diag_scale_).or.& -!!$ & ( prec%prec%iprcparm(mld_smoother_type_) /= mld_bjac_)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$! !$#if defined(HAVE_UMF_) -!!$! !$ prec%prec%iprcparm(mld_sub_solve_) = mld_umf_ -!!$! !$#elif defined(HAVE_SLU_) -!!$! !$ prec%prec%iprcparm(mld_sub_solve_) = mld_slu_ -!!$! !$#else -!!$ prec%prec%iprcparm(mld_sub_solve_) = mld_ilu_n_ -!!$! !$#endif -!!$ end if -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ -!!$ case(mld_umf_, mld_slu_) -!!$ if ((prec%iprcparm(mld_coarse_mat_) /= mld_repl_mat_).or.& -!!$ & (prec%prec%iprcparm(mld_sub_solve_) /= val)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%iprcparm(mld_coarse_mat_) = mld_repl_mat_ -!!$ prec%prec%iprcparm(mld_sub_solve_) = val -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ end if -!!$ case(mld_sludist_) -!!$ if ((prec%iprcparm(mld_coarse_mat_) /= mld_distr_mat_).or.& -!!$ & (prec%prec%iprcparm(mld_sub_solve_) /= val)) then -!!$ if (me == 0) write(debug_unit,*)& -!!$ & 'Warning: inconsistent coarse level specification.' -!!$ if (me == 0) write(debug_unit,*)& -!!$ & ' Resetting according to the value specified for mld_coarse_solve_.' -!!$ prec%iprcparm(mld_coarse_mat_) = mld_distr_mat_ -!!$ prec%prec%iprcparm(mld_sub_solve_) = val -!!$ prec%prec%iprcparm(mld_smoother_type_) = mld_bjac_ -!!$ prec%prec%iprcparm(mld_smoother_sweeps_) = 1 -!!$ end if -!!$ end select - end subroutine check_coarse_lev - end subroutine mld_zmlprec_bld diff --git a/mlprec/impl/mld_zprecinit.F90 b/mlprec/impl/mld_zprecinit.F90 index 1bf0f172..b5cd4a0e 100644 --- a/mlprec/impl/mld_zprecinit.F90 +++ b/mlprec/impl/mld_zprecinit.F90 @@ -125,6 +125,7 @@ subroutine mld_zprecinit(p,ptype,info,nlev) ! Do we want to do something? endif endif + p%coarse_aggr_size = -1 select case(psb_toupper(ptype(1:len_trim(ptype)))) case ('NOPREC','NONE') diff --git a/mlprec/impl/mld_zprecset.F90 b/mlprec/impl/mld_zprecset.F90 index c49244ac..cd1bfb1d 100644 --- a/mlprec/impl/mld_zprecset.F90 +++ b/mlprec/impl/mld_zprecset.F90 @@ -129,6 +129,11 @@ subroutine mld_zprecseti(p,what,val,info,ilev) return endif + if (what == mld_coarse_aggr_size_) then + p%coarse_aggr_size = max(val,-1) + return + end if + ! ! Set preconditioner parameters at level ilev. ! diff --git a/mlprec/mld_base_prec_type.F90 b/mlprec/mld_base_prec_type.F90 index 83de58d5..f5d2caa3 100644 --- a/mlprec/mld_base_prec_type.F90 +++ b/mlprec/mld_base_prec_type.F90 @@ -104,6 +104,7 @@ module mld_base_prec_type procedure, pass(pm) :: descr => ml_parms_descr procedure, pass(pm) :: mldescr => ml_parms_mldescr procedure, pass(pm) :: coarsedescr => ml_parms_coarsedescr + procedure, pass(pm) :: printout => ml_parms_printout end type mld_ml_parms @@ -111,12 +112,14 @@ module mld_base_prec_type real(psb_spk_) :: aggr_omega_val, aggr_thresh contains procedure, pass(pm) :: descr => s_ml_parms_descr + procedure, pass(pm) :: printout => s_ml_parms_printout end type mld_sml_parms type, extends(mld_ml_parms) :: mld_dml_parms real(psb_dpk_) :: aggr_omega_val, aggr_thresh contains procedure, pass(pm) :: descr => d_ml_parms_descr + procedure, pass(pm) :: printout => d_ml_parms_printout end type mld_dml_parms @@ -437,6 +440,40 @@ contains end subroutine mld_stringval + + subroutine ml_parms_printout(pm,iout) + implicit none + class(mld_ml_parms), intent(in) :: pm + integer, intent(in) :: iout + + 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,*) ' : ',pm%aggr_omega_alg,pm%aggr_eig,pm%aggr_filter + write(iout,*) 'COARSE: ',pm%coarse_mat,pm%coarse_solve + end subroutine ml_parms_printout + + + subroutine s_ml_parms_printout(pm,iout) + implicit none + class(mld_sml_parms), intent(in) :: pm + integer, intent(in) :: iout + + call pm%mld_ml_parms%printout(iout) + write(iout,*) 'REAL : ',pm%aggr_omega_val,pm%aggr_thresh + end subroutine s_ml_parms_printout + + + subroutine d_ml_parms_printout(pm,iout) + implicit none + class(mld_dml_parms), intent(in) :: pm + integer, intent(in) :: iout + + call pm%mld_ml_parms%printout(iout) + write(iout,*) 'REAL : ',pm%aggr_omega_val,pm%aggr_thresh + end subroutine d_ml_parms_printout + + ! ! Routines printing out a description of the preconditioner ! diff --git a/mlprec/mld_c_inner_mod.f90 b/mlprec/mld_c_inner_mod.f90 index 8cf313d0..5ff82c70 100644 --- a/mlprec/mld_c_inner_mod.f90 +++ b/mlprec/mld_c_inner_mod.f90 @@ -46,8 +46,6 @@ ! module mld_c_inner_mod use mld_c_prec_type - use mld_c_move_alloc_mod - interface mld_mlprec_bld subroutine mld_cmlprec_bld(a,desc_a,prec,info, amold, vmold) diff --git a/mlprec/mld_c_move_alloc_mod.f90 b/mlprec/mld_c_move_alloc_mod.f90 index 85f8120b..d936c0fd 100644 --- a/mlprec/mld_c_move_alloc_mod.f90 +++ b/mlprec/mld_c_move_alloc_mod.f90 @@ -62,6 +62,7 @@ contains integer, intent(out) :: info call b%free(info) + b%parms = a%parms call move_alloc(a%sm,b%sm) if (info == psb_success_) call psb_move_alloc(a%ac,b%ac,info) if (info == psb_success_) call psb_move_alloc(a%desc_ac,b%desc_ac,info) diff --git a/mlprec/mld_c_onelev_mod.f90 b/mlprec/mld_c_onelev_mod.f90 index ba329fb7..2f9d30a2 100644 --- a/mlprec/mld_c_onelev_mod.f90 +++ b/mlprec/mld_c_onelev_mod.f90 @@ -141,6 +141,11 @@ module mld_c_onelev_mod procedure, pass(lv) :: get_nzeros => c_base_onelev_get_nzeros end type mld_c_onelev_type + type mld_c_onelev_node + type(mld_c_onelev_type) :: item + type(mld_c_onelev_node), pointer :: prev=>null(), next=>null() + end type mld_c_onelev_node + private :: c_base_onelev_default, c_base_onelev_sizeof, & & c_base_onelev_nullify, c_base_onelev_get_nzeros @@ -234,6 +239,9 @@ module mld_c_onelev_mod end subroutine mld_c_base_onelev_dump end interface + interface mld_move_alloc + module procedure mld_c_onelev_move_alloc + end interface contains ! @@ -312,4 +320,22 @@ contains end subroutine c_base_onelev_default + + subroutine mld_c_onelev_move_alloc(a, b,info) + use psb_base_mod + implicit none + type(mld_c_onelev_type), intent(inout) :: a, b + integer, intent(out) :: info + + call b%free(info) + b%parms = a%parms + call move_alloc(a%sm,b%sm) + if (info == psb_success_) call psb_move_alloc(a%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(a%desc_ac,b%desc_ac,info) + if (info == psb_success_) call psb_move_alloc(a%map,b%map,info) + b%base_a => a%base_a + b%base_desc => a%base_desc + + end subroutine mld_c_onelev_move_alloc + end module mld_c_onelev_mod diff --git a/mlprec/mld_c_prec_mod.f90 b/mlprec/mld_c_prec_mod.f90 index 142a6d93..52f0ad8c 100644 --- a/mlprec/mld_c_prec_mod.f90 +++ b/mlprec/mld_c_prec_mod.f90 @@ -46,7 +46,6 @@ module mld_c_prec_mod use mld_c_prec_type - use mld_c_move_alloc_mod interface mld_precinit subroutine mld_cprecinit(p,ptype,info,nlev) diff --git a/mlprec/mld_c_prec_type.f90 b/mlprec/mld_c_prec_type.f90 index a39eef0b..690def55 100644 --- a/mlprec/mld_c_prec_type.f90 +++ b/mlprec/mld_c_prec_type.f90 @@ -81,6 +81,7 @@ module mld_c_prec_type type, extends(psb_cprec_type) :: mld_cprec_type integer :: ictxt + integer(psb_ipk_) :: coarse_aggr_size real(psb_spk_) :: op_complexity=szero type(mld_c_onelev_type), allocatable :: precv(:) contains @@ -159,6 +160,10 @@ module mld_c_prec_type end subroutine mld_cprecaply1 end interface + interface mld_move_alloc + module procedure mld_cprec_move_alloc + end interface + contains ! ! Function returning the size of the mld_prec_type data structure @@ -565,7 +570,7 @@ contains if (present(istart)) then il1 = max(1,istart) else - il1 = 2 + il1 = min(2,iln) end if if (present(iend)) then iln = min(iln, iend) @@ -577,5 +582,32 @@ contains end do end subroutine mld_c_dump + + subroutine mld_cprec_move_alloc(a, b,info) + use psb_base_mod + implicit none + type(mld_cprec_type), intent(inout) :: a + type(mld_cprec_type), intent(inout), target :: b + integer, intent(out) :: info + integer :: i + + if (allocated(b%precv)) then + ! This might not be required if FINAL procedures are available. + call mld_precfree(b,info) + if (info /= psb_success_) then + ! ????? + !!$ return + endif + end if + + call move_alloc(a%precv,b%precv) + ! Fix the pointers except on level 1. + do i=2, size(b%precv) + b%precv(i)%base_a => b%precv(i)%ac + b%precv(i)%base_desc => b%precv(i)%desc_ac + b%precv(i)%map%p_desc_X => b%precv(i-1)%base_desc + b%precv(i)%map%p_desc_Y => b%precv(i)%base_desc + end do + end subroutine mld_cprec_move_alloc end module mld_c_prec_type diff --git a/mlprec/mld_d_inner_mod.f90 b/mlprec/mld_d_inner_mod.f90 index b07aec6c..5b8a2979 100644 --- a/mlprec/mld_d_inner_mod.f90 +++ b/mlprec/mld_d_inner_mod.f90 @@ -46,8 +46,6 @@ ! module mld_d_inner_mod use mld_d_prec_type - use mld_d_move_alloc_mod - interface mld_mlprec_bld subroutine mld_dmlprec_bld(a,desc_a,prec,info, amold, vmold) diff --git a/mlprec/mld_d_move_alloc_mod.f90 b/mlprec/mld_d_move_alloc_mod.f90 index d7bc7f1e..27f26fee 100644 --- a/mlprec/mld_d_move_alloc_mod.f90 +++ b/mlprec/mld_d_move_alloc_mod.f90 @@ -62,6 +62,7 @@ contains integer, intent(out) :: info call b%free(info) + b%parms = a%parms call move_alloc(a%sm,b%sm) if (info == psb_success_) call psb_move_alloc(a%ac,b%ac,info) if (info == psb_success_) call psb_move_alloc(a%desc_ac,b%desc_ac,info) diff --git a/mlprec/mld_d_onelev_mod.f90 b/mlprec/mld_d_onelev_mod.f90 index a17a0333..a5689a69 100644 --- a/mlprec/mld_d_onelev_mod.f90 +++ b/mlprec/mld_d_onelev_mod.f90 @@ -141,6 +141,11 @@ module mld_d_onelev_mod procedure, pass(lv) :: get_nzeros => d_base_onelev_get_nzeros end type mld_d_onelev_type + type mld_d_onelev_node + type(mld_d_onelev_type) :: item + type(mld_d_onelev_node), pointer :: prev=>null(), next=>null() + end type mld_d_onelev_node + private :: d_base_onelev_default, d_base_onelev_sizeof, & & d_base_onelev_nullify, d_base_onelev_get_nzeros @@ -234,6 +239,9 @@ module mld_d_onelev_mod end subroutine mld_d_base_onelev_dump end interface + interface mld_move_alloc + module procedure mld_d_onelev_move_alloc + end interface contains ! @@ -312,4 +320,22 @@ contains end subroutine d_base_onelev_default + + subroutine mld_d_onelev_move_alloc(a, b,info) + use psb_base_mod + implicit none + type(mld_d_onelev_type), intent(inout) :: a, b + integer, intent(out) :: info + + call b%free(info) + b%parms = a%parms + call move_alloc(a%sm,b%sm) + if (info == psb_success_) call psb_move_alloc(a%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(a%desc_ac,b%desc_ac,info) + if (info == psb_success_) call psb_move_alloc(a%map,b%map,info) + b%base_a => a%base_a + b%base_desc => a%base_desc + + end subroutine mld_d_onelev_move_alloc + end module mld_d_onelev_mod diff --git a/mlprec/mld_d_prec_mod.f90 b/mlprec/mld_d_prec_mod.f90 index a4006aca..34231e12 100644 --- a/mlprec/mld_d_prec_mod.f90 +++ b/mlprec/mld_d_prec_mod.f90 @@ -46,7 +46,6 @@ module mld_d_prec_mod use mld_d_prec_type - use mld_d_move_alloc_mod interface mld_precinit subroutine mld_dprecinit(p,ptype,info,nlev) diff --git a/mlprec/mld_d_prec_type.f90 b/mlprec/mld_d_prec_type.f90 index 86ba03f2..296751c5 100644 --- a/mlprec/mld_d_prec_type.f90 +++ b/mlprec/mld_d_prec_type.f90 @@ -81,7 +81,7 @@ module mld_d_prec_type type, extends(psb_dprec_type) :: mld_dprec_type integer :: ictxt - integer :: coarse_aggr_size + integer(psb_ipk_) :: coarse_aggr_size real(psb_dpk_) :: op_complexity=dzero type(mld_d_onelev_type), allocatable :: precv(:) contains @@ -160,6 +160,10 @@ module mld_d_prec_type end subroutine mld_dprecaply1 end interface + interface mld_move_alloc + module procedure mld_dprec_move_alloc + end interface + contains ! ! Function returning the size of the mld_prec_type data structure @@ -566,7 +570,7 @@ contains if (present(istart)) then il1 = max(1,istart) else - il1 = 2 + il1 = min(2,iln) end if if (present(iend)) then iln = min(iln, iend) @@ -578,5 +582,32 @@ contains end do end subroutine mld_d_dump + + subroutine mld_dprec_move_alloc(a, b,info) + use psb_base_mod + implicit none + type(mld_dprec_type), intent(inout) :: a + type(mld_dprec_type), intent(inout), target :: b + integer, intent(out) :: info + integer :: i + + if (allocated(b%precv)) then + ! This might not be required if FINAL procedures are available. + call mld_precfree(b,info) + if (info /= psb_success_) then + ! ????? + !!$ return + endif + end if + + call move_alloc(a%precv,b%precv) + ! Fix the pointers except on level 1. + do i=2, size(b%precv) + b%precv(i)%base_a => b%precv(i)%ac + b%precv(i)%base_desc => b%precv(i)%desc_ac + b%precv(i)%map%p_desc_X => b%precv(i-1)%base_desc + b%precv(i)%map%p_desc_Y => b%precv(i)%base_desc + end do + end subroutine mld_dprec_move_alloc end module mld_d_prec_type diff --git a/mlprec/mld_s_inner_mod.f90 b/mlprec/mld_s_inner_mod.f90 index db48e913..b09fa587 100644 --- a/mlprec/mld_s_inner_mod.f90 +++ b/mlprec/mld_s_inner_mod.f90 @@ -46,8 +46,6 @@ ! module mld_s_inner_mod use mld_s_prec_type - use mld_s_move_alloc_mod - interface mld_mlprec_bld subroutine mld_smlprec_bld(a,desc_a,prec,info, amold, vmold) diff --git a/mlprec/mld_s_move_alloc_mod.f90 b/mlprec/mld_s_move_alloc_mod.f90 index d23fefe9..07676612 100644 --- a/mlprec/mld_s_move_alloc_mod.f90 +++ b/mlprec/mld_s_move_alloc_mod.f90 @@ -62,6 +62,7 @@ contains integer, intent(out) :: info call b%free(info) + b%parms = a%parms call move_alloc(a%sm,b%sm) if (info == psb_success_) call psb_move_alloc(a%ac,b%ac,info) if (info == psb_success_) call psb_move_alloc(a%desc_ac,b%desc_ac,info) diff --git a/mlprec/mld_s_onelev_mod.f90 b/mlprec/mld_s_onelev_mod.f90 index 17f69512..890467c2 100644 --- a/mlprec/mld_s_onelev_mod.f90 +++ b/mlprec/mld_s_onelev_mod.f90 @@ -141,6 +141,11 @@ module mld_s_onelev_mod procedure, pass(lv) :: get_nzeros => s_base_onelev_get_nzeros end type mld_s_onelev_type + type mld_s_onelev_node + type(mld_s_onelev_type) :: item + type(mld_s_onelev_node), pointer :: prev=>null(), next=>null() + end type mld_s_onelev_node + private :: s_base_onelev_default, s_base_onelev_sizeof, & & s_base_onelev_nullify, s_base_onelev_get_nzeros @@ -234,6 +239,9 @@ module mld_s_onelev_mod end subroutine mld_s_base_onelev_dump end interface + interface mld_move_alloc + module procedure mld_s_onelev_move_alloc + end interface contains ! @@ -312,4 +320,22 @@ contains end subroutine s_base_onelev_default + + subroutine mld_s_onelev_move_alloc(a, b,info) + use psb_base_mod + implicit none + type(mld_s_onelev_type), intent(inout) :: a, b + integer, intent(out) :: info + + call b%free(info) + b%parms = a%parms + call move_alloc(a%sm,b%sm) + if (info == psb_success_) call psb_move_alloc(a%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(a%desc_ac,b%desc_ac,info) + if (info == psb_success_) call psb_move_alloc(a%map,b%map,info) + b%base_a => a%base_a + b%base_desc => a%base_desc + + end subroutine mld_s_onelev_move_alloc + end module mld_s_onelev_mod diff --git a/mlprec/mld_s_prec_mod.f90 b/mlprec/mld_s_prec_mod.f90 index f67c0df6..22981b22 100644 --- a/mlprec/mld_s_prec_mod.f90 +++ b/mlprec/mld_s_prec_mod.f90 @@ -46,7 +46,6 @@ module mld_s_prec_mod use mld_s_prec_type - use mld_s_move_alloc_mod interface mld_precinit subroutine mld_sprecinit(p,ptype,info,nlev) diff --git a/mlprec/mld_s_prec_type.f90 b/mlprec/mld_s_prec_type.f90 index 5972a517..da4cc788 100644 --- a/mlprec/mld_s_prec_type.f90 +++ b/mlprec/mld_s_prec_type.f90 @@ -81,6 +81,7 @@ module mld_s_prec_type type, extends(psb_sprec_type) :: mld_sprec_type integer :: ictxt + integer(psb_ipk_) :: coarse_aggr_size real(psb_spk_) :: op_complexity=szero type(mld_s_onelev_type), allocatable :: precv(:) contains @@ -159,6 +160,10 @@ module mld_s_prec_type end subroutine mld_sprecaply1 end interface + interface mld_move_alloc + module procedure mld_sprec_move_alloc + end interface + contains ! ! Function returning the size of the mld_prec_type data structure @@ -565,7 +570,7 @@ contains if (present(istart)) then il1 = max(1,istart) else - il1 = 2 + il1 = min(2,iln) end if if (present(iend)) then iln = min(iln, iend) @@ -577,5 +582,32 @@ contains end do end subroutine mld_s_dump + + subroutine mld_sprec_move_alloc(a, b,info) + use psb_base_mod + implicit none + type(mld_sprec_type), intent(inout) :: a + type(mld_sprec_type), intent(inout), target :: b + integer, intent(out) :: info + integer :: i + + if (allocated(b%precv)) then + ! This might not be required if FINAL procedures are available. + call mld_precfree(b,info) + if (info /= psb_success_) then + ! ????? + !!$ return + endif + end if + + call move_alloc(a%precv,b%precv) + ! Fix the pointers except on level 1. + do i=2, size(b%precv) + b%precv(i)%base_a => b%precv(i)%ac + b%precv(i)%base_desc => b%precv(i)%desc_ac + b%precv(i)%map%p_desc_X => b%precv(i-1)%base_desc + b%precv(i)%map%p_desc_Y => b%precv(i)%base_desc + end do + end subroutine mld_sprec_move_alloc end module mld_s_prec_type diff --git a/mlprec/mld_z_inner_mod.f90 b/mlprec/mld_z_inner_mod.f90 index 4b05ceb8..2428a3ec 100644 --- a/mlprec/mld_z_inner_mod.f90 +++ b/mlprec/mld_z_inner_mod.f90 @@ -46,8 +46,6 @@ ! module mld_z_inner_mod use mld_z_prec_type - use mld_z_move_alloc_mod - interface mld_mlprec_bld subroutine mld_zmlprec_bld(a,desc_a,prec,info, amold, vmold) diff --git a/mlprec/mld_z_move_alloc_mod.f90 b/mlprec/mld_z_move_alloc_mod.f90 index 098d763a..59be7ece 100644 --- a/mlprec/mld_z_move_alloc_mod.f90 +++ b/mlprec/mld_z_move_alloc_mod.f90 @@ -62,6 +62,7 @@ contains integer, intent(out) :: info call b%free(info) + b%parms = a%parms call move_alloc(a%sm,b%sm) if (info == psb_success_) call psb_move_alloc(a%ac,b%ac,info) if (info == psb_success_) call psb_move_alloc(a%desc_ac,b%desc_ac,info) diff --git a/mlprec/mld_z_onelev_mod.f90 b/mlprec/mld_z_onelev_mod.f90 index 6dfed8c2..3e2e7da2 100644 --- a/mlprec/mld_z_onelev_mod.f90 +++ b/mlprec/mld_z_onelev_mod.f90 @@ -141,6 +141,11 @@ module mld_z_onelev_mod procedure, pass(lv) :: get_nzeros => z_base_onelev_get_nzeros end type mld_z_onelev_type + type mld_z_onelev_node + type(mld_z_onelev_type) :: item + type(mld_z_onelev_node), pointer :: prev=>null(), next=>null() + end type mld_z_onelev_node + private :: z_base_onelev_default, z_base_onelev_sizeof, & & z_base_onelev_nullify, z_base_onelev_get_nzeros @@ -234,6 +239,9 @@ module mld_z_onelev_mod end subroutine mld_z_base_onelev_dump end interface + interface mld_move_alloc + module procedure mld_z_onelev_move_alloc + end interface contains ! @@ -312,4 +320,22 @@ contains end subroutine z_base_onelev_default + + subroutine mld_z_onelev_move_alloc(a, b,info) + use psb_base_mod + implicit none + type(mld_z_onelev_type), intent(inout) :: a, b + integer, intent(out) :: info + + call b%free(info) + b%parms = a%parms + call move_alloc(a%sm,b%sm) + if (info == psb_success_) call psb_move_alloc(a%ac,b%ac,info) + if (info == psb_success_) call psb_move_alloc(a%desc_ac,b%desc_ac,info) + if (info == psb_success_) call psb_move_alloc(a%map,b%map,info) + b%base_a => a%base_a + b%base_desc => a%base_desc + + end subroutine mld_z_onelev_move_alloc + end module mld_z_onelev_mod diff --git a/mlprec/mld_z_prec_mod.f90 b/mlprec/mld_z_prec_mod.f90 index 0ec94a43..542aa12a 100644 --- a/mlprec/mld_z_prec_mod.f90 +++ b/mlprec/mld_z_prec_mod.f90 @@ -46,7 +46,6 @@ module mld_z_prec_mod use mld_z_prec_type - use mld_z_move_alloc_mod interface mld_precinit subroutine mld_zprecinit(p,ptype,info,nlev) diff --git a/mlprec/mld_z_prec_type.f90 b/mlprec/mld_z_prec_type.f90 index cd29b8af..f910d7ba 100644 --- a/mlprec/mld_z_prec_type.f90 +++ b/mlprec/mld_z_prec_type.f90 @@ -81,6 +81,7 @@ module mld_z_prec_type type, extends(psb_zprec_type) :: mld_zprec_type integer :: ictxt + integer(psb_ipk_) :: coarse_aggr_size real(psb_dpk_) :: op_complexity=dzero type(mld_z_onelev_type), allocatable :: precv(:) contains @@ -159,6 +160,10 @@ module mld_z_prec_type end subroutine mld_zprecaply1 end interface + interface mld_move_alloc + module procedure mld_zprec_move_alloc + end interface + contains ! ! Function returning the size of the mld_prec_type data structure @@ -565,7 +570,7 @@ contains if (present(istart)) then il1 = max(1,istart) else - il1 = 2 + il1 = min(2,iln) end if if (present(iend)) then iln = min(iln, iend) @@ -577,5 +582,32 @@ contains end do end subroutine mld_z_dump + + subroutine mld_zprec_move_alloc(a, b,info) + use psb_base_mod + implicit none + type(mld_zprec_type), intent(inout) :: a + type(mld_zprec_type), intent(inout), target :: b + integer, intent(out) :: info + integer :: i + + if (allocated(b%precv)) then + ! This might not be required if FINAL procedures are available. + call mld_precfree(b,info) + if (info /= psb_success_) then + ! ????? + !!$ return + endif + end if + + call move_alloc(a%precv,b%precv) + ! Fix the pointers except on level 1. + do i=2, size(b%precv) + b%precv(i)%base_a => b%precv(i)%ac + b%precv(i)%base_desc => b%precv(i)%desc_ac + b%precv(i)%map%p_desc_X => b%precv(i-1)%base_desc + b%precv(i)%map%p_desc_Y => b%precv(i)%base_desc + end do + end subroutine mld_zprec_move_alloc end module mld_z_prec_type diff --git a/tests/pdegen/ppde2d.f90 b/tests/pdegen/ppde2d.f90 index 69fc77e6..d3d8e95d 100644 --- a/tests/pdegen/ppde2d.f90 +++ b/tests/pdegen/ppde2d.f90 @@ -110,6 +110,7 @@ program ppde2d character(len=16) :: aggr_alg ! local or global aggregation character(len=16) :: mltype ! additive or multiplicative 2nd level prec character(len=16) :: smthpos ! side: pre, post, both smoothing + integer :: csize ! aggregation size at which to stop. character(len=16) :: cmat ! coarse mat character(len=16) :: csolve ! Coarse solver: bjac, umf, slu, sludist character(len=16) :: csbsolve ! Coarse subsolver: ILU, ILU(T), SuperLU, UMFPACK. @@ -196,6 +197,7 @@ program ppde2d call mld_precset(prec,mld_coarse_fillin_, prectype%cfill, info) call mld_precset(prec,mld_coarse_iluthrs_, prectype%cthres, info) call mld_precset(prec,mld_coarse_sweeps_, prectype%cjswp, info) + call mld_precset(prec,mld_coarse_aggr_size_, prectype%csize, info) else nlv = 1 call mld_precinit(prec,prectype%prec, info, nlev=nlv) @@ -336,6 +338,7 @@ contains call read_data(prectype%cthres,5) ! Threshold for fact. 1 ILU(T) call read_data(prectype%cjswp,5) ! Jacobi sweeps call read_data(prectype%athres,5) ! smoother aggr thresh + call read_data(prectype%csize,5) ! coarse size end if end if @@ -373,6 +376,7 @@ contains call psb_bcast(ictxt,prectype%cthres) ! Threshold for fact. 1 ILU(T) call psb_bcast(ictxt,prectype%cjswp) ! Jacobi sweeps call psb_bcast(ictxt,prectype%athres) ! smoother aggr thresh + call psb_bcast(ictxt,prectype%csize) ! coarse size end if if (iam == psb_root_) then diff --git a/tests/pdegen/ppde3d.f90 b/tests/pdegen/ppde3d.f90 index db96119c..dfeebc3a 100644 --- a/tests/pdegen/ppde3d.f90 +++ b/tests/pdegen/ppde3d.f90 @@ -111,6 +111,7 @@ program ppde3d character(len=16) :: aggr_alg ! local or global aggregation character(len=16) :: mltype ! additive or multiplicative 2nd level prec character(len=16) :: smthpos ! side: pre, post, both smoothing + integer :: csize ! aggregation size at which to stop. character(len=16) :: cmat ! coarse mat character(len=16) :: csolve ! Coarse solver: bjac, umf, slu, sludist character(len=16) :: csbsolve ! Coarse subsolver: ILU, ILU(T), SuperLU, UMFPACK. @@ -200,6 +201,7 @@ program ppde3d call mld_precset(prec,mld_coarse_fillin_, prectype%cfill, info) call mld_precset(prec,mld_coarse_iluthrs_, prectype%cthres, info) call mld_precset(prec,mld_coarse_sweeps_, prectype%cjswp, info) + call mld_precset(prec,mld_coarse_aggr_size_, prectype%csize, info) else nlv = 1 call mld_precinit(prec,prectype%prec, info, nlev=nlv) @@ -340,6 +342,7 @@ contains call read_data(prectype%cthres,5) ! Threshold for fact. 1 ILU(T) call read_data(prectype%cjswp,5) ! Jacobi sweeps call read_data(prectype%athres,5) ! smoother aggr thresh + call read_data(prectype%csize,5) ! coarse size end if end if @@ -377,6 +380,7 @@ contains call psb_bcast(ictxt,prectype%cthres) ! Threshold for fact. 1 ILU(T) call psb_bcast(ictxt,prectype%cjswp) ! Jacobi sweeps call psb_bcast(ictxt,prectype%athres) ! smoother aggr thresh + call psb_bcast(ictxt,prectype%csize) ! coarse size end if if (iam == psb_root_) then diff --git a/tests/pdegen/runs/ppde.inp b/tests/pdegen/runs/ppde.inp index 4e674318..464d1c02 100644 --- a/tests/pdegen/runs/ppde.inp +++ b/tests/pdegen/runs/ppde.inp @@ -1,6 +1,6 @@ BICGSTAB ! Iterative method: BiCGSTAB BiCG CGS RGMRES BiCGSTABL CG CSR ! Storage format CSR COO JAD -040 ! IDIM; domain size is idim**3 +100 ! IDIM; domain size is idim**3 2 ! ISTOPC 0100 ! ITMAX -1 ! ITRACE @@ -21,10 +21,11 @@ NONSMOOTHED ! Kind of aggregation: SMOOTHED, NONSMOOTHED DEC ! Type of aggregation DEC SYMDEC GLB MULT ! Type of multilevel correction: ADD MULT TWOSIDE ! Side of correction PRE POST TWOSIDE (ignored for ADD) -DIST ! Coarse level: matrix distribution DIST REPL +REPL ! Coarse level: matrix distribution DIST REPL BJAC ! Coarse level: solver JACOBI BJAC UMF SLU SLUDIST ILU ! Coarse level: subsolver DSCALE ILU UMF SLU SLUDIST 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 -0.10d0 ! Smoother Aggregation Threshold: >= 0.0 default if <0 +100 ! Coarse size limit to determine levels. If <0, then use NL diff --git a/tests/pdegen/spde2d.f90 b/tests/pdegen/spde2d.f90 index aca14652..b0e5147d 100644 --- a/tests/pdegen/spde2d.f90 +++ b/tests/pdegen/spde2d.f90 @@ -110,6 +110,7 @@ program spde2d character(len=16) :: aggr_alg ! local or global aggregation character(len=16) :: mltype ! additive or multiplicative 2nd level prec character(len=16) :: smthpos ! side: pre, post, both smoothing + integer :: csize ! aggregation size at which to stop. character(len=16) :: cmat ! coarse mat character(len=16) :: csolve ! Coarse solver: bjac, umf, slu, sludist character(len=16) :: csbsolve ! Coarse subsolver: ILU, ILU(T), SuperLU, UMFPACK. @@ -196,6 +197,7 @@ program spde2d call mld_precset(prec,mld_coarse_fillin_, prectype%cfill, info) call mld_precset(prec,mld_coarse_iluthrs_, prectype%cthres, info) call mld_precset(prec,mld_coarse_sweeps_, prectype%cjswp, info) + call mld_precset(prec,mld_coarse_aggr_size_, prectype%csize, info) else nlv = 1 call mld_precinit(prec,prectype%prec, info, nlev=nlv) @@ -336,6 +338,7 @@ contains call read_data(prectype%cthres,5) ! Threshold for fact. 1 ILU(T) call read_data(prectype%cjswp,5) ! Jacobi sweeps call read_data(prectype%athres,5) ! smoother aggr thresh + call read_data(prectype%csize,5) ! coarse size end if end if @@ -373,6 +376,7 @@ contains call psb_bcast(ictxt,prectype%cthres) ! Threshold for fact. 1 ILU(T) call psb_bcast(ictxt,prectype%cjswp) ! Jacobi sweeps call psb_bcast(ictxt,prectype%athres) ! smoother aggr thresh + call psb_bcast(ictxt,prectype%csize) ! coarse size end if if (iam == psb_root_) then diff --git a/tests/pdegen/spde3d.f90 b/tests/pdegen/spde3d.f90 index 32b12f88..8d2ab5c4 100644 --- a/tests/pdegen/spde3d.f90 +++ b/tests/pdegen/spde3d.f90 @@ -111,6 +111,7 @@ program spde3d character(len=16) :: aggr_alg ! local or global aggregation character(len=16) :: mltype ! additive or multiplicative 2nd level prec character(len=16) :: smthpos ! side: pre, post, both smoothing + integer :: csize ! aggregation size at which to stop. character(len=16) :: cmat ! coarse mat character(len=16) :: csolve ! Coarse solver: bjac, umf, slu, sludist character(len=16) :: csbsolve ! Coarse subsolver: ILU, ILU(T), SuperLU, UMFPACK. @@ -200,6 +201,7 @@ program spde3d call mld_precset(prec,mld_coarse_fillin_, prectype%cfill, info) call mld_precset(prec,mld_coarse_iluthrs_, prectype%cthres, info) call mld_precset(prec,mld_coarse_sweeps_, prectype%cjswp, info) + call mld_precset(prec,mld_coarse_aggr_size_, prectype%csize, info) else nlv = 1 call mld_precinit(prec,prectype%prec, info, nlev=nlv) @@ -340,6 +342,7 @@ contains call read_data(prectype%cthres,5) ! Threshold for fact. 1 ILU(T) call read_data(prectype%cjswp,5) ! Jacobi sweeps call read_data(prectype%athres,5) ! smoother aggr thresh + call read_data(prectype%csize,5) ! coarse size end if end if @@ -377,6 +380,7 @@ contains call psb_bcast(ictxt,prectype%cthres) ! Threshold for fact. 1 ILU(T) call psb_bcast(ictxt,prectype%cjswp) ! Jacobi sweeps call psb_bcast(ictxt,prectype%athres) ! smoother aggr thresh + call psb_bcast(ictxt,prectype%csize) ! coarse size end if if (iam == psb_root_) then