From 1e484f6a189d7befb77695c1406f24d5b783ed7c Mon Sep 17 00:00:00 2001 From: Salvatore Filippone Date: Tue, 14 May 2019 16:02:20 +0100 Subject: [PATCH] Aligned new aggregator object internals. --- mlprec/impl/aggregator/Makefile | 22 +- .../mld_c_dec_aggregator_mat_asb.f90 | 229 +++++++++--------- .../mld_c_dec_aggregator_mat_bld.f90 | 216 +++++++++++++++++ ...t_biz_asb.f90 => mld_caggrmat_biz_bld.f90} | 14 +- ...rg_asb.f90 => mld_caggrmat_minnrg_bld.f90} | 12 +- ...th_asb.f90 => mld_caggrmat_nosmth_bld.f90} | 14 +- ...smth_asb.f90 => mld_caggrmat_smth_bld.f90} | 14 +- .../mld_d_dec_aggregator_mat_asb.f90 | 229 +++++++++--------- .../mld_d_dec_aggregator_mat_bld.f90 | 216 +++++++++++++++++ ...t_biz_asb.f90 => mld_daggrmat_biz_bld.f90} | 14 +- ...rg_asb.f90 => mld_daggrmat_minnrg_bld.f90} | 12 +- ...th_asb.f90 => mld_daggrmat_nosmth_bld.f90} | 14 +- ...smth_asb.f90 => mld_daggrmat_smth_bld.f90} | 14 +- .../mld_s_dec_aggregator_mat_asb.f90 | 229 +++++++++--------- .../mld_s_dec_aggregator_mat_bld.f90 | 216 +++++++++++++++++ ...t_biz_asb.f90 => mld_saggrmat_biz_bld.f90} | 14 +- ...rg_asb.f90 => mld_saggrmat_minnrg_bld.f90} | 12 +- ...th_asb.f90 => mld_saggrmat_nosmth_bld.f90} | 14 +- ...smth_asb.f90 => mld_saggrmat_smth_bld.f90} | 14 +- .../mld_z_dec_aggregator_mat_asb.f90 | 229 +++++++++--------- .../mld_z_dec_aggregator_mat_bld.f90 | 216 +++++++++++++++++ ...t_biz_asb.f90 => mld_zaggrmat_biz_bld.f90} | 14 +- ...rg_asb.f90 => mld_zaggrmat_minnrg_bld.f90} | 12 +- ...th_asb.f90 => mld_zaggrmat_nosmth_bld.f90} | 14 +- ...smth_asb.f90 => mld_zaggrmat_smth_bld.f90} | 14 +- .../impl/level/mld_c_base_onelev_mat_asb.f90 | 121 ++------- .../impl/level/mld_d_base_onelev_mat_asb.f90 | 121 ++------- .../impl/level/mld_s_base_onelev_mat_asb.f90 | 121 ++------- .../impl/level/mld_z_base_onelev_mat_asb.f90 | 121 ++------- mlprec/mld_c_base_aggregator_mod.f90 | 63 ++++- mlprec/mld_c_dec_aggregator_mod.f90 | 25 +- mlprec/mld_c_inner_mod.f90 | 38 +-- mlprec/mld_c_onelev_mod.f90 | 4 +- mlprec/mld_c_symdec_aggregator_mod.f90 | 4 +- mlprec/mld_d_base_aggregator_mod.f90 | 63 ++++- mlprec/mld_d_dec_aggregator_mod.f90 | 25 +- mlprec/mld_d_inner_mod.f90 | 38 +-- mlprec/mld_d_onelev_mod.f90 | 4 +- mlprec/mld_d_symdec_aggregator_mod.f90 | 4 +- mlprec/mld_s_base_aggregator_mod.f90 | 63 ++++- mlprec/mld_s_dec_aggregator_mod.f90 | 25 +- mlprec/mld_s_inner_mod.f90 | 38 +-- mlprec/mld_s_onelev_mod.f90 | 4 +- mlprec/mld_s_symdec_aggregator_mod.f90 | 4 +- mlprec/mld_z_base_aggregator_mod.f90 | 63 ++++- mlprec/mld_z_dec_aggregator_mod.f90 | 25 +- mlprec/mld_z_inner_mod.f90 | 38 +-- mlprec/mld_z_onelev_mod.f90 | 4 +- mlprec/mld_z_symdec_aggregator_mod.f90 | 4 +- 49 files changed, 1849 insertions(+), 1189 deletions(-) create mode 100644 mlprec/impl/aggregator/mld_c_dec_aggregator_mat_bld.f90 rename mlprec/impl/aggregator/{mld_caggrmat_biz_asb.f90 => mld_caggrmat_biz_bld.f90} (97%) rename mlprec/impl/aggregator/{mld_caggrmat_minnrg_asb.f90 => mld_caggrmat_minnrg_bld.f90} (98%) rename mlprec/impl/aggregator/{mld_caggrmat_nosmth_asb.f90 => mld_caggrmat_nosmth_bld.f90} (96%) rename mlprec/impl/aggregator/{mld_caggrmat_smth_asb.f90 => mld_caggrmat_smth_bld.f90} (97%) create mode 100644 mlprec/impl/aggregator/mld_d_dec_aggregator_mat_bld.f90 rename mlprec/impl/aggregator/{mld_daggrmat_biz_asb.f90 => mld_daggrmat_biz_bld.f90} (97%) rename mlprec/impl/aggregator/{mld_daggrmat_minnrg_asb.f90 => mld_daggrmat_minnrg_bld.f90} (98%) rename mlprec/impl/aggregator/{mld_daggrmat_nosmth_asb.f90 => mld_daggrmat_nosmth_bld.f90} (96%) rename mlprec/impl/aggregator/{mld_daggrmat_smth_asb.f90 => mld_daggrmat_smth_bld.f90} (97%) create mode 100644 mlprec/impl/aggregator/mld_s_dec_aggregator_mat_bld.f90 rename mlprec/impl/aggregator/{mld_saggrmat_biz_asb.f90 => mld_saggrmat_biz_bld.f90} (97%) rename mlprec/impl/aggregator/{mld_saggrmat_minnrg_asb.f90 => mld_saggrmat_minnrg_bld.f90} (98%) rename mlprec/impl/aggregator/{mld_saggrmat_nosmth_asb.f90 => mld_saggrmat_nosmth_bld.f90} (96%) rename mlprec/impl/aggregator/{mld_saggrmat_smth_asb.f90 => mld_saggrmat_smth_bld.f90} (97%) create mode 100644 mlprec/impl/aggregator/mld_z_dec_aggregator_mat_bld.f90 rename mlprec/impl/aggregator/{mld_zaggrmat_biz_asb.f90 => mld_zaggrmat_biz_bld.f90} (97%) rename mlprec/impl/aggregator/{mld_zaggrmat_minnrg_asb.f90 => mld_zaggrmat_minnrg_bld.f90} (98%) rename mlprec/impl/aggregator/{mld_zaggrmat_nosmth_asb.f90 => mld_zaggrmat_nosmth_bld.f90} (96%) rename mlprec/impl/aggregator/{mld_zaggrmat_smth_asb.f90 => mld_zaggrmat_smth_bld.f90} (97%) diff --git a/mlprec/impl/aggregator/Makefile b/mlprec/impl/aggregator/Makefile index 81a51e87..ebaead0f 100644 --- a/mlprec/impl/aggregator/Makefile +++ b/mlprec/impl/aggregator/Makefile @@ -10,29 +10,33 @@ FINCLUDES=$(FMFLAG)$(HERE) $(FMFLAG)$(MODDIR) $(FMFLAG)$(INCDIR) $(PSBLAS_INCLUD OBJS= \ mld_s_dec_aggregator_mat_asb.o \ +mld_s_dec_aggregator_mat_bld.o \ mld_s_dec_aggregator_tprol.o \ mld_s_symdec_aggregator_tprol.o \ mld_s_map_to_tprol.o mld_s_soc1_map_bld.o mld_s_soc2_map_bld.o\ -mld_saggrmat_biz_asb.o mld_saggrmat_minnrg_asb.o\ -mld_saggrmat_nosmth_asb.o mld_saggrmat_smth_asb.o \ +mld_saggrmat_biz_bld.o mld_saggrmat_minnrg_bld.o\ +mld_saggrmat_nosmth_bld.o mld_saggrmat_smth_bld.o \ mld_d_dec_aggregator_mat_asb.o \ +mld_d_dec_aggregator_mat_bld.o \ mld_d_dec_aggregator_tprol.o \ mld_d_symdec_aggregator_tprol.o \ mld_d_map_to_tprol.o mld_d_soc1_map_bld.o mld_d_soc2_map_bld.o\ -mld_daggrmat_biz_asb.o mld_daggrmat_minnrg_asb.o\ -mld_daggrmat_nosmth_asb.o mld_daggrmat_smth_asb.o \ +mld_daggrmat_biz_bld.o mld_daggrmat_minnrg_bld.o\ +mld_daggrmat_nosmth_bld.o mld_daggrmat_smth_bld.o \ mld_c_dec_aggregator_mat_asb.o \ +mld_c_dec_aggregator_mat_bld.o \ mld_c_dec_aggregator_tprol.o \ mld_c_symdec_aggregator_tprol.o \ mld_c_map_to_tprol.o mld_c_soc1_map_bld.o mld_c_soc2_map_bld.o\ -mld_caggrmat_biz_asb.o mld_caggrmat_minnrg_asb.o\ -mld_caggrmat_nosmth_asb.o mld_caggrmat_smth_asb.o \ +mld_caggrmat_biz_bld.o mld_caggrmat_minnrg_bld.o\ +mld_caggrmat_nosmth_bld.o mld_caggrmat_smth_bld.o \ mld_z_dec_aggregator_mat_asb.o \ +mld_z_dec_aggregator_mat_bld.o \ mld_z_dec_aggregator_tprol.o \ mld_z_symdec_aggregator_tprol.o \ mld_z_map_to_tprol.o mld_z_soc1_map_bld.o mld_z_soc2_map_bld.o\ -mld_zaggrmat_biz_asb.o mld_zaggrmat_minnrg_asb.o\ -mld_zaggrmat_nosmth_asb.o mld_zaggrmat_smth_asb.o +mld_zaggrmat_biz_bld.o mld_zaggrmat_minnrg_bld.o\ +mld_zaggrmat_nosmth_bld.o mld_zaggrmat_smth_bld.o #mld_s_hybrid_aggregator_tprol.o \ #mld_d_hybrid_aggregator_tprol.o \ @@ -40,7 +44,7 @@ mld_zaggrmat_nosmth_asb.o mld_zaggrmat_smth_asb.o #mld_z_hybrid_aggregator_tprol.o \ #bootCMatch_interface.o mld_d_bcmatch_aggregator_tprol.o\ -#mld_d_bcmatch_map_to_tprol.o mld_d_bcmatch_aggregator_mat_asb.o \ +#mld_d_bcmatch_map_to_tprol.o mld_d_bcmatch_aggregator_mat_bld.o \ LIBNAME=libmld_prec.a diff --git a/mlprec/impl/aggregator/mld_c_dec_aggregator_mat_asb.f90 b/mlprec/impl/aggregator/mld_c_dec_aggregator_mat_asb.f90 index f6f0717f..7335ce8f 100644 --- a/mlprec/impl/aggregator/mld_c_dec_aggregator_mat_asb.f90 +++ b/mlprec/impl/aggregator/mld_c_dec_aggregator_mat_asb.f90 @@ -40,63 +40,8 @@ ! Subroutine: mld_c_dec_aggregator_mat_asb ! Version: complex ! -! This routine builds the matrix associated to the current level of the -! multilevel preconditioner from the matrix associated to the previous level, -! by using the user-specified aggregation technique (therefore, it also builds the -! prolongation and restriction operators mapping the current level to the -! previous one and vice versa). -! The current level is regarded as the coarse one, while the previous as -! the fine one. This is in agreement with the fact that the routine is called, -! by mld_mlprec_bld, only on levels >=2. -! The coarse-level matrix A_C is built from a fine-level matrix A -! by using the Galerkin approach, i.e. -! -! A_C = P_C^T A P_C, -! -! where P_C is a prolongator from the coarse level to the fine one. -! -! A mapping from the nodes of the adjacency graph of A to the nodes of the -! adjacency graph of A_C has been computed by the mld_aggrmap_bld subroutine. -! The prolongator P_C is built here from this mapping, according to the -! value of p%iprcparm(mld_aggr_kind_), specified by the user through -! mld_cprecinit and mld_zprecset. -! On output from this routine the entries of AC, op_prol, op_restr -! are still in "global numbering" mode; this is fixed in the calling routine -! mld_c_lev_aggrmat_asb. -! -! Currently four different prolongators are implemented, corresponding to -! four aggregation algorithms: -! 1. un-smoothed aggregation, -! 2. smoothed aggregation, -! 3. "bizarre" aggregation. -! 4. minimum energy -! 1. The non-smoothed aggregation uses as prolongator the piecewise constant -! interpolation operator corresponding to the fine-to-coarse level mapping built -! by p%aggr%bld_tprol. This is called tentative prolongator. -! 2. The smoothed aggregation uses as prolongator the operator obtained by applying -! a damped Jacobi smoother to the tentative prolongator. -! 3. The "bizarre" aggregation uses a prolongator proposed by the authors of MLD2P4. -! This prolongator still requires a deep analysis and testing and its use is -! not recommended. -! 4. Minimum energy aggregation -! -! For more details see -! M. Brezina and P. Vanek, A black-box iterative solver based on a two-level -! Schwarz method, Computing, 63 (1999), 233-263. -! P. Vanek, J. Mandel and M. Brezina, Algebraic Multigrid by Smoothed -! Aggregation for Second and Fourth Order Elliptic Problems, Computing, 56 -! (1996), 179-196. -! P. D'Ambra, D. di Serafino and S. Filippone, On the development of PSBLAS-based -! parallel two-level Schwarz preconditioners, Appl. Num. Math., 57 (2007), -! 1181-1196. -! M. Sala, R. Tuminaro: A new Petrov-Galerkin smoothed aggregation preconditioner -! for nonsymmetric linear systems, SIAM J. Sci. Comput., 31(1):143-166 (2008) -! -! -! The main structure is: -! 1. Perform sanity checks; -! 2. Compute prolongator/restrictor/AC ! +! From a given AC to final format, generating DESC_AC ! ! Arguments: ! ag - type(mld_c_dec_aggregator_type), input/output. @@ -121,89 +66,143 @@ ! the various processes do not overlap. ! nlaggr - integer, dimension(:) input ! nlaggr(i) contains the aggregates held by process i. -! ac - type(psb_cspmat_type), output -! The coarse matrix on output +! ac - type(psb_cspmat_type), inout +! The coarse matrix +! desc_ac - type(psb_desc_type), output. +! The communication descriptor of the fine-level matrix. +! The 'one-level' data structure that will contain the local +! part of the matrix to be built as well as the information +! concerning the prolongator and its transpose. ! ! op_prol - type(psb_cspmat_type), input/output ! The tentative prolongator on input, the computed prolongator on output ! -! op_restr - type(psb_cspmat_type), output +! op_restr - type(psb_cspmat_type), input/output ! The restrictor operator; normally, it is the transpose of the prolongator. ! ! info - integer, output. ! Error code. ! -subroutine mld_c_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_prol,op_restr,info) +subroutine mld_c_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac, op_prol,op_restr,info) use psb_base_mod - use mld_c_prec_type, mld_protect_name => mld_c_dec_aggregator_mat_asb - use mld_c_inner_mod + use mld_base_prec_type + use mld_c_dec_aggregator_mod, mld_protect_name => mld_c_dec_aggregator_mat_asb implicit none - class(mld_c_dec_aggregator_type), target, intent(inout) :: ag - type(mld_sml_parms), intent(inout) :: parms + type(mld_sml_parms), intent(inout) :: parms type(psb_cspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) - type(psb_lcspmat_type), intent(inout) :: op_prol - type(psb_lcspmat_type), intent(out) :: ac,op_restr - integer(psb_ipk_), intent(out) :: info - - ! Local variables - character(len=20) :: name - integer(psb_mpk_) :: ictxt, np, me - type(psb_lc_coo_sparse_mat) :: acoo, bcoo - type(psb_lc_csr_sparse_mat) :: acsr1 - integer(psb_lpk_) :: nzl,ntaggr - integer(psb_ipk_) :: err_act - integer(psb_ipk_) :: debug_level, debug_unit - - name='mld_c_dec_aggregator_mat_asb' + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lcspmat_type), intent(inout) :: op_prol, ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + ! + integer(psb_mpk_) :: ictxt, np, me + type(psb_lc_coo_sparse_mat) :: acoo, bcoo + type(psb_lc_csr_sparse_mat) :: acsr1 + type(psb_lcspmat_type) :: lac, lac1 + type(psb_cspmat_type) :: tmp_ac + integer(psb_ipk_) :: i_nr, i_nc, i_nl, nzl + integer(psb_lpk_) :: ntaggr, l_nr, l_nc + integer(psb_ipk_) :: err_act, debug_level, debug_unit + character(len=20) :: name='c_dec_aggregator_mat_asb' + + + if (psb_get_errstatus().ne.0) return call psb_erractionsave(err_act) - if (psb_errstatus_fatal()) then - info = psb_err_internal_error_; goto 9999 - end if debug_unit = psb_get_debug_unit() debug_level = psb_get_debug_level() info = psb_success_ ictxt = desc_a%get_context() call psb_info(ictxt,me,np) - ! - ! Build the coarse-level matrix from the fine-level one, starting from - ! the mapping defined by mld_aggrmap_bld and applying the aggregation - ! algorithm specified by - ! - select case (parms%aggr_prol) - case (mld_no_smooth_) - - call mld_caggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,& - & parms,ac,op_prol,op_restr,info) - - case(mld_smooth_prol_) - - call mld_caggrmat_smth_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - - case(mld_biz_prol_) - call mld_caggrmat_biz_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - case(mld_min_energy_) - - call mld_caggrmat_minnrg_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - - case default + ntaggr = sum(nlaggr) + + select case(parms%coarse_mat) + + case(mld_distr_mat_) + + call ac%mv_to(bcoo) + nzl = bcoo%get_nzeros() + i_nl = nlaggr(me+1) + if (info == psb_success_) call psb_cdall(ictxt,desc_ac,info,nl=i_nl) + if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,desc_ac,info) + if (info == psb_success_) call psb_cdasb(desc_ac,info) + if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),desc_ac,info,iact='I') + if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),desc_ac,info,iact='I') + if (info /= psb_success_) then + call psb_errpush(psb_err_internal_error_,name,& + & a_err='Creating desc_ac and converting ac') + goto 9999 + end if + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),& + & 'Assembld aux descr. distr.' + call ac%mv_from(bcoo) + l_nr = desc_ac%get_local_rows() + l_nc = desc_ac%get_local_cols() + call ac%set_nrows(l_nr) + call ac%set_ncols(l_nc) + call ac%set_asb() + + if (info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free') + goto 9999 + end if + + if (np>1) then + call op_prol%mv_to(acsr1) + nzl = acsr1%get_nzeros() + call psb_glob_to_loc(acsr1%ja(1:nzl),desc_ac,info,'I') + if(info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc') + goto 9999 + end if + call op_prol%mv_from(acsr1) + endif + call op_prol%set_ncols(l_nc) + + if (np>1) then + call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_) + call op_restr%mv_to(acoo) + nzl = acoo%get_nzeros() + if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),desc_ac,info,'I') + call acoo%set_dupl(psb_dupl_add_) + if (info == psb_success_) call op_restr%mv_from(acoo) + if (info == psb_success_) call op_restr%cscnv(info,type='csr') + if(info /= psb_success_) then + call psb_errpush(psb_err_internal_error_,name,& + & a_err='Converting op_restr to local') + goto 9999 + end if + end if + ! + ! Clip to local rows. + ! + call op_restr%set_nrows(l_nr) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),& + & 'Done ac ' + + case(mld_repl_mat_) + ! + ! + call psb_cdall(ictxt,desc_ac,info,mg=ntaggr,repl=.true.) + if (info == psb_success_) call psb_cdasb(desc_ac,info) + if (info == psb_success_) call tmp_ac%mv_from_l(ac) + if (info == psb_success_) & + & call psb_gather(ac,tmp_ac,desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) + if (info /= psb_success_) goto 9999 + + case default info = psb_err_internal_error_ - call psb_errpush(info,name,a_err='Invalid aggr kind') + call psb_errpush(info,name,a_err='invalid mld_coarse_mat_') goto 9999 - end select - if (info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='Inner aggrmat asb') - goto 9999 - end if call psb_erractionrestore(err_act) @@ -212,5 +211,5 @@ subroutine mld_c_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_p 9999 call psb_error_handler(err_act) return - + end subroutine mld_c_dec_aggregator_mat_asb diff --git a/mlprec/impl/aggregator/mld_c_dec_aggregator_mat_bld.f90 b/mlprec/impl/aggregator/mld_c_dec_aggregator_mat_bld.f90 new file mode 100644 index 00000000..1a73cca2 --- /dev/null +++ b/mlprec/impl/aggregator/mld_c_dec_aggregator_mat_bld.f90 @@ -0,0 +1,216 @@ +! +! +! MLD2P4 version 2.2 +! MultiLevel Domain Decomposition Parallel Preconditioners Package +! based on PSBLAS (Parallel Sparse BLAS version 3.5) +! +! (C) Copyright 2008-2018 +! +! Salvatore Filippone +! Pasqua D'Ambra +! Daniela di Serafino +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the MLD2P4 group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE MLD2P4 GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: mld_c_dec_aggregator_mat_bld.f90 +! +! Subroutine: mld_c_dec_aggregator_mat_bld +! Version: complex +! +! This routine builds the matrix associated to the current level of the +! multilevel preconditioner from the matrix associated to the previous level, +! by using the user-specified aggregation technique (therefore, it also builds the +! prolongation and restriction operators mapping the current level to the +! previous one and vice versa). +! The current level is regarded as the coarse one, while the previous as +! the fine one. This is in agreement with the fact that the routine is called, +! by mld_mlprec_bld, only on levels >=2. +! The coarse-level matrix A_C is built from a fine-level matrix A +! by using the Galerkin approach, i.e. +! +! A_C = P_C^T A P_C, +! +! where P_C is a prolongator from the coarse level to the fine one. +! +! A mapping from the nodes of the adjacency graph of A to the nodes of the +! adjacency graph of A_C has been computed by the mld_aggrmap_bld subroutine. +! The prolongator P_C is built here from this mapping, according to the +! value of p%iprcparm(mld_aggr_kind_), specified by the user through +! mld_cprecinit and mld_zprecset. +! On output from this routine the entries of AC, op_prol, op_restr +! are still in "global numbering" mode; this is fixed in the calling routine +! mld_c_lev_aggrmat_bld. +! +! Currently four different prolongators are implemented, corresponding to +! four aggregation algorithms: +! 1. un-smoothed aggregation, +! 2. smoothed aggregation, +! 3. "bizarre" aggregation. +! 4. minimum energy +! 1. The non-smoothed aggregation uses as prolongator the piecewise constant +! interpolation operator corresponding to the fine-to-coarse level mapping built +! by p%aggr%bld_tprol. This is called tentative prolongator. +! 2. The smoothed aggregation uses as prolongator the operator obtained by applying +! a damped Jacobi smoother to the tentative prolongator. +! 3. The "bizarre" aggregation uses a prolongator proposed by the authors of MLD2P4. +! This prolongator still requires a deep analysis and testing and its use is +! not recommended. +! 4. Minimum energy aggregation +! +! For more details see +! M. Brezina and P. Vanek, A black-box iterative solver based on a two-level +! Schwarz method, Computing, 63 (1999), 233-263. +! P. Vanek, J. Mandel and M. Brezina, Algebraic Multigrid by Smoothed +! Aggregation for Second and Fourth Order Elliptic Problems, Computing, 56 +! (1996), 179-196. +! P. D'Ambra, D. di Serafino and S. Filippone, On the development of PSBLAS-based +! parallel two-level Schwarz preconditioners, Appl. Num. Math., 57 (2007), +! 1181-1196. +! M. Sala, R. Tuminaro: A new Petrov-Galerkin smoothed aggregation preconditioner +! for nonsymmetric linear systems, SIAM J. Sci. Comput., 31(1):143-166 (2008) +! +! +! The main structure is: +! 1. Perform sanity checks; +! 2. Compute prolongator/restrictor/AC +! +! +! Arguments: +! ag - type(mld_c_dec_aggregator_type), input/output. +! The aggregator object +! parms - type(mld_sml_parms), input +! The aggregation parameters +! a - type(psb_cspmat_type), input. +! The sparse matrix structure containing the local part of +! the fine-level matrix. +! desc_a - type(psb_desc_type), input. +! The communication descriptor of the fine-level matrix. +! The 'one-level' data structure that will contain the local +! part of the matrix to be built as well as the information +! concerning the prolongator and its transpose. +! ilaggr - integer, dimension(:), input +! The mapping between the row indices of the coarse-level +! matrix and the row indices of the fine-level matrix. +! ilaggr(i)=j means that node i in the adjacency graph +! of the fine-level matrix is mapped onto node j in the +! adjacency graph of the coarse-level matrix. Note that the indices +! are assumed to be shifted so as to make sure the ranges on +! the various processes do not overlap. +! nlaggr - integer, dimension(:) input +! nlaggr(i) contains the aggregates held by process i. +! ac - type(psb_cspmat_type), output +! The coarse matrix on output +! +! op_prol - type(psb_cspmat_type), input/output +! The tentative prolongator on input, the computed prolongator on output +! +! op_restr - type(psb_cspmat_type), output +! The restrictor operator; normally, it is the transpose of the prolongator. +! +! info - integer, output. +! Error code. +! +subroutine mld_c_dec_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_prol,op_restr,info) + use psb_base_mod + use mld_c_prec_type, mld_protect_name => mld_c_dec_aggregator_mat_bld + use mld_c_inner_mod + implicit none + + class(mld_c_dec_aggregator_type), target, intent(inout) :: ag + type(mld_sml_parms), intent(inout) :: parms + type(psb_cspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lcspmat_type), intent(inout) :: op_prol + type(psb_lcspmat_type), intent(out) :: ac,op_restr + integer(psb_ipk_), intent(out) :: info + + ! Local variables + character(len=20) :: name + integer(psb_mpk_) :: ictxt, np, me + type(psb_lc_coo_sparse_mat) :: acoo, bcoo + type(psb_lc_csr_sparse_mat) :: acsr1 + integer(psb_lpk_) :: nzl,ntaggr + integer(psb_ipk_) :: err_act + integer(psb_ipk_) :: debug_level, debug_unit + + name='mld_c_dec_aggregator_mat_bld' + call psb_erractionsave(err_act) + if (psb_errstatus_fatal()) then + info = psb_err_internal_error_; goto 9999 + end if + debug_unit = psb_get_debug_unit() + debug_level = psb_get_debug_level() + info = psb_success_ + ictxt = desc_a%get_context() + call psb_info(ictxt,me,np) + + ! + ! Build the coarse-level matrix from the fine-level one, starting from + ! the mapping defined by mld_aggrmap_bld and applying the aggregation + ! algorithm specified by + ! + select case (parms%aggr_prol) + case (mld_no_smooth_) + + call mld_caggrmat_nosmth_bld(a,desc_a,ilaggr,nlaggr,& + & parms,ac,op_prol,op_restr,info) + + case(mld_smooth_prol_) + + call mld_caggrmat_smth_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case(mld_biz_prol_) + + call mld_caggrmat_biz_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case(mld_min_energy_) + + call mld_caggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case default + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Invalid aggr kind') + goto 9999 + + end select + if (info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='Inner aggrmat bld') + goto 9999 + end if + + + call psb_erractionrestore(err_act) + return + +9999 call psb_error_handler(err_act) + return + + +end subroutine mld_c_dec_aggregator_mat_bld diff --git a/mlprec/impl/aggregator/mld_caggrmat_biz_asb.f90 b/mlprec/impl/aggregator/mld_caggrmat_biz_bld.f90 similarity index 97% rename from mlprec/impl/aggregator/mld_caggrmat_biz_asb.f90 rename to mlprec/impl/aggregator/mld_caggrmat_biz_bld.f90 index 3d27c449..10f1d0cd 100644 --- a/mlprec/impl/aggregator/mld_caggrmat_biz_asb.f90 +++ b/mlprec/impl/aggregator/mld_caggrmat_biz_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_caggrmat_biz_asb.F90 +! File: mld_caggrmat_biz_bld.F90 ! -! Subroutine: mld_caggrmat_biz_asb +! Subroutine: mld_caggrmat_biz_bld ! Version: complex ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -57,7 +57,7 @@ ! specified by the user through mld_cprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! mld_c_lev_aggrmat_asb. +! mld_c_lev_aggrmat_bld. ! ! Arguments: ! a - type(psb_cspmat_type), input. @@ -80,10 +80,10 @@ ! info - integer, output. ! Error code. ! -subroutine mld_caggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_caggrmat_biz_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_c_inner_mod, mld_protect_name => mld_caggrmat_biz_asb + use mld_c_inner_mod, mld_protect_name => mld_caggrmat_biz_bld implicit none @@ -111,7 +111,7 @@ subroutine mld_caggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr integer(psb_ipk_), parameter :: ncmax=16 real(psb_spk_) :: anorm, omega, tmp, dg, theta - name='mld_aggrmat_biz_asb' + name='mld_aggrmat_biz_bld' info=psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -372,4 +372,4 @@ subroutine mld_caggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr call psb_error_handler(err_act) return -end subroutine mld_caggrmat_biz_asb +end subroutine mld_caggrmat_biz_bld diff --git a/mlprec/impl/aggregator/mld_caggrmat_minnrg_asb.f90 b/mlprec/impl/aggregator/mld_caggrmat_minnrg_bld.f90 similarity index 98% rename from mlprec/impl/aggregator/mld_caggrmat_minnrg_asb.f90 rename to mlprec/impl/aggregator/mld_caggrmat_minnrg_bld.f90 index 12eb9379..21365f6f 100644 --- a/mlprec/impl/aggregator/mld_caggrmat_minnrg_asb.f90 +++ b/mlprec/impl/aggregator/mld_caggrmat_minnrg_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_caggrmat_minnrg_asb.F90 +! File: mld_caggrmat_minnrg_bld.F90 ! -! Subroutine: mld_caggrmat_minnrg_asb +! Subroutine: mld_caggrmat_minnrg_bld ! Version: complex ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -65,7 +65,7 @@ ! ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! ! Arguments: @@ -104,10 +104,10 @@ ! Error code. ! ! -subroutine mld_caggrmat_minnrg_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_caggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_c_inner_mod, mld_protect_name => mld_caggrmat_minnrg_asb + use mld_c_inner_mod, mld_protect_name => mld_caggrmat_minnrg_bld implicit none @@ -651,4 +651,4 @@ contains close(20+me) end subroutine local_dump -end subroutine mld_caggrmat_minnrg_asb +end subroutine mld_caggrmat_minnrg_bld diff --git a/mlprec/impl/aggregator/mld_caggrmat_nosmth_asb.f90 b/mlprec/impl/aggregator/mld_caggrmat_nosmth_bld.f90 similarity index 96% rename from mlprec/impl/aggregator/mld_caggrmat_nosmth_asb.f90 rename to mlprec/impl/aggregator/mld_caggrmat_nosmth_bld.f90 index 7eb18dcc..a8f5dc6e 100644 --- a/mlprec/impl/aggregator/mld_caggrmat_nosmth_asb.f90 +++ b/mlprec/impl/aggregator/mld_caggrmat_nosmth_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_caggrmat_nosmth_asb.F90 +! File: mld_caggrmat_nosmth_bld.F90 ! -! Subroutine: mld_caggrmat_nosmth_asb +! Subroutine: mld_caggrmat_nosmth_bld ! Version: complex ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -53,7 +53,7 @@ ! specified by the user through mld_cprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! For details see ! P. D'Ambra, D. di Serafino and S. Filippone, On the development of @@ -96,10 +96,10 @@ ! Error code. ! ! -subroutine mld_caggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_caggrmat_nosmth_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_c_inner_mod, mld_protect_name => mld_caggrmat_nosmth_asb + use mld_c_inner_mod, mld_protect_name => mld_caggrmat_nosmth_bld implicit none @@ -124,7 +124,7 @@ subroutine mld_caggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_re integer(psb_lpk_) :: nrow, nglob, ncol, ntaggr, nzl, ip, & & naggr, nzt, naggrm1, naggrp1, i, k - name = 'mld_aggrmat_nosmth_asb' + name = 'mld_aggrmat_nosmth_bld' info = psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -198,4 +198,4 @@ subroutine mld_caggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_re return -end subroutine mld_caggrmat_nosmth_asb +end subroutine mld_caggrmat_nosmth_bld diff --git a/mlprec/impl/aggregator/mld_caggrmat_smth_asb.f90 b/mlprec/impl/aggregator/mld_caggrmat_smth_bld.f90 similarity index 97% rename from mlprec/impl/aggregator/mld_caggrmat_smth_asb.f90 rename to mlprec/impl/aggregator/mld_caggrmat_smth_bld.f90 index a8564eb6..b63f7e06 100644 --- a/mlprec/impl/aggregator/mld_caggrmat_smth_asb.f90 +++ b/mlprec/impl/aggregator/mld_caggrmat_smth_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_caggrmat_smth_asb.F90 +! File: mld_caggrmat_smth_bld.F90 ! -! Subroutine: mld_caggrmat_smth_asb +! Subroutine: mld_caggrmat_smth_bld ! Version: complex ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -65,7 +65,7 @@ ! specified by the user through mld_cprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! ! Arguments: @@ -102,10 +102,10 @@ ! info - integer, output. ! Error code. ! -subroutine mld_caggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_caggrmat_smth_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_c_inner_mod, mld_protect_name => mld_caggrmat_smth_asb + use mld_c_inner_mod, mld_protect_name => mld_caggrmat_smth_bld implicit none @@ -133,7 +133,7 @@ subroutine mld_caggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_rest integer(psb_ipk_), parameter :: ncmax=16 real(psb_spk_) :: anorm, omega, tmp, dg, theta - name='mld_aggrmat_smth_asb' + name='mld_aggrmat_smth_bld' info=psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -411,4 +411,4 @@ subroutine mld_caggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_rest call psb_error_handler(err_act) return -end subroutine mld_caggrmat_smth_asb +end subroutine mld_caggrmat_smth_bld diff --git a/mlprec/impl/aggregator/mld_d_dec_aggregator_mat_asb.f90 b/mlprec/impl/aggregator/mld_d_dec_aggregator_mat_asb.f90 index 4004c189..82f9ee9e 100644 --- a/mlprec/impl/aggregator/mld_d_dec_aggregator_mat_asb.f90 +++ b/mlprec/impl/aggregator/mld_d_dec_aggregator_mat_asb.f90 @@ -40,63 +40,8 @@ ! Subroutine: mld_d_dec_aggregator_mat_asb ! Version: real ! -! This routine builds the matrix associated to the current level of the -! multilevel preconditioner from the matrix associated to the previous level, -! by using the user-specified aggregation technique (therefore, it also builds the -! prolongation and restriction operators mapping the current level to the -! previous one and vice versa). -! The current level is regarded as the coarse one, while the previous as -! the fine one. This is in agreement with the fact that the routine is called, -! by mld_mlprec_bld, only on levels >=2. -! The coarse-level matrix A_C is built from a fine-level matrix A -! by using the Galerkin approach, i.e. -! -! A_C = P_C^T A P_C, -! -! where P_C is a prolongator from the coarse level to the fine one. -! -! A mapping from the nodes of the adjacency graph of A to the nodes of the -! adjacency graph of A_C has been computed by the mld_aggrmap_bld subroutine. -! The prolongator P_C is built here from this mapping, according to the -! value of p%iprcparm(mld_aggr_kind_), specified by the user through -! mld_dprecinit and mld_zprecset. -! On output from this routine the entries of AC, op_prol, op_restr -! are still in "global numbering" mode; this is fixed in the calling routine -! mld_d_lev_aggrmat_asb. -! -! Currently four different prolongators are implemented, corresponding to -! four aggregation algorithms: -! 1. un-smoothed aggregation, -! 2. smoothed aggregation, -! 3. "bizarre" aggregation. -! 4. minimum energy -! 1. The non-smoothed aggregation uses as prolongator the piecewise constant -! interpolation operator corresponding to the fine-to-coarse level mapping built -! by p%aggr%bld_tprol. This is called tentative prolongator. -! 2. The smoothed aggregation uses as prolongator the operator obtained by applying -! a damped Jacobi smoother to the tentative prolongator. -! 3. The "bizarre" aggregation uses a prolongator proposed by the authors of MLD2P4. -! This prolongator still requires a deep analysis and testing and its use is -! not recommended. -! 4. Minimum energy aggregation -! -! For more details see -! M. Brezina and P. Vanek, A black-box iterative solver based on a two-level -! Schwarz method, Computing, 63 (1999), 233-263. -! P. Vanek, J. Mandel and M. Brezina, Algebraic Multigrid by Smoothed -! Aggregation for Second and Fourth Order Elliptic Problems, Computing, 56 -! (1996), 179-196. -! P. D'Ambra, D. di Serafino and S. Filippone, On the development of PSBLAS-based -! parallel two-level Schwarz preconditioners, Appl. Num. Math., 57 (2007), -! 1181-1196. -! M. Sala, R. Tuminaro: A new Petrov-Galerkin smoothed aggregation preconditioner -! for nonsymmetric linear systems, SIAM J. Sci. Comput., 31(1):143-166 (2008) -! -! -! The main structure is: -! 1. Perform sanity checks; -! 2. Compute prolongator/restrictor/AC ! +! From a given AC to final format, generating DESC_AC ! ! Arguments: ! ag - type(mld_d_dec_aggregator_type), input/output. @@ -121,89 +66,143 @@ ! the various processes do not overlap. ! nlaggr - integer, dimension(:) input ! nlaggr(i) contains the aggregates held by process i. -! ac - type(psb_dspmat_type), output -! The coarse matrix on output +! ac - type(psb_dspmat_type), inout +! The coarse matrix +! desc_ac - type(psb_desc_type), output. +! The communication descriptor of the fine-level matrix. +! The 'one-level' data structure that will contain the local +! part of the matrix to be built as well as the information +! concerning the prolongator and its transpose. ! ! op_prol - type(psb_dspmat_type), input/output ! The tentative prolongator on input, the computed prolongator on output ! -! op_restr - type(psb_dspmat_type), output +! op_restr - type(psb_dspmat_type), input/output ! The restrictor operator; normally, it is the transpose of the prolongator. ! ! info - integer, output. ! Error code. ! -subroutine mld_d_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_prol,op_restr,info) +subroutine mld_d_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac, op_prol,op_restr,info) use psb_base_mod - use mld_d_prec_type, mld_protect_name => mld_d_dec_aggregator_mat_asb - use mld_d_inner_mod + use mld_base_prec_type + use mld_d_dec_aggregator_mod, mld_protect_name => mld_d_dec_aggregator_mat_asb implicit none - class(mld_d_dec_aggregator_type), target, intent(inout) :: ag - type(mld_dml_parms), intent(inout) :: parms + type(mld_dml_parms), intent(inout) :: parms type(psb_dspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) - type(psb_ldspmat_type), intent(inout) :: op_prol - type(psb_ldspmat_type), intent(out) :: ac,op_restr - integer(psb_ipk_), intent(out) :: info - - ! Local variables - character(len=20) :: name - integer(psb_mpk_) :: ictxt, np, me - type(psb_ld_coo_sparse_mat) :: acoo, bcoo - type(psb_ld_csr_sparse_mat) :: acsr1 - integer(psb_lpk_) :: nzl,ntaggr - integer(psb_ipk_) :: err_act - integer(psb_ipk_) :: debug_level, debug_unit - - name='mld_d_dec_aggregator_mat_asb' + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_ldspmat_type), intent(inout) :: op_prol, ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + ! + integer(psb_mpk_) :: ictxt, np, me + type(psb_ld_coo_sparse_mat) :: acoo, bcoo + type(psb_ld_csr_sparse_mat) :: acsr1 + type(psb_ldspmat_type) :: lac, lac1 + type(psb_dspmat_type) :: tmp_ac + integer(psb_ipk_) :: i_nr, i_nc, i_nl, nzl + integer(psb_lpk_) :: ntaggr, l_nr, l_nc + integer(psb_ipk_) :: err_act, debug_level, debug_unit + character(len=20) :: name='d_dec_aggregator_mat_asb' + + + if (psb_get_errstatus().ne.0) return call psb_erractionsave(err_act) - if (psb_errstatus_fatal()) then - info = psb_err_internal_error_; goto 9999 - end if debug_unit = psb_get_debug_unit() debug_level = psb_get_debug_level() info = psb_success_ ictxt = desc_a%get_context() call psb_info(ictxt,me,np) - ! - ! Build the coarse-level matrix from the fine-level one, starting from - ! the mapping defined by mld_aggrmap_bld and applying the aggregation - ! algorithm specified by - ! - select case (parms%aggr_prol) - case (mld_no_smooth_) - - call mld_daggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,& - & parms,ac,op_prol,op_restr,info) - - case(mld_smooth_prol_) - - call mld_daggrmat_smth_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - - case(mld_biz_prol_) - call mld_daggrmat_biz_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - case(mld_min_energy_) - - call mld_daggrmat_minnrg_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - - case default + ntaggr = sum(nlaggr) + + select case(parms%coarse_mat) + + case(mld_distr_mat_) + + call ac%mv_to(bcoo) + nzl = bcoo%get_nzeros() + i_nl = nlaggr(me+1) + if (info == psb_success_) call psb_cdall(ictxt,desc_ac,info,nl=i_nl) + if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,desc_ac,info) + if (info == psb_success_) call psb_cdasb(desc_ac,info) + if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),desc_ac,info,iact='I') + if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),desc_ac,info,iact='I') + if (info /= psb_success_) then + call psb_errpush(psb_err_internal_error_,name,& + & a_err='Creating desc_ac and converting ac') + goto 9999 + end if + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),& + & 'Assembld aux descr. distr.' + call ac%mv_from(bcoo) + l_nr = desc_ac%get_local_rows() + l_nc = desc_ac%get_local_cols() + call ac%set_nrows(l_nr) + call ac%set_ncols(l_nc) + call ac%set_asb() + + if (info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free') + goto 9999 + end if + + if (np>1) then + call op_prol%mv_to(acsr1) + nzl = acsr1%get_nzeros() + call psb_glob_to_loc(acsr1%ja(1:nzl),desc_ac,info,'I') + if(info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc') + goto 9999 + end if + call op_prol%mv_from(acsr1) + endif + call op_prol%set_ncols(l_nc) + + if (np>1) then + call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_) + call op_restr%mv_to(acoo) + nzl = acoo%get_nzeros() + if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),desc_ac,info,'I') + call acoo%set_dupl(psb_dupl_add_) + if (info == psb_success_) call op_restr%mv_from(acoo) + if (info == psb_success_) call op_restr%cscnv(info,type='csr') + if(info /= psb_success_) then + call psb_errpush(psb_err_internal_error_,name,& + & a_err='Converting op_restr to local') + goto 9999 + end if + end if + ! + ! Clip to local rows. + ! + call op_restr%set_nrows(l_nr) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),& + & 'Done ac ' + + case(mld_repl_mat_) + ! + ! + call psb_cdall(ictxt,desc_ac,info,mg=ntaggr,repl=.true.) + if (info == psb_success_) call psb_cdasb(desc_ac,info) + if (info == psb_success_) call tmp_ac%mv_from_l(ac) + if (info == psb_success_) & + & call psb_gather(ac,tmp_ac,desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) + if (info /= psb_success_) goto 9999 + + case default info = psb_err_internal_error_ - call psb_errpush(info,name,a_err='Invalid aggr kind') + call psb_errpush(info,name,a_err='invalid mld_coarse_mat_') goto 9999 - end select - if (info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='Inner aggrmat asb') - goto 9999 - end if call psb_erractionrestore(err_act) @@ -212,5 +211,5 @@ subroutine mld_d_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_p 9999 call psb_error_handler(err_act) return - + end subroutine mld_d_dec_aggregator_mat_asb diff --git a/mlprec/impl/aggregator/mld_d_dec_aggregator_mat_bld.f90 b/mlprec/impl/aggregator/mld_d_dec_aggregator_mat_bld.f90 new file mode 100644 index 00000000..b0317d76 --- /dev/null +++ b/mlprec/impl/aggregator/mld_d_dec_aggregator_mat_bld.f90 @@ -0,0 +1,216 @@ +! +! +! MLD2P4 version 2.2 +! MultiLevel Domain Decomposition Parallel Preconditioners Package +! based on PSBLAS (Parallel Sparse BLAS version 3.5) +! +! (C) Copyright 2008-2018 +! +! Salvatore Filippone +! Pasqua D'Ambra +! Daniela di Serafino +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the MLD2P4 group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE MLD2P4 GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: mld_d_dec_aggregator_mat_bld.f90 +! +! Subroutine: mld_d_dec_aggregator_mat_bld +! Version: real +! +! This routine builds the matrix associated to the current level of the +! multilevel preconditioner from the matrix associated to the previous level, +! by using the user-specified aggregation technique (therefore, it also builds the +! prolongation and restriction operators mapping the current level to the +! previous one and vice versa). +! The current level is regarded as the coarse one, while the previous as +! the fine one. This is in agreement with the fact that the routine is called, +! by mld_mlprec_bld, only on levels >=2. +! The coarse-level matrix A_C is built from a fine-level matrix A +! by using the Galerkin approach, i.e. +! +! A_C = P_C^T A P_C, +! +! where P_C is a prolongator from the coarse level to the fine one. +! +! A mapping from the nodes of the adjacency graph of A to the nodes of the +! adjacency graph of A_C has been computed by the mld_aggrmap_bld subroutine. +! The prolongator P_C is built here from this mapping, according to the +! value of p%iprcparm(mld_aggr_kind_), specified by the user through +! mld_dprecinit and mld_zprecset. +! On output from this routine the entries of AC, op_prol, op_restr +! are still in "global numbering" mode; this is fixed in the calling routine +! mld_d_lev_aggrmat_bld. +! +! Currently four different prolongators are implemented, corresponding to +! four aggregation algorithms: +! 1. un-smoothed aggregation, +! 2. smoothed aggregation, +! 3. "bizarre" aggregation. +! 4. minimum energy +! 1. The non-smoothed aggregation uses as prolongator the piecewise constant +! interpolation operator corresponding to the fine-to-coarse level mapping built +! by p%aggr%bld_tprol. This is called tentative prolongator. +! 2. The smoothed aggregation uses as prolongator the operator obtained by applying +! a damped Jacobi smoother to the tentative prolongator. +! 3. The "bizarre" aggregation uses a prolongator proposed by the authors of MLD2P4. +! This prolongator still requires a deep analysis and testing and its use is +! not recommended. +! 4. Minimum energy aggregation +! +! For more details see +! M. Brezina and P. Vanek, A black-box iterative solver based on a two-level +! Schwarz method, Computing, 63 (1999), 233-263. +! P. Vanek, J. Mandel and M. Brezina, Algebraic Multigrid by Smoothed +! Aggregation for Second and Fourth Order Elliptic Problems, Computing, 56 +! (1996), 179-196. +! P. D'Ambra, D. di Serafino and S. Filippone, On the development of PSBLAS-based +! parallel two-level Schwarz preconditioners, Appl. Num. Math., 57 (2007), +! 1181-1196. +! M. Sala, R. Tuminaro: A new Petrov-Galerkin smoothed aggregation preconditioner +! for nonsymmetric linear systems, SIAM J. Sci. Comput., 31(1):143-166 (2008) +! +! +! The main structure is: +! 1. Perform sanity checks; +! 2. Compute prolongator/restrictor/AC +! +! +! Arguments: +! ag - type(mld_d_dec_aggregator_type), input/output. +! The aggregator object +! parms - type(mld_dml_parms), input +! The aggregation parameters +! a - type(psb_dspmat_type), input. +! The sparse matrix structure containing the local part of +! the fine-level matrix. +! desc_a - type(psb_desc_type), input. +! The communication descriptor of the fine-level matrix. +! The 'one-level' data structure that will contain the local +! part of the matrix to be built as well as the information +! concerning the prolongator and its transpose. +! ilaggr - integer, dimension(:), input +! The mapping between the row indices of the coarse-level +! matrix and the row indices of the fine-level matrix. +! ilaggr(i)=j means that node i in the adjacency graph +! of the fine-level matrix is mapped onto node j in the +! adjacency graph of the coarse-level matrix. Note that the indices +! are assumed to be shifted so as to make sure the ranges on +! the various processes do not overlap. +! nlaggr - integer, dimension(:) input +! nlaggr(i) contains the aggregates held by process i. +! ac - type(psb_dspmat_type), output +! The coarse matrix on output +! +! op_prol - type(psb_dspmat_type), input/output +! The tentative prolongator on input, the computed prolongator on output +! +! op_restr - type(psb_dspmat_type), output +! The restrictor operator; normally, it is the transpose of the prolongator. +! +! info - integer, output. +! Error code. +! +subroutine mld_d_dec_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_prol,op_restr,info) + use psb_base_mod + use mld_d_prec_type, mld_protect_name => mld_d_dec_aggregator_mat_bld + use mld_d_inner_mod + implicit none + + class(mld_d_dec_aggregator_type), target, intent(inout) :: ag + type(mld_dml_parms), intent(inout) :: parms + type(psb_dspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_ldspmat_type), intent(inout) :: op_prol + type(psb_ldspmat_type), intent(out) :: ac,op_restr + integer(psb_ipk_), intent(out) :: info + + ! Local variables + character(len=20) :: name + integer(psb_mpk_) :: ictxt, np, me + type(psb_ld_coo_sparse_mat) :: acoo, bcoo + type(psb_ld_csr_sparse_mat) :: acsr1 + integer(psb_lpk_) :: nzl,ntaggr + integer(psb_ipk_) :: err_act + integer(psb_ipk_) :: debug_level, debug_unit + + name='mld_d_dec_aggregator_mat_bld' + call psb_erractionsave(err_act) + if (psb_errstatus_fatal()) then + info = psb_err_internal_error_; goto 9999 + end if + debug_unit = psb_get_debug_unit() + debug_level = psb_get_debug_level() + info = psb_success_ + ictxt = desc_a%get_context() + call psb_info(ictxt,me,np) + + ! + ! Build the coarse-level matrix from the fine-level one, starting from + ! the mapping defined by mld_aggrmap_bld and applying the aggregation + ! algorithm specified by + ! + select case (parms%aggr_prol) + case (mld_no_smooth_) + + call mld_daggrmat_nosmth_bld(a,desc_a,ilaggr,nlaggr,& + & parms,ac,op_prol,op_restr,info) + + case(mld_smooth_prol_) + + call mld_daggrmat_smth_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case(mld_biz_prol_) + + call mld_daggrmat_biz_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case(mld_min_energy_) + + call mld_daggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case default + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Invalid aggr kind') + goto 9999 + + end select + if (info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='Inner aggrmat bld') + goto 9999 + end if + + + call psb_erractionrestore(err_act) + return + +9999 call psb_error_handler(err_act) + return + + +end subroutine mld_d_dec_aggregator_mat_bld diff --git a/mlprec/impl/aggregator/mld_daggrmat_biz_asb.f90 b/mlprec/impl/aggregator/mld_daggrmat_biz_bld.f90 similarity index 97% rename from mlprec/impl/aggregator/mld_daggrmat_biz_asb.f90 rename to mlprec/impl/aggregator/mld_daggrmat_biz_bld.f90 index d1bb1233..30aa8ed4 100644 --- a/mlprec/impl/aggregator/mld_daggrmat_biz_asb.f90 +++ b/mlprec/impl/aggregator/mld_daggrmat_biz_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_daggrmat_biz_asb.F90 +! File: mld_daggrmat_biz_bld.F90 ! -! Subroutine: mld_daggrmat_biz_asb +! Subroutine: mld_daggrmat_biz_bld ! Version: real ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -57,7 +57,7 @@ ! specified by the user through mld_dprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! mld_d_lev_aggrmat_asb. +! mld_d_lev_aggrmat_bld. ! ! Arguments: ! a - type(psb_dspmat_type), input. @@ -80,10 +80,10 @@ ! info - integer, output. ! Error code. ! -subroutine mld_daggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_daggrmat_biz_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_d_inner_mod, mld_protect_name => mld_daggrmat_biz_asb + use mld_d_inner_mod, mld_protect_name => mld_daggrmat_biz_bld implicit none @@ -111,7 +111,7 @@ subroutine mld_daggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr integer(psb_ipk_), parameter :: ncmax=16 real(psb_dpk_) :: anorm, omega, tmp, dg, theta - name='mld_aggrmat_biz_asb' + name='mld_aggrmat_biz_bld' info=psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -372,4 +372,4 @@ subroutine mld_daggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr call psb_error_handler(err_act) return -end subroutine mld_daggrmat_biz_asb +end subroutine mld_daggrmat_biz_bld diff --git a/mlprec/impl/aggregator/mld_daggrmat_minnrg_asb.f90 b/mlprec/impl/aggregator/mld_daggrmat_minnrg_bld.f90 similarity index 98% rename from mlprec/impl/aggregator/mld_daggrmat_minnrg_asb.f90 rename to mlprec/impl/aggregator/mld_daggrmat_minnrg_bld.f90 index e6bf12af..dac54e8f 100644 --- a/mlprec/impl/aggregator/mld_daggrmat_minnrg_asb.f90 +++ b/mlprec/impl/aggregator/mld_daggrmat_minnrg_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_daggrmat_minnrg_asb.F90 +! File: mld_daggrmat_minnrg_bld.F90 ! -! Subroutine: mld_daggrmat_minnrg_asb +! Subroutine: mld_daggrmat_minnrg_bld ! Version: real ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -65,7 +65,7 @@ ! ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! ! Arguments: @@ -104,10 +104,10 @@ ! Error code. ! ! -subroutine mld_daggrmat_minnrg_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_daggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_d_inner_mod, mld_protect_name => mld_daggrmat_minnrg_asb + use mld_d_inner_mod, mld_protect_name => mld_daggrmat_minnrg_bld implicit none @@ -651,4 +651,4 @@ contains close(20+me) end subroutine local_dump -end subroutine mld_daggrmat_minnrg_asb +end subroutine mld_daggrmat_minnrg_bld diff --git a/mlprec/impl/aggregator/mld_daggrmat_nosmth_asb.f90 b/mlprec/impl/aggregator/mld_daggrmat_nosmth_bld.f90 similarity index 96% rename from mlprec/impl/aggregator/mld_daggrmat_nosmth_asb.f90 rename to mlprec/impl/aggregator/mld_daggrmat_nosmth_bld.f90 index 15c8cb4b..705fd1ca 100644 --- a/mlprec/impl/aggregator/mld_daggrmat_nosmth_asb.f90 +++ b/mlprec/impl/aggregator/mld_daggrmat_nosmth_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_daggrmat_nosmth_asb.F90 +! File: mld_daggrmat_nosmth_bld.F90 ! -! Subroutine: mld_daggrmat_nosmth_asb +! Subroutine: mld_daggrmat_nosmth_bld ! Version: real ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -53,7 +53,7 @@ ! specified by the user through mld_dprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! For details see ! P. D'Ambra, D. di Serafino and S. Filippone, On the development of @@ -96,10 +96,10 @@ ! Error code. ! ! -subroutine mld_daggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_daggrmat_nosmth_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_d_inner_mod, mld_protect_name => mld_daggrmat_nosmth_asb + use mld_d_inner_mod, mld_protect_name => mld_daggrmat_nosmth_bld implicit none @@ -124,7 +124,7 @@ subroutine mld_daggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_re integer(psb_lpk_) :: nrow, nglob, ncol, ntaggr, nzl, ip, & & naggr, nzt, naggrm1, naggrp1, i, k - name = 'mld_aggrmat_nosmth_asb' + name = 'mld_aggrmat_nosmth_bld' info = psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -198,4 +198,4 @@ subroutine mld_daggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_re return -end subroutine mld_daggrmat_nosmth_asb +end subroutine mld_daggrmat_nosmth_bld diff --git a/mlprec/impl/aggregator/mld_daggrmat_smth_asb.f90 b/mlprec/impl/aggregator/mld_daggrmat_smth_bld.f90 similarity index 97% rename from mlprec/impl/aggregator/mld_daggrmat_smth_asb.f90 rename to mlprec/impl/aggregator/mld_daggrmat_smth_bld.f90 index e6d0c58c..381dc03e 100644 --- a/mlprec/impl/aggregator/mld_daggrmat_smth_asb.f90 +++ b/mlprec/impl/aggregator/mld_daggrmat_smth_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_daggrmat_smth_asb.F90 +! File: mld_daggrmat_smth_bld.F90 ! -! Subroutine: mld_daggrmat_smth_asb +! Subroutine: mld_daggrmat_smth_bld ! Version: real ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -65,7 +65,7 @@ ! specified by the user through mld_dprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! ! Arguments: @@ -102,10 +102,10 @@ ! info - integer, output. ! Error code. ! -subroutine mld_daggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_daggrmat_smth_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_d_inner_mod, mld_protect_name => mld_daggrmat_smth_asb + use mld_d_inner_mod, mld_protect_name => mld_daggrmat_smth_bld implicit none @@ -133,7 +133,7 @@ subroutine mld_daggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_rest integer(psb_ipk_), parameter :: ncmax=16 real(psb_dpk_) :: anorm, omega, tmp, dg, theta - name='mld_aggrmat_smth_asb' + name='mld_aggrmat_smth_bld' info=psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -411,4 +411,4 @@ subroutine mld_daggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_rest call psb_error_handler(err_act) return -end subroutine mld_daggrmat_smth_asb +end subroutine mld_daggrmat_smth_bld diff --git a/mlprec/impl/aggregator/mld_s_dec_aggregator_mat_asb.f90 b/mlprec/impl/aggregator/mld_s_dec_aggregator_mat_asb.f90 index 38160c08..cee1207f 100644 --- a/mlprec/impl/aggregator/mld_s_dec_aggregator_mat_asb.f90 +++ b/mlprec/impl/aggregator/mld_s_dec_aggregator_mat_asb.f90 @@ -40,63 +40,8 @@ ! Subroutine: mld_s_dec_aggregator_mat_asb ! Version: real ! -! This routine builds the matrix associated to the current level of the -! multilevel preconditioner from the matrix associated to the previous level, -! by using the user-specified aggregation technique (therefore, it also builds the -! prolongation and restriction operators mapping the current level to the -! previous one and vice versa). -! The current level is regarded as the coarse one, while the previous as -! the fine one. This is in agreement with the fact that the routine is called, -! by mld_mlprec_bld, only on levels >=2. -! The coarse-level matrix A_C is built from a fine-level matrix A -! by using the Galerkin approach, i.e. -! -! A_C = P_C^T A P_C, -! -! where P_C is a prolongator from the coarse level to the fine one. -! -! A mapping from the nodes of the adjacency graph of A to the nodes of the -! adjacency graph of A_C has been computed by the mld_aggrmap_bld subroutine. -! The prolongator P_C is built here from this mapping, according to the -! value of p%iprcparm(mld_aggr_kind_), specified by the user through -! mld_sprecinit and mld_zprecset. -! On output from this routine the entries of AC, op_prol, op_restr -! are still in "global numbering" mode; this is fixed in the calling routine -! mld_s_lev_aggrmat_asb. -! -! Currently four different prolongators are implemented, corresponding to -! four aggregation algorithms: -! 1. un-smoothed aggregation, -! 2. smoothed aggregation, -! 3. "bizarre" aggregation. -! 4. minimum energy -! 1. The non-smoothed aggregation uses as prolongator the piecewise constant -! interpolation operator corresponding to the fine-to-coarse level mapping built -! by p%aggr%bld_tprol. This is called tentative prolongator. -! 2. The smoothed aggregation uses as prolongator the operator obtained by applying -! a damped Jacobi smoother to the tentative prolongator. -! 3. The "bizarre" aggregation uses a prolongator proposed by the authors of MLD2P4. -! This prolongator still requires a deep analysis and testing and its use is -! not recommended. -! 4. Minimum energy aggregation -! -! For more details see -! M. Brezina and P. Vanek, A black-box iterative solver based on a two-level -! Schwarz method, Computing, 63 (1999), 233-263. -! P. Vanek, J. Mandel and M. Brezina, Algebraic Multigrid by Smoothed -! Aggregation for Second and Fourth Order Elliptic Problems, Computing, 56 -! (1996), 179-196. -! P. D'Ambra, D. di Serafino and S. Filippone, On the development of PSBLAS-based -! parallel two-level Schwarz preconditioners, Appl. Num. Math., 57 (2007), -! 1181-1196. -! M. Sala, R. Tuminaro: A new Petrov-Galerkin smoothed aggregation preconditioner -! for nonsymmetric linear systems, SIAM J. Sci. Comput., 31(1):143-166 (2008) -! -! -! The main structure is: -! 1. Perform sanity checks; -! 2. Compute prolongator/restrictor/AC ! +! From a given AC to final format, generating DESC_AC ! ! Arguments: ! ag - type(mld_s_dec_aggregator_type), input/output. @@ -121,89 +66,143 @@ ! the various processes do not overlap. ! nlaggr - integer, dimension(:) input ! nlaggr(i) contains the aggregates held by process i. -! ac - type(psb_sspmat_type), output -! The coarse matrix on output +! ac - type(psb_sspmat_type), inout +! The coarse matrix +! desc_ac - type(psb_desc_type), output. +! The communication descriptor of the fine-level matrix. +! The 'one-level' data structure that will contain the local +! part of the matrix to be built as well as the information +! concerning the prolongator and its transpose. ! ! op_prol - type(psb_sspmat_type), input/output ! The tentative prolongator on input, the computed prolongator on output ! -! op_restr - type(psb_sspmat_type), output +! op_restr - type(psb_sspmat_type), input/output ! The restrictor operator; normally, it is the transpose of the prolongator. ! ! info - integer, output. ! Error code. ! -subroutine mld_s_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_prol,op_restr,info) +subroutine mld_s_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac, op_prol,op_restr,info) use psb_base_mod - use mld_s_prec_type, mld_protect_name => mld_s_dec_aggregator_mat_asb - use mld_s_inner_mod + use mld_base_prec_type + use mld_s_dec_aggregator_mod, mld_protect_name => mld_s_dec_aggregator_mat_asb implicit none - class(mld_s_dec_aggregator_type), target, intent(inout) :: ag - type(mld_sml_parms), intent(inout) :: parms + type(mld_sml_parms), intent(inout) :: parms type(psb_sspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) - type(psb_lsspmat_type), intent(inout) :: op_prol - type(psb_lsspmat_type), intent(out) :: ac,op_restr - integer(psb_ipk_), intent(out) :: info - - ! Local variables - character(len=20) :: name - integer(psb_mpk_) :: ictxt, np, me - type(psb_ls_coo_sparse_mat) :: acoo, bcoo - type(psb_ls_csr_sparse_mat) :: acsr1 - integer(psb_lpk_) :: nzl,ntaggr - integer(psb_ipk_) :: err_act - integer(psb_ipk_) :: debug_level, debug_unit - - name='mld_s_dec_aggregator_mat_asb' + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lsspmat_type), intent(inout) :: op_prol, ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + ! + integer(psb_mpk_) :: ictxt, np, me + type(psb_ls_coo_sparse_mat) :: acoo, bcoo + type(psb_ls_csr_sparse_mat) :: acsr1 + type(psb_lsspmat_type) :: lac, lac1 + type(psb_sspmat_type) :: tmp_ac + integer(psb_ipk_) :: i_nr, i_nc, i_nl, nzl + integer(psb_lpk_) :: ntaggr, l_nr, l_nc + integer(psb_ipk_) :: err_act, debug_level, debug_unit + character(len=20) :: name='s_dec_aggregator_mat_asb' + + + if (psb_get_errstatus().ne.0) return call psb_erractionsave(err_act) - if (psb_errstatus_fatal()) then - info = psb_err_internal_error_; goto 9999 - end if debug_unit = psb_get_debug_unit() debug_level = psb_get_debug_level() info = psb_success_ ictxt = desc_a%get_context() call psb_info(ictxt,me,np) - ! - ! Build the coarse-level matrix from the fine-level one, starting from - ! the mapping defined by mld_aggrmap_bld and applying the aggregation - ! algorithm specified by - ! - select case (parms%aggr_prol) - case (mld_no_smooth_) - - call mld_saggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,& - & parms,ac,op_prol,op_restr,info) - - case(mld_smooth_prol_) - - call mld_saggrmat_smth_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - - case(mld_biz_prol_) - call mld_saggrmat_biz_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - case(mld_min_energy_) - - call mld_saggrmat_minnrg_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - - case default + ntaggr = sum(nlaggr) + + select case(parms%coarse_mat) + + case(mld_distr_mat_) + + call ac%mv_to(bcoo) + nzl = bcoo%get_nzeros() + i_nl = nlaggr(me+1) + if (info == psb_success_) call psb_cdall(ictxt,desc_ac,info,nl=i_nl) + if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,desc_ac,info) + if (info == psb_success_) call psb_cdasb(desc_ac,info) + if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),desc_ac,info,iact='I') + if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),desc_ac,info,iact='I') + if (info /= psb_success_) then + call psb_errpush(psb_err_internal_error_,name,& + & a_err='Creating desc_ac and converting ac') + goto 9999 + end if + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),& + & 'Assembld aux descr. distr.' + call ac%mv_from(bcoo) + l_nr = desc_ac%get_local_rows() + l_nc = desc_ac%get_local_cols() + call ac%set_nrows(l_nr) + call ac%set_ncols(l_nc) + call ac%set_asb() + + if (info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free') + goto 9999 + end if + + if (np>1) then + call op_prol%mv_to(acsr1) + nzl = acsr1%get_nzeros() + call psb_glob_to_loc(acsr1%ja(1:nzl),desc_ac,info,'I') + if(info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc') + goto 9999 + end if + call op_prol%mv_from(acsr1) + endif + call op_prol%set_ncols(l_nc) + + if (np>1) then + call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_) + call op_restr%mv_to(acoo) + nzl = acoo%get_nzeros() + if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),desc_ac,info,'I') + call acoo%set_dupl(psb_dupl_add_) + if (info == psb_success_) call op_restr%mv_from(acoo) + if (info == psb_success_) call op_restr%cscnv(info,type='csr') + if(info /= psb_success_) then + call psb_errpush(psb_err_internal_error_,name,& + & a_err='Converting op_restr to local') + goto 9999 + end if + end if + ! + ! Clip to local rows. + ! + call op_restr%set_nrows(l_nr) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),& + & 'Done ac ' + + case(mld_repl_mat_) + ! + ! + call psb_cdall(ictxt,desc_ac,info,mg=ntaggr,repl=.true.) + if (info == psb_success_) call psb_cdasb(desc_ac,info) + if (info == psb_success_) call tmp_ac%mv_from_l(ac) + if (info == psb_success_) & + & call psb_gather(ac,tmp_ac,desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) + if (info /= psb_success_) goto 9999 + + case default info = psb_err_internal_error_ - call psb_errpush(info,name,a_err='Invalid aggr kind') + call psb_errpush(info,name,a_err='invalid mld_coarse_mat_') goto 9999 - end select - if (info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='Inner aggrmat asb') - goto 9999 - end if call psb_erractionrestore(err_act) @@ -212,5 +211,5 @@ subroutine mld_s_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_p 9999 call psb_error_handler(err_act) return - + end subroutine mld_s_dec_aggregator_mat_asb diff --git a/mlprec/impl/aggregator/mld_s_dec_aggregator_mat_bld.f90 b/mlprec/impl/aggregator/mld_s_dec_aggregator_mat_bld.f90 new file mode 100644 index 00000000..e1008f91 --- /dev/null +++ b/mlprec/impl/aggregator/mld_s_dec_aggregator_mat_bld.f90 @@ -0,0 +1,216 @@ +! +! +! MLD2P4 version 2.2 +! MultiLevel Domain Decomposition Parallel Preconditioners Package +! based on PSBLAS (Parallel Sparse BLAS version 3.5) +! +! (C) Copyright 2008-2018 +! +! Salvatore Filippone +! Pasqua D'Ambra +! Daniela di Serafino +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the MLD2P4 group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE MLD2P4 GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: mld_s_dec_aggregator_mat_bld.f90 +! +! Subroutine: mld_s_dec_aggregator_mat_bld +! Version: real +! +! This routine builds the matrix associated to the current level of the +! multilevel preconditioner from the matrix associated to the previous level, +! by using the user-specified aggregation technique (therefore, it also builds the +! prolongation and restriction operators mapping the current level to the +! previous one and vice versa). +! The current level is regarded as the coarse one, while the previous as +! the fine one. This is in agreement with the fact that the routine is called, +! by mld_mlprec_bld, only on levels >=2. +! The coarse-level matrix A_C is built from a fine-level matrix A +! by using the Galerkin approach, i.e. +! +! A_C = P_C^T A P_C, +! +! where P_C is a prolongator from the coarse level to the fine one. +! +! A mapping from the nodes of the adjacency graph of A to the nodes of the +! adjacency graph of A_C has been computed by the mld_aggrmap_bld subroutine. +! The prolongator P_C is built here from this mapping, according to the +! value of p%iprcparm(mld_aggr_kind_), specified by the user through +! mld_sprecinit and mld_zprecset. +! On output from this routine the entries of AC, op_prol, op_restr +! are still in "global numbering" mode; this is fixed in the calling routine +! mld_s_lev_aggrmat_bld. +! +! Currently four different prolongators are implemented, corresponding to +! four aggregation algorithms: +! 1. un-smoothed aggregation, +! 2. smoothed aggregation, +! 3. "bizarre" aggregation. +! 4. minimum energy +! 1. The non-smoothed aggregation uses as prolongator the piecewise constant +! interpolation operator corresponding to the fine-to-coarse level mapping built +! by p%aggr%bld_tprol. This is called tentative prolongator. +! 2. The smoothed aggregation uses as prolongator the operator obtained by applying +! a damped Jacobi smoother to the tentative prolongator. +! 3. The "bizarre" aggregation uses a prolongator proposed by the authors of MLD2P4. +! This prolongator still requires a deep analysis and testing and its use is +! not recommended. +! 4. Minimum energy aggregation +! +! For more details see +! M. Brezina and P. Vanek, A black-box iterative solver based on a two-level +! Schwarz method, Computing, 63 (1999), 233-263. +! P. Vanek, J. Mandel and M. Brezina, Algebraic Multigrid by Smoothed +! Aggregation for Second and Fourth Order Elliptic Problems, Computing, 56 +! (1996), 179-196. +! P. D'Ambra, D. di Serafino and S. Filippone, On the development of PSBLAS-based +! parallel two-level Schwarz preconditioners, Appl. Num. Math., 57 (2007), +! 1181-1196. +! M. Sala, R. Tuminaro: A new Petrov-Galerkin smoothed aggregation preconditioner +! for nonsymmetric linear systems, SIAM J. Sci. Comput., 31(1):143-166 (2008) +! +! +! The main structure is: +! 1. Perform sanity checks; +! 2. Compute prolongator/restrictor/AC +! +! +! Arguments: +! ag - type(mld_s_dec_aggregator_type), input/output. +! The aggregator object +! parms - type(mld_sml_parms), input +! The aggregation parameters +! a - type(psb_sspmat_type), input. +! The sparse matrix structure containing the local part of +! the fine-level matrix. +! desc_a - type(psb_desc_type), input. +! The communication descriptor of the fine-level matrix. +! The 'one-level' data structure that will contain the local +! part of the matrix to be built as well as the information +! concerning the prolongator and its transpose. +! ilaggr - integer, dimension(:), input +! The mapping between the row indices of the coarse-level +! matrix and the row indices of the fine-level matrix. +! ilaggr(i)=j means that node i in the adjacency graph +! of the fine-level matrix is mapped onto node j in the +! adjacency graph of the coarse-level matrix. Note that the indices +! are assumed to be shifted so as to make sure the ranges on +! the various processes do not overlap. +! nlaggr - integer, dimension(:) input +! nlaggr(i) contains the aggregates held by process i. +! ac - type(psb_sspmat_type), output +! The coarse matrix on output +! +! op_prol - type(psb_sspmat_type), input/output +! The tentative prolongator on input, the computed prolongator on output +! +! op_restr - type(psb_sspmat_type), output +! The restrictor operator; normally, it is the transpose of the prolongator. +! +! info - integer, output. +! Error code. +! +subroutine mld_s_dec_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_prol,op_restr,info) + use psb_base_mod + use mld_s_prec_type, mld_protect_name => mld_s_dec_aggregator_mat_bld + use mld_s_inner_mod + implicit none + + class(mld_s_dec_aggregator_type), target, intent(inout) :: ag + type(mld_sml_parms), intent(inout) :: parms + type(psb_sspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lsspmat_type), intent(inout) :: op_prol + type(psb_lsspmat_type), intent(out) :: ac,op_restr + integer(psb_ipk_), intent(out) :: info + + ! Local variables + character(len=20) :: name + integer(psb_mpk_) :: ictxt, np, me + type(psb_ls_coo_sparse_mat) :: acoo, bcoo + type(psb_ls_csr_sparse_mat) :: acsr1 + integer(psb_lpk_) :: nzl,ntaggr + integer(psb_ipk_) :: err_act + integer(psb_ipk_) :: debug_level, debug_unit + + name='mld_s_dec_aggregator_mat_bld' + call psb_erractionsave(err_act) + if (psb_errstatus_fatal()) then + info = psb_err_internal_error_; goto 9999 + end if + debug_unit = psb_get_debug_unit() + debug_level = psb_get_debug_level() + info = psb_success_ + ictxt = desc_a%get_context() + call psb_info(ictxt,me,np) + + ! + ! Build the coarse-level matrix from the fine-level one, starting from + ! the mapping defined by mld_aggrmap_bld and applying the aggregation + ! algorithm specified by + ! + select case (parms%aggr_prol) + case (mld_no_smooth_) + + call mld_saggrmat_nosmth_bld(a,desc_a,ilaggr,nlaggr,& + & parms,ac,op_prol,op_restr,info) + + case(mld_smooth_prol_) + + call mld_saggrmat_smth_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case(mld_biz_prol_) + + call mld_saggrmat_biz_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case(mld_min_energy_) + + call mld_saggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case default + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Invalid aggr kind') + goto 9999 + + end select + if (info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='Inner aggrmat bld') + goto 9999 + end if + + + call psb_erractionrestore(err_act) + return + +9999 call psb_error_handler(err_act) + return + + +end subroutine mld_s_dec_aggregator_mat_bld diff --git a/mlprec/impl/aggregator/mld_saggrmat_biz_asb.f90 b/mlprec/impl/aggregator/mld_saggrmat_biz_bld.f90 similarity index 97% rename from mlprec/impl/aggregator/mld_saggrmat_biz_asb.f90 rename to mlprec/impl/aggregator/mld_saggrmat_biz_bld.f90 index d7027abf..de50f467 100644 --- a/mlprec/impl/aggregator/mld_saggrmat_biz_asb.f90 +++ b/mlprec/impl/aggregator/mld_saggrmat_biz_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_saggrmat_biz_asb.F90 +! File: mld_saggrmat_biz_bld.F90 ! -! Subroutine: mld_saggrmat_biz_asb +! Subroutine: mld_saggrmat_biz_bld ! Version: real ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -57,7 +57,7 @@ ! specified by the user through mld_sprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! mld_s_lev_aggrmat_asb. +! mld_s_lev_aggrmat_bld. ! ! Arguments: ! a - type(psb_sspmat_type), input. @@ -80,10 +80,10 @@ ! info - integer, output. ! Error code. ! -subroutine mld_saggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_saggrmat_biz_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_s_inner_mod, mld_protect_name => mld_saggrmat_biz_asb + use mld_s_inner_mod, mld_protect_name => mld_saggrmat_biz_bld implicit none @@ -111,7 +111,7 @@ subroutine mld_saggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr integer(psb_ipk_), parameter :: ncmax=16 real(psb_spk_) :: anorm, omega, tmp, dg, theta - name='mld_aggrmat_biz_asb' + name='mld_aggrmat_biz_bld' info=psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -372,4 +372,4 @@ subroutine mld_saggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr call psb_error_handler(err_act) return -end subroutine mld_saggrmat_biz_asb +end subroutine mld_saggrmat_biz_bld diff --git a/mlprec/impl/aggregator/mld_saggrmat_minnrg_asb.f90 b/mlprec/impl/aggregator/mld_saggrmat_minnrg_bld.f90 similarity index 98% rename from mlprec/impl/aggregator/mld_saggrmat_minnrg_asb.f90 rename to mlprec/impl/aggregator/mld_saggrmat_minnrg_bld.f90 index 56da4362..4f765e85 100644 --- a/mlprec/impl/aggregator/mld_saggrmat_minnrg_asb.f90 +++ b/mlprec/impl/aggregator/mld_saggrmat_minnrg_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_saggrmat_minnrg_asb.F90 +! File: mld_saggrmat_minnrg_bld.F90 ! -! Subroutine: mld_saggrmat_minnrg_asb +! Subroutine: mld_saggrmat_minnrg_bld ! Version: real ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -65,7 +65,7 @@ ! ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! ! Arguments: @@ -104,10 +104,10 @@ ! Error code. ! ! -subroutine mld_saggrmat_minnrg_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_saggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_s_inner_mod, mld_protect_name => mld_saggrmat_minnrg_asb + use mld_s_inner_mod, mld_protect_name => mld_saggrmat_minnrg_bld implicit none @@ -651,4 +651,4 @@ contains close(20+me) end subroutine local_dump -end subroutine mld_saggrmat_minnrg_asb +end subroutine mld_saggrmat_minnrg_bld diff --git a/mlprec/impl/aggregator/mld_saggrmat_nosmth_asb.f90 b/mlprec/impl/aggregator/mld_saggrmat_nosmth_bld.f90 similarity index 96% rename from mlprec/impl/aggregator/mld_saggrmat_nosmth_asb.f90 rename to mlprec/impl/aggregator/mld_saggrmat_nosmth_bld.f90 index 9b8d99bf..c4c3ab8f 100644 --- a/mlprec/impl/aggregator/mld_saggrmat_nosmth_asb.f90 +++ b/mlprec/impl/aggregator/mld_saggrmat_nosmth_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_saggrmat_nosmth_asb.F90 +! File: mld_saggrmat_nosmth_bld.F90 ! -! Subroutine: mld_saggrmat_nosmth_asb +! Subroutine: mld_saggrmat_nosmth_bld ! Version: real ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -53,7 +53,7 @@ ! specified by the user through mld_sprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! For details see ! P. D'Ambra, D. di Serafino and S. Filippone, On the development of @@ -96,10 +96,10 @@ ! Error code. ! ! -subroutine mld_saggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_saggrmat_nosmth_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_s_inner_mod, mld_protect_name => mld_saggrmat_nosmth_asb + use mld_s_inner_mod, mld_protect_name => mld_saggrmat_nosmth_bld implicit none @@ -124,7 +124,7 @@ subroutine mld_saggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_re integer(psb_lpk_) :: nrow, nglob, ncol, ntaggr, nzl, ip, & & naggr, nzt, naggrm1, naggrp1, i, k - name = 'mld_aggrmat_nosmth_asb' + name = 'mld_aggrmat_nosmth_bld' info = psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -198,4 +198,4 @@ subroutine mld_saggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_re return -end subroutine mld_saggrmat_nosmth_asb +end subroutine mld_saggrmat_nosmth_bld diff --git a/mlprec/impl/aggregator/mld_saggrmat_smth_asb.f90 b/mlprec/impl/aggregator/mld_saggrmat_smth_bld.f90 similarity index 97% rename from mlprec/impl/aggregator/mld_saggrmat_smth_asb.f90 rename to mlprec/impl/aggregator/mld_saggrmat_smth_bld.f90 index 2557076e..fcf159d3 100644 --- a/mlprec/impl/aggregator/mld_saggrmat_smth_asb.f90 +++ b/mlprec/impl/aggregator/mld_saggrmat_smth_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_saggrmat_smth_asb.F90 +! File: mld_saggrmat_smth_bld.F90 ! -! Subroutine: mld_saggrmat_smth_asb +! Subroutine: mld_saggrmat_smth_bld ! Version: real ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -65,7 +65,7 @@ ! specified by the user through mld_sprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! ! Arguments: @@ -102,10 +102,10 @@ ! info - integer, output. ! Error code. ! -subroutine mld_saggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_saggrmat_smth_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_s_inner_mod, mld_protect_name => mld_saggrmat_smth_asb + use mld_s_inner_mod, mld_protect_name => mld_saggrmat_smth_bld implicit none @@ -133,7 +133,7 @@ subroutine mld_saggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_rest integer(psb_ipk_), parameter :: ncmax=16 real(psb_spk_) :: anorm, omega, tmp, dg, theta - name='mld_aggrmat_smth_asb' + name='mld_aggrmat_smth_bld' info=psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -411,4 +411,4 @@ subroutine mld_saggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_rest call psb_error_handler(err_act) return -end subroutine mld_saggrmat_smth_asb +end subroutine mld_saggrmat_smth_bld diff --git a/mlprec/impl/aggregator/mld_z_dec_aggregator_mat_asb.f90 b/mlprec/impl/aggregator/mld_z_dec_aggregator_mat_asb.f90 index a709179d..bdbf6cd2 100644 --- a/mlprec/impl/aggregator/mld_z_dec_aggregator_mat_asb.f90 +++ b/mlprec/impl/aggregator/mld_z_dec_aggregator_mat_asb.f90 @@ -40,63 +40,8 @@ ! Subroutine: mld_z_dec_aggregator_mat_asb ! Version: complex ! -! This routine builds the matrix associated to the current level of the -! multilevel preconditioner from the matrix associated to the previous level, -! by using the user-specified aggregation technique (therefore, it also builds the -! prolongation and restriction operators mapping the current level to the -! previous one and vice versa). -! The current level is regarded as the coarse one, while the previous as -! the fine one. This is in agreement with the fact that the routine is called, -! by mld_mlprec_bld, only on levels >=2. -! The coarse-level matrix A_C is built from a fine-level matrix A -! by using the Galerkin approach, i.e. -! -! A_C = P_C^T A P_C, -! -! where P_C is a prolongator from the coarse level to the fine one. -! -! A mapping from the nodes of the adjacency graph of A to the nodes of the -! adjacency graph of A_C has been computed by the mld_aggrmap_bld subroutine. -! The prolongator P_C is built here from this mapping, according to the -! value of p%iprcparm(mld_aggr_kind_), specified by the user through -! mld_zprecinit and mld_zprecset. -! On output from this routine the entries of AC, op_prol, op_restr -! are still in "global numbering" mode; this is fixed in the calling routine -! mld_z_lev_aggrmat_asb. -! -! Currently four different prolongators are implemented, corresponding to -! four aggregation algorithms: -! 1. un-smoothed aggregation, -! 2. smoothed aggregation, -! 3. "bizarre" aggregation. -! 4. minimum energy -! 1. The non-smoothed aggregation uses as prolongator the piecewise constant -! interpolation operator corresponding to the fine-to-coarse level mapping built -! by p%aggr%bld_tprol. This is called tentative prolongator. -! 2. The smoothed aggregation uses as prolongator the operator obtained by applying -! a damped Jacobi smoother to the tentative prolongator. -! 3. The "bizarre" aggregation uses a prolongator proposed by the authors of MLD2P4. -! This prolongator still requires a deep analysis and testing and its use is -! not recommended. -! 4. Minimum energy aggregation -! -! For more details see -! M. Brezina and P. Vanek, A black-box iterative solver based on a two-level -! Schwarz method, Computing, 63 (1999), 233-263. -! P. Vanek, J. Mandel and M. Brezina, Algebraic Multigrid by Smoothed -! Aggregation for Second and Fourth Order Elliptic Problems, Computing, 56 -! (1996), 179-196. -! P. D'Ambra, D. di Serafino and S. Filippone, On the development of PSBLAS-based -! parallel two-level Schwarz preconditioners, Appl. Num. Math., 57 (2007), -! 1181-1196. -! M. Sala, R. Tuminaro: A new Petrov-Galerkin smoothed aggregation preconditioner -! for nonsymmetric linear systems, SIAM J. Sci. Comput., 31(1):143-166 (2008) -! -! -! The main structure is: -! 1. Perform sanity checks; -! 2. Compute prolongator/restrictor/AC ! +! From a given AC to final format, generating DESC_AC ! ! Arguments: ! ag - type(mld_z_dec_aggregator_type), input/output. @@ -121,89 +66,143 @@ ! the various processes do not overlap. ! nlaggr - integer, dimension(:) input ! nlaggr(i) contains the aggregates held by process i. -! ac - type(psb_zspmat_type), output -! The coarse matrix on output +! ac - type(psb_zspmat_type), inout +! The coarse matrix +! desc_ac - type(psb_desc_type), output. +! The communication descriptor of the fine-level matrix. +! The 'one-level' data structure that will contain the local +! part of the matrix to be built as well as the information +! concerning the prolongator and its transpose. ! ! op_prol - type(psb_zspmat_type), input/output ! The tentative prolongator on input, the computed prolongator on output ! -! op_restr - type(psb_zspmat_type), output +! op_restr - type(psb_zspmat_type), input/output ! The restrictor operator; normally, it is the transpose of the prolongator. ! ! info - integer, output. ! Error code. ! -subroutine mld_z_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_prol,op_restr,info) +subroutine mld_z_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac, op_prol,op_restr,info) use psb_base_mod - use mld_z_prec_type, mld_protect_name => mld_z_dec_aggregator_mat_asb - use mld_z_inner_mod + use mld_base_prec_type + use mld_z_dec_aggregator_mod, mld_protect_name => mld_z_dec_aggregator_mat_asb implicit none - class(mld_z_dec_aggregator_type), target, intent(inout) :: ag - type(mld_dml_parms), intent(inout) :: parms + type(mld_dml_parms), intent(inout) :: parms type(psb_zspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) - type(psb_lzspmat_type), intent(inout) :: op_prol - type(psb_lzspmat_type), intent(out) :: ac,op_restr - integer(psb_ipk_), intent(out) :: info - - ! Local variables - character(len=20) :: name - integer(psb_mpk_) :: ictxt, np, me - type(psb_lz_coo_sparse_mat) :: acoo, bcoo - type(psb_lz_csr_sparse_mat) :: acsr1 - integer(psb_lpk_) :: nzl,ntaggr - integer(psb_ipk_) :: err_act - integer(psb_ipk_) :: debug_level, debug_unit - - name='mld_z_dec_aggregator_mat_asb' + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lzspmat_type), intent(inout) :: op_prol, ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + ! + integer(psb_mpk_) :: ictxt, np, me + type(psb_lz_coo_sparse_mat) :: acoo, bcoo + type(psb_lz_csr_sparse_mat) :: acsr1 + type(psb_lzspmat_type) :: lac, lac1 + type(psb_zspmat_type) :: tmp_ac + integer(psb_ipk_) :: i_nr, i_nc, i_nl, nzl + integer(psb_lpk_) :: ntaggr, l_nr, l_nc + integer(psb_ipk_) :: err_act, debug_level, debug_unit + character(len=20) :: name='z_dec_aggregator_mat_asb' + + + if (psb_get_errstatus().ne.0) return call psb_erractionsave(err_act) - if (psb_errstatus_fatal()) then - info = psb_err_internal_error_; goto 9999 - end if debug_unit = psb_get_debug_unit() debug_level = psb_get_debug_level() info = psb_success_ ictxt = desc_a%get_context() call psb_info(ictxt,me,np) - ! - ! Build the coarse-level matrix from the fine-level one, starting from - ! the mapping defined by mld_aggrmap_bld and applying the aggregation - ! algorithm specified by - ! - select case (parms%aggr_prol) - case (mld_no_smooth_) - - call mld_zaggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,& - & parms,ac,op_prol,op_restr,info) - - case(mld_smooth_prol_) - - call mld_zaggrmat_smth_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - - case(mld_biz_prol_) - call mld_zaggrmat_biz_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - case(mld_min_energy_) - - call mld_zaggrmat_minnrg_asb(a,desc_a,ilaggr,nlaggr, & - & parms,ac,op_prol,op_restr,info) - - case default + ntaggr = sum(nlaggr) + + select case(parms%coarse_mat) + + case(mld_distr_mat_) + + call ac%mv_to(bcoo) + nzl = bcoo%get_nzeros() + i_nl = nlaggr(me+1) + if (info == psb_success_) call psb_cdall(ictxt,desc_ac,info,nl=i_nl) + if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,desc_ac,info) + if (info == psb_success_) call psb_cdasb(desc_ac,info) + if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),desc_ac,info,iact='I') + if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),desc_ac,info,iact='I') + if (info /= psb_success_) then + call psb_errpush(psb_err_internal_error_,name,& + & a_err='Creating desc_ac and converting ac') + goto 9999 + end if + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),& + & 'Assembld aux descr. distr.' + call ac%mv_from(bcoo) + l_nr = desc_ac%get_local_rows() + l_nc = desc_ac%get_local_cols() + call ac%set_nrows(l_nr) + call ac%set_ncols(l_nc) + call ac%set_asb() + + if (info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free') + goto 9999 + end if + + if (np>1) then + call op_prol%mv_to(acsr1) + nzl = acsr1%get_nzeros() + call psb_glob_to_loc(acsr1%ja(1:nzl),desc_ac,info,'I') + if(info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc') + goto 9999 + end if + call op_prol%mv_from(acsr1) + endif + call op_prol%set_ncols(l_nc) + + if (np>1) then + call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_) + call op_restr%mv_to(acoo) + nzl = acoo%get_nzeros() + if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),desc_ac,info,'I') + call acoo%set_dupl(psb_dupl_add_) + if (info == psb_success_) call op_restr%mv_from(acoo) + if (info == psb_success_) call op_restr%cscnv(info,type='csr') + if(info /= psb_success_) then + call psb_errpush(psb_err_internal_error_,name,& + & a_err='Converting op_restr to local') + goto 9999 + end if + end if + ! + ! Clip to local rows. + ! + call op_restr%set_nrows(l_nr) + + if (debug_level >= psb_debug_outer_) & + & write(debug_unit,*) me,' ',trim(name),& + & 'Done ac ' + + case(mld_repl_mat_) + ! + ! + call psb_cdall(ictxt,desc_ac,info,mg=ntaggr,repl=.true.) + if (info == psb_success_) call psb_cdasb(desc_ac,info) + if (info == psb_success_) call tmp_ac%mv_from_l(ac) + if (info == psb_success_) & + & call psb_gather(ac,tmp_ac,desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) + if (info /= psb_success_) goto 9999 + + case default info = psb_err_internal_error_ - call psb_errpush(info,name,a_err='Invalid aggr kind') + call psb_errpush(info,name,a_err='invalid mld_coarse_mat_') goto 9999 - end select - if (info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='Inner aggrmat asb') - goto 9999 - end if call psb_erractionrestore(err_act) @@ -212,5 +211,5 @@ subroutine mld_z_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_p 9999 call psb_error_handler(err_act) return - + end subroutine mld_z_dec_aggregator_mat_asb diff --git a/mlprec/impl/aggregator/mld_z_dec_aggregator_mat_bld.f90 b/mlprec/impl/aggregator/mld_z_dec_aggregator_mat_bld.f90 new file mode 100644 index 00000000..219049eb --- /dev/null +++ b/mlprec/impl/aggregator/mld_z_dec_aggregator_mat_bld.f90 @@ -0,0 +1,216 @@ +! +! +! MLD2P4 version 2.2 +! MultiLevel Domain Decomposition Parallel Preconditioners Package +! based on PSBLAS (Parallel Sparse BLAS version 3.5) +! +! (C) Copyright 2008-2018 +! +! Salvatore Filippone +! Pasqua D'Ambra +! Daniela di Serafino +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the MLD2P4 group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE MLD2P4 GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: mld_z_dec_aggregator_mat_bld.f90 +! +! Subroutine: mld_z_dec_aggregator_mat_bld +! Version: complex +! +! This routine builds the matrix associated to the current level of the +! multilevel preconditioner from the matrix associated to the previous level, +! by using the user-specified aggregation technique (therefore, it also builds the +! prolongation and restriction operators mapping the current level to the +! previous one and vice versa). +! The current level is regarded as the coarse one, while the previous as +! the fine one. This is in agreement with the fact that the routine is called, +! by mld_mlprec_bld, only on levels >=2. +! The coarse-level matrix A_C is built from a fine-level matrix A +! by using the Galerkin approach, i.e. +! +! A_C = P_C^T A P_C, +! +! where P_C is a prolongator from the coarse level to the fine one. +! +! A mapping from the nodes of the adjacency graph of A to the nodes of the +! adjacency graph of A_C has been computed by the mld_aggrmap_bld subroutine. +! The prolongator P_C is built here from this mapping, according to the +! value of p%iprcparm(mld_aggr_kind_), specified by the user through +! mld_zprecinit and mld_zprecset. +! On output from this routine the entries of AC, op_prol, op_restr +! are still in "global numbering" mode; this is fixed in the calling routine +! mld_z_lev_aggrmat_bld. +! +! Currently four different prolongators are implemented, corresponding to +! four aggregation algorithms: +! 1. un-smoothed aggregation, +! 2. smoothed aggregation, +! 3. "bizarre" aggregation. +! 4. minimum energy +! 1. The non-smoothed aggregation uses as prolongator the piecewise constant +! interpolation operator corresponding to the fine-to-coarse level mapping built +! by p%aggr%bld_tprol. This is called tentative prolongator. +! 2. The smoothed aggregation uses as prolongator the operator obtained by applying +! a damped Jacobi smoother to the tentative prolongator. +! 3. The "bizarre" aggregation uses a prolongator proposed by the authors of MLD2P4. +! This prolongator still requires a deep analysis and testing and its use is +! not recommended. +! 4. Minimum energy aggregation +! +! For more details see +! M. Brezina and P. Vanek, A black-box iterative solver based on a two-level +! Schwarz method, Computing, 63 (1999), 233-263. +! P. Vanek, J. Mandel and M. Brezina, Algebraic Multigrid by Smoothed +! Aggregation for Second and Fourth Order Elliptic Problems, Computing, 56 +! (1996), 179-196. +! P. D'Ambra, D. di Serafino and S. Filippone, On the development of PSBLAS-based +! parallel two-level Schwarz preconditioners, Appl. Num. Math., 57 (2007), +! 1181-1196. +! M. Sala, R. Tuminaro: A new Petrov-Galerkin smoothed aggregation preconditioner +! for nonsymmetric linear systems, SIAM J. Sci. Comput., 31(1):143-166 (2008) +! +! +! The main structure is: +! 1. Perform sanity checks; +! 2. Compute prolongator/restrictor/AC +! +! +! Arguments: +! ag - type(mld_z_dec_aggregator_type), input/output. +! The aggregator object +! parms - type(mld_dml_parms), input +! The aggregation parameters +! a - type(psb_zspmat_type), input. +! The sparse matrix structure containing the local part of +! the fine-level matrix. +! desc_a - type(psb_desc_type), input. +! The communication descriptor of the fine-level matrix. +! The 'one-level' data structure that will contain the local +! part of the matrix to be built as well as the information +! concerning the prolongator and its transpose. +! ilaggr - integer, dimension(:), input +! The mapping between the row indices of the coarse-level +! matrix and the row indices of the fine-level matrix. +! ilaggr(i)=j means that node i in the adjacency graph +! of the fine-level matrix is mapped onto node j in the +! adjacency graph of the coarse-level matrix. Note that the indices +! are assumed to be shifted so as to make sure the ranges on +! the various processes do not overlap. +! nlaggr - integer, dimension(:) input +! nlaggr(i) contains the aggregates held by process i. +! ac - type(psb_zspmat_type), output +! The coarse matrix on output +! +! op_prol - type(psb_zspmat_type), input/output +! The tentative prolongator on input, the computed prolongator on output +! +! op_restr - type(psb_zspmat_type), output +! The restrictor operator; normally, it is the transpose of the prolongator. +! +! info - integer, output. +! Error code. +! +subroutine mld_z_dec_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,op_prol,op_restr,info) + use psb_base_mod + use mld_z_prec_type, mld_protect_name => mld_z_dec_aggregator_mat_bld + use mld_z_inner_mod + implicit none + + class(mld_z_dec_aggregator_type), target, intent(inout) :: ag + type(mld_dml_parms), intent(inout) :: parms + type(psb_zspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lzspmat_type), intent(inout) :: op_prol + type(psb_lzspmat_type), intent(out) :: ac,op_restr + integer(psb_ipk_), intent(out) :: info + + ! Local variables + character(len=20) :: name + integer(psb_mpk_) :: ictxt, np, me + type(psb_lz_coo_sparse_mat) :: acoo, bcoo + type(psb_lz_csr_sparse_mat) :: acsr1 + integer(psb_lpk_) :: nzl,ntaggr + integer(psb_ipk_) :: err_act + integer(psb_ipk_) :: debug_level, debug_unit + + name='mld_z_dec_aggregator_mat_bld' + call psb_erractionsave(err_act) + if (psb_errstatus_fatal()) then + info = psb_err_internal_error_; goto 9999 + end if + debug_unit = psb_get_debug_unit() + debug_level = psb_get_debug_level() + info = psb_success_ + ictxt = desc_a%get_context() + call psb_info(ictxt,me,np) + + ! + ! Build the coarse-level matrix from the fine-level one, starting from + ! the mapping defined by mld_aggrmap_bld and applying the aggregation + ! algorithm specified by + ! + select case (parms%aggr_prol) + case (mld_no_smooth_) + + call mld_zaggrmat_nosmth_bld(a,desc_a,ilaggr,nlaggr,& + & parms,ac,op_prol,op_restr,info) + + case(mld_smooth_prol_) + + call mld_zaggrmat_smth_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case(mld_biz_prol_) + + call mld_zaggrmat_biz_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case(mld_min_energy_) + + call mld_zaggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr, & + & parms,ac,op_prol,op_restr,info) + + case default + info = psb_err_internal_error_ + call psb_errpush(info,name,a_err='Invalid aggr kind') + goto 9999 + + end select + if (info /= psb_success_) then + call psb_errpush(psb_err_from_subroutine_,name,a_err='Inner aggrmat bld') + goto 9999 + end if + + + call psb_erractionrestore(err_act) + return + +9999 call psb_error_handler(err_act) + return + + +end subroutine mld_z_dec_aggregator_mat_bld diff --git a/mlprec/impl/aggregator/mld_zaggrmat_biz_asb.f90 b/mlprec/impl/aggregator/mld_zaggrmat_biz_bld.f90 similarity index 97% rename from mlprec/impl/aggregator/mld_zaggrmat_biz_asb.f90 rename to mlprec/impl/aggregator/mld_zaggrmat_biz_bld.f90 index d39b9c26..2d32ecfd 100644 --- a/mlprec/impl/aggregator/mld_zaggrmat_biz_asb.f90 +++ b/mlprec/impl/aggregator/mld_zaggrmat_biz_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_zaggrmat_biz_asb.F90 +! File: mld_zaggrmat_biz_bld.F90 ! -! Subroutine: mld_zaggrmat_biz_asb +! Subroutine: mld_zaggrmat_biz_bld ! Version: complex ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -57,7 +57,7 @@ ! specified by the user through mld_zprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! mld_z_lev_aggrmat_asb. +! mld_z_lev_aggrmat_bld. ! ! Arguments: ! a - type(psb_zspmat_type), input. @@ -80,10 +80,10 @@ ! info - integer, output. ! Error code. ! -subroutine mld_zaggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_zaggrmat_biz_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_z_inner_mod, mld_protect_name => mld_zaggrmat_biz_asb + use mld_z_inner_mod, mld_protect_name => mld_zaggrmat_biz_bld implicit none @@ -111,7 +111,7 @@ subroutine mld_zaggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr integer(psb_ipk_), parameter :: ncmax=16 real(psb_dpk_) :: anorm, omega, tmp, dg, theta - name='mld_aggrmat_biz_asb' + name='mld_aggrmat_biz_bld' info=psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -372,4 +372,4 @@ subroutine mld_zaggrmat_biz_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr call psb_error_handler(err_act) return -end subroutine mld_zaggrmat_biz_asb +end subroutine mld_zaggrmat_biz_bld diff --git a/mlprec/impl/aggregator/mld_zaggrmat_minnrg_asb.f90 b/mlprec/impl/aggregator/mld_zaggrmat_minnrg_bld.f90 similarity index 98% rename from mlprec/impl/aggregator/mld_zaggrmat_minnrg_asb.f90 rename to mlprec/impl/aggregator/mld_zaggrmat_minnrg_bld.f90 index 043855aa..76bdcc7c 100644 --- a/mlprec/impl/aggregator/mld_zaggrmat_minnrg_asb.f90 +++ b/mlprec/impl/aggregator/mld_zaggrmat_minnrg_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_zaggrmat_minnrg_asb.F90 +! File: mld_zaggrmat_minnrg_bld.F90 ! -! Subroutine: mld_zaggrmat_minnrg_asb +! Subroutine: mld_zaggrmat_minnrg_bld ! Version: complex ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -65,7 +65,7 @@ ! ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! ! Arguments: @@ -104,10 +104,10 @@ ! Error code. ! ! -subroutine mld_zaggrmat_minnrg_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_zaggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_z_inner_mod, mld_protect_name => mld_zaggrmat_minnrg_asb + use mld_z_inner_mod, mld_protect_name => mld_zaggrmat_minnrg_bld implicit none @@ -651,4 +651,4 @@ contains close(20+me) end subroutine local_dump -end subroutine mld_zaggrmat_minnrg_asb +end subroutine mld_zaggrmat_minnrg_bld diff --git a/mlprec/impl/aggregator/mld_zaggrmat_nosmth_asb.f90 b/mlprec/impl/aggregator/mld_zaggrmat_nosmth_bld.f90 similarity index 96% rename from mlprec/impl/aggregator/mld_zaggrmat_nosmth_asb.f90 rename to mlprec/impl/aggregator/mld_zaggrmat_nosmth_bld.f90 index 152b2676..13eb1423 100644 --- a/mlprec/impl/aggregator/mld_zaggrmat_nosmth_asb.f90 +++ b/mlprec/impl/aggregator/mld_zaggrmat_nosmth_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_zaggrmat_nosmth_asb.F90 +! File: mld_zaggrmat_nosmth_bld.F90 ! -! Subroutine: mld_zaggrmat_nosmth_asb +! Subroutine: mld_zaggrmat_nosmth_bld ! Version: complex ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -53,7 +53,7 @@ ! specified by the user through mld_zprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! For details see ! P. D'Ambra, D. di Serafino and S. Filippone, On the development of @@ -96,10 +96,10 @@ ! Error code. ! ! -subroutine mld_zaggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_zaggrmat_nosmth_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_z_inner_mod, mld_protect_name => mld_zaggrmat_nosmth_asb + use mld_z_inner_mod, mld_protect_name => mld_zaggrmat_nosmth_bld implicit none @@ -124,7 +124,7 @@ subroutine mld_zaggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_re integer(psb_lpk_) :: nrow, nglob, ncol, ntaggr, nzl, ip, & & naggr, nzt, naggrm1, naggrp1, i, k - name = 'mld_aggrmat_nosmth_asb' + name = 'mld_aggrmat_nosmth_bld' info = psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -198,4 +198,4 @@ subroutine mld_zaggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_re return -end subroutine mld_zaggrmat_nosmth_asb +end subroutine mld_zaggrmat_nosmth_bld diff --git a/mlprec/impl/aggregator/mld_zaggrmat_smth_asb.f90 b/mlprec/impl/aggregator/mld_zaggrmat_smth_bld.f90 similarity index 97% rename from mlprec/impl/aggregator/mld_zaggrmat_smth_asb.f90 rename to mlprec/impl/aggregator/mld_zaggrmat_smth_bld.f90 index 80376ca5..c96b5d12 100644 --- a/mlprec/impl/aggregator/mld_zaggrmat_smth_asb.f90 +++ b/mlprec/impl/aggregator/mld_zaggrmat_smth_bld.f90 @@ -35,9 +35,9 @@ ! POSSIBILITY OF SUCH DAMAGE. ! ! -! File: mld_zaggrmat_smth_asb.F90 +! File: mld_zaggrmat_smth_bld.F90 ! -! Subroutine: mld_zaggrmat_smth_asb +! Subroutine: mld_zaggrmat_smth_bld ! Version: complex ! ! This routine builds a coarse-level matrix A_C from a fine-level matrix A @@ -65,7 +65,7 @@ ! specified by the user through mld_zprecinit and mld_zprecset. ! On output from this routine the entries of AC, op_prol, op_restr ! are still in "global numbering" mode; this is fixed in the calling routine -! aggregator%mat_asb. +! aggregator%mat_bld. ! ! ! Arguments: @@ -102,10 +102,10 @@ ! info - integer, output. ! Error code. ! -subroutine mld_zaggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) +subroutine mld_zaggrmat_smth_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) use psb_base_mod use mld_base_prec_type - use mld_z_inner_mod, mld_protect_name => mld_zaggrmat_smth_asb + use mld_z_inner_mod, mld_protect_name => mld_zaggrmat_smth_bld implicit none @@ -133,7 +133,7 @@ subroutine mld_zaggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_rest integer(psb_ipk_), parameter :: ncmax=16 real(psb_dpk_) :: anorm, omega, tmp, dg, theta - name='mld_aggrmat_smth_asb' + name='mld_aggrmat_smth_bld' info=psb_success_ call psb_erractionsave(err_act) if (psb_errstatus_fatal()) then @@ -411,4 +411,4 @@ subroutine mld_zaggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_rest call psb_error_handler(err_act) return -end subroutine mld_zaggrmat_smth_asb +end subroutine mld_zaggrmat_smth_bld diff --git a/mlprec/impl/level/mld_c_base_onelev_mat_asb.f90 b/mlprec/impl/level/mld_c_base_onelev_mat_asb.f90 index e5643570..3e0794bb 100644 --- a/mlprec/impl/level/mld_c_base_onelev_mat_asb.f90 +++ b/mlprec/impl/level/mld_c_base_onelev_mat_asb.f90 @@ -142,119 +142,28 @@ subroutine mld_c_base_onelev_mat_asb(lv,a,desc_a,ilaggr,nlaggr,op_prol,info) ! the mapping defined by mld_aggrmap_bld and applying the aggregation ! algorithm specified by lv%iprcparm(mld_aggr_prol_) ! - call lv%aggr%mat_asb(lv%parms,a,desc_a,ilaggr,nlaggr,lac,op_prol,op_restr,info) + call lv%aggr%mat_bld(lv%parms,a,desc_a,ilaggr,nlaggr,lac,op_prol,op_restr,info) if(info /= psb_success_) then call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_aggrmat_asb') goto 9999 end if - - ! Common code refactored here. - - ntaggr = sum(nlaggr) - - select case(lv%parms%coarse_mat) - - case(mld_distr_mat_) - - call lac%mv_to(bcoo) - nzl = bcoo%get_nzeros() - inl = nlaggr(me+1) - if (inl < nlaggr(me+1)) then - info = psb_err_bad_int_cnv_ - call psb_errpush(info,name,& - & l_err=(/nlaggr(me+1),inl*1_psb_lpk_/)) - goto 9999 - end if - if (info == psb_success_) call psb_cdall(ictxt,lv%desc_ac,info,nl=inl) - if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,lv%desc_ac,info) - if (info == psb_success_) call psb_cdasb(lv%desc_ac,info) - if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),lv%desc_ac,info,iact='I') - if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),lv%desc_ac,info,iact='I') - if (info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,& - & a_err='Creating lv%desc_ac and converting ac') - goto 9999 - end if - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),& - & 'Assembld aux descr. distr.' - call lac%mv_from(bcoo) - call lv%ac%mv_from_l(lac) - - call lv%ac%set_nrows(lv%desc_ac%get_local_rows()) - call lv%ac%set_ncols(lv%desc_ac%get_local_cols()) - call lv%ac%set_asb() - - if (info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free') - goto 9999 - end if - - if (np>1) then - call op_prol%mv_to(acsr1) - nzl = acsr1%get_nzeros() - call psb_glob_to_loc(acsr1%ja(1:nzl),lv%desc_ac,info,'I') - if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc') - goto 9999 - end if - call op_prol%mv_from(acsr1) - endif - nc = lv%desc_ac%get_local_cols() - call op_prol%set_ncols(nc) - - if (np>1) then - call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_) - call op_restr%mv_to(acoo) - nzl = acoo%get_nzeros() - if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),lv%desc_ac,info,'I') - call acoo%set_dupl(psb_dupl_add_) - if (info == psb_success_) call op_restr%mv_from(acoo) - if (info == psb_success_) call op_restr%cscnv(info,type='csr') - if(info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,& - & a_err='Converting op_restr to local') - goto 9999 - end if - end if - ! - ! Clip to local rows. - ! - nr = lv%desc_ac%get_local_rows() - call op_restr%set_nrows(nr) - - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),& - & 'Done ac ' - - case(mld_repl_mat_) - ! - ! - - call psb_cdall(ictxt,lv%desc_ac,info,mg=ntaggr,repl=.true.) - if (info == psb_success_) call psb_cdasb(lv%desc_ac,info) - if (info == psb_success_) & - & call psb_gather(lac1,lac,lv%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) - call lv%ac%mv_from_l(lac1) - if (info /= psb_success_) goto 9999 - - case default - info = psb_err_internal_error_ - call psb_errpush(info,name,a_err='invalid mld_coarse_mat_') - goto 9999 - end select - - call lv%ac%cscnv(info,type='csr',dupl=psb_dupl_add_) - if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='spcnv') - goto 9999 - end if - - call lv%aggr%bld_map(desc_a, lv%desc_ac,ilaggr,nlaggr,op_restr,op_prol,lv%map,info) + ! + ! Now build its descriptor and convert global indices for + ! ac, op_restr and op_prol + ! + if (info == psb_success_) & + & call lv%aggr%mat_asb(lv%parms,a,desc_a,ilaggr,nlaggr,& + & lac,lv%desc_ac,op_prol,op_restr,info) + + call lv%ac%mv_from_l(lac) + if (info == psb_success_) call lv%ac%cscnv(info,type='csr',dupl=psb_dupl_add_) + + if (info == psb_success_) call lv%aggr%bld_map(desc_a, lv%desc_ac,& + & ilaggr,nlaggr,op_restr,op_prol,lv%map,info) if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='bld_map') + call psb_errpush(psb_err_from_subroutine_,name,a_err='mat_asb/map_bld') goto 9999 end if ! diff --git a/mlprec/impl/level/mld_d_base_onelev_mat_asb.f90 b/mlprec/impl/level/mld_d_base_onelev_mat_asb.f90 index 81b32433..2a192b9e 100644 --- a/mlprec/impl/level/mld_d_base_onelev_mat_asb.f90 +++ b/mlprec/impl/level/mld_d_base_onelev_mat_asb.f90 @@ -142,119 +142,28 @@ subroutine mld_d_base_onelev_mat_asb(lv,a,desc_a,ilaggr,nlaggr,op_prol,info) ! the mapping defined by mld_aggrmap_bld and applying the aggregation ! algorithm specified by lv%iprcparm(mld_aggr_prol_) ! - call lv%aggr%mat_asb(lv%parms,a,desc_a,ilaggr,nlaggr,lac,op_prol,op_restr,info) + call lv%aggr%mat_bld(lv%parms,a,desc_a,ilaggr,nlaggr,lac,op_prol,op_restr,info) if(info /= psb_success_) then call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_aggrmat_asb') goto 9999 end if - - ! Common code refactored here. - - ntaggr = sum(nlaggr) - - select case(lv%parms%coarse_mat) - - case(mld_distr_mat_) - - call lac%mv_to(bcoo) - nzl = bcoo%get_nzeros() - inl = nlaggr(me+1) - if (inl < nlaggr(me+1)) then - info = psb_err_bad_int_cnv_ - call psb_errpush(info,name,& - & l_err=(/nlaggr(me+1),inl*1_psb_lpk_/)) - goto 9999 - end if - if (info == psb_success_) call psb_cdall(ictxt,lv%desc_ac,info,nl=inl) - if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,lv%desc_ac,info) - if (info == psb_success_) call psb_cdasb(lv%desc_ac,info) - if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),lv%desc_ac,info,iact='I') - if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),lv%desc_ac,info,iact='I') - if (info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,& - & a_err='Creating lv%desc_ac and converting ac') - goto 9999 - end if - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),& - & 'Assembld aux descr. distr.' - call lac%mv_from(bcoo) - call lv%ac%mv_from_l(lac) - - call lv%ac%set_nrows(lv%desc_ac%get_local_rows()) - call lv%ac%set_ncols(lv%desc_ac%get_local_cols()) - call lv%ac%set_asb() - - if (info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free') - goto 9999 - end if - - if (np>1) then - call op_prol%mv_to(acsr1) - nzl = acsr1%get_nzeros() - call psb_glob_to_loc(acsr1%ja(1:nzl),lv%desc_ac,info,'I') - if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc') - goto 9999 - end if - call op_prol%mv_from(acsr1) - endif - nc = lv%desc_ac%get_local_cols() - call op_prol%set_ncols(nc) - - if (np>1) then - call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_) - call op_restr%mv_to(acoo) - nzl = acoo%get_nzeros() - if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),lv%desc_ac,info,'I') - call acoo%set_dupl(psb_dupl_add_) - if (info == psb_success_) call op_restr%mv_from(acoo) - if (info == psb_success_) call op_restr%cscnv(info,type='csr') - if(info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,& - & a_err='Converting op_restr to local') - goto 9999 - end if - end if - ! - ! Clip to local rows. - ! - nr = lv%desc_ac%get_local_rows() - call op_restr%set_nrows(nr) - - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),& - & 'Done ac ' - - case(mld_repl_mat_) - ! - ! - - call psb_cdall(ictxt,lv%desc_ac,info,mg=ntaggr,repl=.true.) - if (info == psb_success_) call psb_cdasb(lv%desc_ac,info) - if (info == psb_success_) & - & call psb_gather(lac1,lac,lv%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) - call lv%ac%mv_from_l(lac1) - if (info /= psb_success_) goto 9999 - - case default - info = psb_err_internal_error_ - call psb_errpush(info,name,a_err='invalid mld_coarse_mat_') - goto 9999 - end select - - call lv%ac%cscnv(info,type='csr',dupl=psb_dupl_add_) - if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='spcnv') - goto 9999 - end if - - call lv%aggr%bld_map(desc_a, lv%desc_ac,ilaggr,nlaggr,op_restr,op_prol,lv%map,info) + ! + ! Now build its descriptor and convert global indices for + ! ac, op_restr and op_prol + ! + if (info == psb_success_) & + & call lv%aggr%mat_asb(lv%parms,a,desc_a,ilaggr,nlaggr,& + & lac,lv%desc_ac,op_prol,op_restr,info) + + call lv%ac%mv_from_l(lac) + if (info == psb_success_) call lv%ac%cscnv(info,type='csr',dupl=psb_dupl_add_) + + if (info == psb_success_) call lv%aggr%bld_map(desc_a, lv%desc_ac,& + & ilaggr,nlaggr,op_restr,op_prol,lv%map,info) if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='bld_map') + call psb_errpush(psb_err_from_subroutine_,name,a_err='mat_asb/map_bld') goto 9999 end if ! diff --git a/mlprec/impl/level/mld_s_base_onelev_mat_asb.f90 b/mlprec/impl/level/mld_s_base_onelev_mat_asb.f90 index c6ed12ca..cdbfa918 100644 --- a/mlprec/impl/level/mld_s_base_onelev_mat_asb.f90 +++ b/mlprec/impl/level/mld_s_base_onelev_mat_asb.f90 @@ -142,119 +142,28 @@ subroutine mld_s_base_onelev_mat_asb(lv,a,desc_a,ilaggr,nlaggr,op_prol,info) ! the mapping defined by mld_aggrmap_bld and applying the aggregation ! algorithm specified by lv%iprcparm(mld_aggr_prol_) ! - call lv%aggr%mat_asb(lv%parms,a,desc_a,ilaggr,nlaggr,lac,op_prol,op_restr,info) + call lv%aggr%mat_bld(lv%parms,a,desc_a,ilaggr,nlaggr,lac,op_prol,op_restr,info) if(info /= psb_success_) then call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_aggrmat_asb') goto 9999 end if - - ! Common code refactored here. - - ntaggr = sum(nlaggr) - - select case(lv%parms%coarse_mat) - - case(mld_distr_mat_) - - call lac%mv_to(bcoo) - nzl = bcoo%get_nzeros() - inl = nlaggr(me+1) - if (inl < nlaggr(me+1)) then - info = psb_err_bad_int_cnv_ - call psb_errpush(info,name,& - & l_err=(/nlaggr(me+1),inl*1_psb_lpk_/)) - goto 9999 - end if - if (info == psb_success_) call psb_cdall(ictxt,lv%desc_ac,info,nl=inl) - if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,lv%desc_ac,info) - if (info == psb_success_) call psb_cdasb(lv%desc_ac,info) - if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),lv%desc_ac,info,iact='I') - if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),lv%desc_ac,info,iact='I') - if (info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,& - & a_err='Creating lv%desc_ac and converting ac') - goto 9999 - end if - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),& - & 'Assembld aux descr. distr.' - call lac%mv_from(bcoo) - call lv%ac%mv_from_l(lac) - - call lv%ac%set_nrows(lv%desc_ac%get_local_rows()) - call lv%ac%set_ncols(lv%desc_ac%get_local_cols()) - call lv%ac%set_asb() - - if (info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free') - goto 9999 - end if - - if (np>1) then - call op_prol%mv_to(acsr1) - nzl = acsr1%get_nzeros() - call psb_glob_to_loc(acsr1%ja(1:nzl),lv%desc_ac,info,'I') - if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc') - goto 9999 - end if - call op_prol%mv_from(acsr1) - endif - nc = lv%desc_ac%get_local_cols() - call op_prol%set_ncols(nc) - - if (np>1) then - call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_) - call op_restr%mv_to(acoo) - nzl = acoo%get_nzeros() - if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),lv%desc_ac,info,'I') - call acoo%set_dupl(psb_dupl_add_) - if (info == psb_success_) call op_restr%mv_from(acoo) - if (info == psb_success_) call op_restr%cscnv(info,type='csr') - if(info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,& - & a_err='Converting op_restr to local') - goto 9999 - end if - end if - ! - ! Clip to local rows. - ! - nr = lv%desc_ac%get_local_rows() - call op_restr%set_nrows(nr) - - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),& - & 'Done ac ' - - case(mld_repl_mat_) - ! - ! - - call psb_cdall(ictxt,lv%desc_ac,info,mg=ntaggr,repl=.true.) - if (info == psb_success_) call psb_cdasb(lv%desc_ac,info) - if (info == psb_success_) & - & call psb_gather(lac1,lac,lv%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) - call lv%ac%mv_from_l(lac1) - if (info /= psb_success_) goto 9999 - - case default - info = psb_err_internal_error_ - call psb_errpush(info,name,a_err='invalid mld_coarse_mat_') - goto 9999 - end select - - call lv%ac%cscnv(info,type='csr',dupl=psb_dupl_add_) - if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='spcnv') - goto 9999 - end if - - call lv%aggr%bld_map(desc_a, lv%desc_ac,ilaggr,nlaggr,op_restr,op_prol,lv%map,info) + ! + ! Now build its descriptor and convert global indices for + ! ac, op_restr and op_prol + ! + if (info == psb_success_) & + & call lv%aggr%mat_asb(lv%parms,a,desc_a,ilaggr,nlaggr,& + & lac,lv%desc_ac,op_prol,op_restr,info) + + call lv%ac%mv_from_l(lac) + if (info == psb_success_) call lv%ac%cscnv(info,type='csr',dupl=psb_dupl_add_) + + if (info == psb_success_) call lv%aggr%bld_map(desc_a, lv%desc_ac,& + & ilaggr,nlaggr,op_restr,op_prol,lv%map,info) if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='bld_map') + call psb_errpush(psb_err_from_subroutine_,name,a_err='mat_asb/map_bld') goto 9999 end if ! diff --git a/mlprec/impl/level/mld_z_base_onelev_mat_asb.f90 b/mlprec/impl/level/mld_z_base_onelev_mat_asb.f90 index 7e57e3b2..449922e9 100644 --- a/mlprec/impl/level/mld_z_base_onelev_mat_asb.f90 +++ b/mlprec/impl/level/mld_z_base_onelev_mat_asb.f90 @@ -142,119 +142,28 @@ subroutine mld_z_base_onelev_mat_asb(lv,a,desc_a,ilaggr,nlaggr,op_prol,info) ! the mapping defined by mld_aggrmap_bld and applying the aggregation ! algorithm specified by lv%iprcparm(mld_aggr_prol_) ! - call lv%aggr%mat_asb(lv%parms,a,desc_a,ilaggr,nlaggr,lac,op_prol,op_restr,info) + call lv%aggr%mat_bld(lv%parms,a,desc_a,ilaggr,nlaggr,lac,op_prol,op_restr,info) if(info /= psb_success_) then call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_aggrmat_asb') goto 9999 end if - - ! Common code refactored here. - - ntaggr = sum(nlaggr) - - select case(lv%parms%coarse_mat) - - case(mld_distr_mat_) - - call lac%mv_to(bcoo) - nzl = bcoo%get_nzeros() - inl = nlaggr(me+1) - if (inl < nlaggr(me+1)) then - info = psb_err_bad_int_cnv_ - call psb_errpush(info,name,& - & l_err=(/nlaggr(me+1),inl*1_psb_lpk_/)) - goto 9999 - end if - if (info == psb_success_) call psb_cdall(ictxt,lv%desc_ac,info,nl=inl) - if (info == psb_success_) call psb_cdins(nzl,bcoo%ia,bcoo%ja,lv%desc_ac,info) - if (info == psb_success_) call psb_cdasb(lv%desc_ac,info) - if (info == psb_success_) call psb_glob_to_loc(bcoo%ia(1:nzl),lv%desc_ac,info,iact='I') - if (info == psb_success_) call psb_glob_to_loc(bcoo%ja(1:nzl),lv%desc_ac,info,iact='I') - if (info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,& - & a_err='Creating lv%desc_ac and converting ac') - goto 9999 - end if - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),& - & 'Assembld aux descr. distr.' - call lac%mv_from(bcoo) - call lv%ac%mv_from_l(lac) - - call lv%ac%set_nrows(lv%desc_ac%get_local_rows()) - call lv%ac%set_ncols(lv%desc_ac%get_local_cols()) - call lv%ac%set_asb() - - if (info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_sp_free') - goto 9999 - end if - - if (np>1) then - call op_prol%mv_to(acsr1) - nzl = acsr1%get_nzeros() - call psb_glob_to_loc(acsr1%ja(1:nzl),lv%desc_ac,info,'I') - if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_glob_to_loc') - goto 9999 - end if - call op_prol%mv_from(acsr1) - endif - nc = lv%desc_ac%get_local_cols() - call op_prol%set_ncols(nc) - - if (np>1) then - call op_restr%cscnv(info,type='coo',dupl=psb_dupl_add_) - call op_restr%mv_to(acoo) - nzl = acoo%get_nzeros() - if (info == psb_success_) call psb_glob_to_loc(acoo%ia(1:nzl),lv%desc_ac,info,'I') - call acoo%set_dupl(psb_dupl_add_) - if (info == psb_success_) call op_restr%mv_from(acoo) - if (info == psb_success_) call op_restr%cscnv(info,type='csr') - if(info /= psb_success_) then - call psb_errpush(psb_err_internal_error_,name,& - & a_err='Converting op_restr to local') - goto 9999 - end if - end if - ! - ! Clip to local rows. - ! - nr = lv%desc_ac%get_local_rows() - call op_restr%set_nrows(nr) - - if (debug_level >= psb_debug_outer_) & - & write(debug_unit,*) me,' ',trim(name),& - & 'Done ac ' - - case(mld_repl_mat_) - ! - ! - - call psb_cdall(ictxt,lv%desc_ac,info,mg=ntaggr,repl=.true.) - if (info == psb_success_) call psb_cdasb(lv%desc_ac,info) - if (info == psb_success_) & - & call psb_gather(lac1,lac,lv%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) - call lv%ac%mv_from_l(lac1) - if (info /= psb_success_) goto 9999 - - case default - info = psb_err_internal_error_ - call psb_errpush(info,name,a_err='invalid mld_coarse_mat_') - goto 9999 - end select - - call lv%ac%cscnv(info,type='csr',dupl=psb_dupl_add_) - if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='spcnv') - goto 9999 - end if - - call lv%aggr%bld_map(desc_a, lv%desc_ac,ilaggr,nlaggr,op_restr,op_prol,lv%map,info) + ! + ! Now build its descriptor and convert global indices for + ! ac, op_restr and op_prol + ! + if (info == psb_success_) & + & call lv%aggr%mat_asb(lv%parms,a,desc_a,ilaggr,nlaggr,& + & lac,lv%desc_ac,op_prol,op_restr,info) + + call lv%ac%mv_from_l(lac) + if (info == psb_success_) call lv%ac%cscnv(info,type='csr',dupl=psb_dupl_add_) + + if (info == psb_success_) call lv%aggr%bld_map(desc_a, lv%desc_ac,& + & ilaggr,nlaggr,op_restr,op_prol,lv%map,info) if(info /= psb_success_) then - call psb_errpush(psb_err_from_subroutine_,name,a_err='bld_map') + call psb_errpush(psb_err_from_subroutine_,name,a_err='mat_asb/map_bld') goto 9999 end if ! diff --git a/mlprec/mld_c_base_aggregator_mod.f90 b/mlprec/mld_c_base_aggregator_mod.f90 index 79d7d0e3..dcaee1d1 100644 --- a/mlprec/mld_c_base_aggregator_mod.f90 +++ b/mlprec/mld_c_base_aggregator_mod.f90 @@ -64,8 +64,10 @@ module mld_c_base_aggregator_mod !! !! bld_tprol - Build a tentative prolongator !! - !! mat_asb - Build the final prolongator/restrictor and the - !! coarse matrix ac + !! mat_bld - Build prolongator/restrictor and coarse matrix ac + !! + !! mat_asb - Convert prolongator/restrictor/coarse matrix + !! and fix their descriptor(s) !! !! update_next - Transfer information to the next level; default is !! to do nothing, i.e. aggregators at different @@ -83,6 +85,7 @@ module mld_c_base_aggregator_mod logical :: do_clean_zeros contains procedure, pass(ag) :: bld_tprol => mld_c_base_aggregator_build_tprol + procedure, pass(ag) :: mat_bld => mld_c_base_aggregator_mat_bld procedure, pass(ag) :: mat_asb => mld_c_base_aggregator_mat_asb procedure, pass(ag) :: bld_map => mld_c_base_aggregator_bld_map procedure, pass(ag) :: update_next => mld_c_base_aggregator_update_next @@ -273,8 +276,8 @@ contains class(mld_c_base_aggregator_type), target, intent(inout) :: ag type(mld_sml_parms), intent(inout) :: parms type(mld_saggr_data), intent(in) :: ag_data - type(psb_cspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_cspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lcspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info @@ -298,7 +301,7 @@ contains end subroutine mld_c_base_aggregator_build_tprol ! - !> Function mat_asb + !> Function mat_bld !! \memberof mld_c_base_aggregator_type !! \brief Build prolongator/restrictor/coarse matrix. !! @@ -316,7 +319,7 @@ contains !! in many cases it is the transpose of the prolongator. !! \param info Return code !! - subroutine mld_c_base_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& + subroutine mld_c_base_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& & op_prol,op_restr,info) use psb_base_mod implicit none @@ -329,6 +332,54 @@ contains type(psb_lcspmat_type), intent(out) :: ac,op_restr integer(psb_ipk_), intent(out) :: info integer(psb_ipk_) :: err_act + character(len=20) :: name='c_base_aggregator_mat_bld' + + call psb_erractionsave(err_act) + + info = psb_err_missing_override_method_ + call psb_errpush(info,name) + goto 9999 + + call psb_erractionrestore(err_act) + return + +9999 call psb_error_handler(err_act) + + return + end subroutine mld_c_base_aggregator_mat_bld + + ! + !> Function mat_asb + !! \memberof mld_c_base_aggregator_type + !! \brief Build prolongator/restrictor/coarse matrix. + !! + !! + !! \param ag The input aggregator object + !! \param parms The auxiliary parameters object + !! \param a The local matrix part + !! \param desc_a The descriptor + !! \param ilaggr Aggregation map + !! \param nlaggr Sizes of ilaggr on all processes + !! \param ac On output the coarse matrix + !! \param op_prol On input, the tentative prolongator operator, on output + !! the final prolongator + !! \param op_restr On output, the restrictor operator; + !! in many cases it is the transpose of the prolongator. + !! \param info Return code + !! + subroutine mld_c_base_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac, op_prol,op_restr,info) + use psb_base_mod + implicit none + class(mld_c_base_aggregator_type), target, intent(inout) :: ag + type(mld_sml_parms), intent(inout) :: parms + type(psb_cspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lcspmat_type), intent(inout) :: op_prol,ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_) :: err_act character(len=20) :: name='c_base_aggregator_mat_asb' call psb_erractionsave(err_act) diff --git a/mlprec/mld_c_dec_aggregator_mod.f90 b/mlprec/mld_c_dec_aggregator_mod.f90 index 40070838..30620c1e 100644 --- a/mlprec/mld_c_dec_aggregator_mod.f90 +++ b/mlprec/mld_c_dec_aggregator_mod.f90 @@ -80,6 +80,7 @@ module mld_c_dec_aggregator_mod contains procedure, pass(ag) :: bld_tprol => mld_c_dec_aggregator_build_tprol + procedure, pass(ag) :: mat_bld => mld_c_dec_aggregator_mat_bld procedure, pass(ag) :: mat_asb => mld_c_dec_aggregator_mat_asb procedure, pass(ag) :: default => mld_c_dec_aggregator_default procedure, pass(ag) :: set_aggr_type => mld_c_dec_aggregator_set_aggr_type @@ -99,8 +100,8 @@ module mld_c_dec_aggregator_mod class(mld_c_dec_aggregator_type), target, intent(inout) :: ag type(mld_sml_parms), intent(inout) :: parms type(mld_saggr_data), intent(in) :: ag_data - type(psb_cspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_cspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lcspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info @@ -108,7 +109,7 @@ module mld_c_dec_aggregator_mod end interface interface - subroutine mld_c_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& + subroutine mld_c_dec_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& & op_prol,op_restr,info) import :: mld_c_dec_aggregator_type, psb_desc_type, psb_cspmat_type, psb_spk_, & & psb_ipk_, psb_lpk_, psb_lcspmat_type, mld_sml_parms @@ -121,9 +122,25 @@ module mld_c_dec_aggregator_mod type(psb_lcspmat_type), intent(inout) :: op_prol type(psb_lcspmat_type), intent(out) :: ac,op_restr integer(psb_ipk_), intent(out) :: info - end subroutine mld_c_dec_aggregator_mat_asb + end subroutine mld_c_dec_aggregator_mat_bld end interface + interface + subroutine mld_c_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac,op_prol,op_restr,info) + import :: mld_c_dec_aggregator_type, psb_desc_type, psb_cspmat_type, psb_spk_, & + & psb_ipk_, psb_lpk_, psb_lcspmat_type, mld_sml_parms + implicit none + class(mld_c_dec_aggregator_type), target, intent(inout) :: ag + type(mld_sml_parms), intent(inout) :: parms + type(psb_cspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lcspmat_type), intent(inout) :: op_prol,ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + end subroutine mld_c_dec_aggregator_mat_asb + end interface contains diff --git a/mlprec/mld_c_inner_mod.f90 b/mlprec/mld_c_inner_mod.f90 index 9a141310..d5b7b50a 100644 --- a/mlprec/mld_c_inner_mod.f90 +++ b/mlprec/mld_c_inner_mod.f90 @@ -108,37 +108,9 @@ module mld_c_inner_mod integer(psb_ipk_), intent(out) :: info end subroutine mld_c_map_to_tprol end interface mld_map_to_tprol - - interface mld_lev_mat_asb - subroutine mld_c_lev_aggrmat_asb(p,a,desc_a,ilaggr,nlaggr,op_prol,info) - import :: psb_cspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lcspmat_type - import :: mld_c_onelev_type - implicit none - type(mld_c_onelev_type), intent(inout), target :: p - type(psb_cspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:),nlaggr(:) - type(psb_lcspmat_type), intent(inout) :: op_prol - integer(psb_ipk_), intent(out) :: info - end subroutine mld_c_lev_aggrmat_asb - end interface mld_lev_mat_asb - - interface mld_aggrmat_asb - subroutine mld_caggrmat_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) - import :: psb_cspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lcspmat_type - import :: mld_sml_parms - implicit none - type(psb_cspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) - type(mld_sml_parms), intent(inout) :: parms - type(psb_lcspmat_type), intent(out) :: ac,op_prol,op_restr - integer(psb_ipk_), intent(out) :: info - end subroutine mld_caggrmat_asb - end interface mld_aggrmat_asb abstract interface - subroutine mld_caggrmat_var_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) + subroutine mld_caggrmat_var_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) import :: psb_cspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lcspmat_type import :: mld_c_onelev_type, mld_sml_parms implicit none @@ -148,11 +120,11 @@ module mld_c_inner_mod type(mld_sml_parms), intent(inout) :: parms type(psb_lcspmat_type), intent(out) :: ac,op_prol,op_restr integer(psb_ipk_), intent(out) :: info - end subroutine mld_caggrmat_var_asb + end subroutine mld_caggrmat_var_bld end interface - procedure(mld_caggrmat_var_asb) :: mld_caggrmat_nosmth_asb, & - & mld_caggrmat_smth_asb, mld_caggrmat_minnrg_asb, & - & mld_caggrmat_biz_asb + procedure(mld_caggrmat_var_bld) :: mld_caggrmat_nosmth_bld, & + & mld_caggrmat_smth_bld, mld_caggrmat_minnrg_bld, & + & mld_caggrmat_biz_bld end module mld_c_inner_mod diff --git a/mlprec/mld_c_onelev_mod.f90 b/mlprec/mld_c_onelev_mod.f90 index 2ab4253b..e8223d2b 100644 --- a/mlprec/mld_c_onelev_mod.f90 +++ b/mlprec/mld_c_onelev_mod.f90 @@ -488,8 +488,8 @@ contains & ilaggr,nlaggr,op_prol,ag_data,info) implicit none class(mld_c_onelev_type), intent(inout), target :: lv - type(psb_cspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_cspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lcspmat_type), intent(out) :: op_prol type(mld_saggr_data), intent(in) :: ag_data diff --git a/mlprec/mld_c_symdec_aggregator_mod.f90 b/mlprec/mld_c_symdec_aggregator_mod.f90 index b96dfa8d..99e31af6 100644 --- a/mlprec/mld_c_symdec_aggregator_mod.f90 +++ b/mlprec/mld_c_symdec_aggregator_mod.f90 @@ -70,8 +70,8 @@ module mld_c_symdec_aggregator_mod class(mld_c_symdec_aggregator_type), target, intent(inout) :: ag type(mld_sml_parms), intent(inout) :: parms type(mld_saggr_data), intent(in) :: ag_data - type(psb_cspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_cspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lcspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info diff --git a/mlprec/mld_d_base_aggregator_mod.f90 b/mlprec/mld_d_base_aggregator_mod.f90 index 04b58f48..49307c4f 100644 --- a/mlprec/mld_d_base_aggregator_mod.f90 +++ b/mlprec/mld_d_base_aggregator_mod.f90 @@ -64,8 +64,10 @@ module mld_d_base_aggregator_mod !! !! bld_tprol - Build a tentative prolongator !! - !! mat_asb - Build the final prolongator/restrictor and the - !! coarse matrix ac + !! mat_bld - Build prolongator/restrictor and coarse matrix ac + !! + !! mat_asb - Convert prolongator/restrictor/coarse matrix + !! and fix their descriptor(s) !! !! update_next - Transfer information to the next level; default is !! to do nothing, i.e. aggregators at different @@ -83,6 +85,7 @@ module mld_d_base_aggregator_mod logical :: do_clean_zeros contains procedure, pass(ag) :: bld_tprol => mld_d_base_aggregator_build_tprol + procedure, pass(ag) :: mat_bld => mld_d_base_aggregator_mat_bld procedure, pass(ag) :: mat_asb => mld_d_base_aggregator_mat_asb procedure, pass(ag) :: bld_map => mld_d_base_aggregator_bld_map procedure, pass(ag) :: update_next => mld_d_base_aggregator_update_next @@ -273,8 +276,8 @@ contains class(mld_d_base_aggregator_type), target, intent(inout) :: ag type(mld_dml_parms), intent(inout) :: parms type(mld_daggr_data), intent(in) :: ag_data - type(psb_dspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_dspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_ldspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info @@ -298,7 +301,7 @@ contains end subroutine mld_d_base_aggregator_build_tprol ! - !> Function mat_asb + !> Function mat_bld !! \memberof mld_d_base_aggregator_type !! \brief Build prolongator/restrictor/coarse matrix. !! @@ -316,7 +319,7 @@ contains !! in many cases it is the transpose of the prolongator. !! \param info Return code !! - subroutine mld_d_base_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& + subroutine mld_d_base_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& & op_prol,op_restr,info) use psb_base_mod implicit none @@ -329,6 +332,54 @@ contains type(psb_ldspmat_type), intent(out) :: ac,op_restr integer(psb_ipk_), intent(out) :: info integer(psb_ipk_) :: err_act + character(len=20) :: name='d_base_aggregator_mat_bld' + + call psb_erractionsave(err_act) + + info = psb_err_missing_override_method_ + call psb_errpush(info,name) + goto 9999 + + call psb_erractionrestore(err_act) + return + +9999 call psb_error_handler(err_act) + + return + end subroutine mld_d_base_aggregator_mat_bld + + ! + !> Function mat_asb + !! \memberof mld_d_base_aggregator_type + !! \brief Build prolongator/restrictor/coarse matrix. + !! + !! + !! \param ag The input aggregator object + !! \param parms The auxiliary parameters object + !! \param a The local matrix part + !! \param desc_a The descriptor + !! \param ilaggr Aggregation map + !! \param nlaggr Sizes of ilaggr on all processes + !! \param ac On output the coarse matrix + !! \param op_prol On input, the tentative prolongator operator, on output + !! the final prolongator + !! \param op_restr On output, the restrictor operator; + !! in many cases it is the transpose of the prolongator. + !! \param info Return code + !! + subroutine mld_d_base_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac, op_prol,op_restr,info) + use psb_base_mod + implicit none + class(mld_d_base_aggregator_type), target, intent(inout) :: ag + type(mld_dml_parms), intent(inout) :: parms + type(psb_dspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_ldspmat_type), intent(inout) :: op_prol,ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_) :: err_act character(len=20) :: name='d_base_aggregator_mat_asb' call psb_erractionsave(err_act) diff --git a/mlprec/mld_d_dec_aggregator_mod.f90 b/mlprec/mld_d_dec_aggregator_mod.f90 index 71dfaf0a..2a7b2336 100644 --- a/mlprec/mld_d_dec_aggregator_mod.f90 +++ b/mlprec/mld_d_dec_aggregator_mod.f90 @@ -80,6 +80,7 @@ module mld_d_dec_aggregator_mod contains procedure, pass(ag) :: bld_tprol => mld_d_dec_aggregator_build_tprol + procedure, pass(ag) :: mat_bld => mld_d_dec_aggregator_mat_bld procedure, pass(ag) :: mat_asb => mld_d_dec_aggregator_mat_asb procedure, pass(ag) :: default => mld_d_dec_aggregator_default procedure, pass(ag) :: set_aggr_type => mld_d_dec_aggregator_set_aggr_type @@ -99,8 +100,8 @@ module mld_d_dec_aggregator_mod class(mld_d_dec_aggregator_type), target, intent(inout) :: ag type(mld_dml_parms), intent(inout) :: parms type(mld_daggr_data), intent(in) :: ag_data - type(psb_dspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_dspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_ldspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info @@ -108,7 +109,7 @@ module mld_d_dec_aggregator_mod end interface interface - subroutine mld_d_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& + subroutine mld_d_dec_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& & op_prol,op_restr,info) import :: mld_d_dec_aggregator_type, psb_desc_type, psb_dspmat_type, psb_dpk_, & & psb_ipk_, psb_lpk_, psb_ldspmat_type, mld_dml_parms @@ -121,9 +122,25 @@ module mld_d_dec_aggregator_mod type(psb_ldspmat_type), intent(inout) :: op_prol type(psb_ldspmat_type), intent(out) :: ac,op_restr integer(psb_ipk_), intent(out) :: info - end subroutine mld_d_dec_aggregator_mat_asb + end subroutine mld_d_dec_aggregator_mat_bld end interface + interface + subroutine mld_d_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac,op_prol,op_restr,info) + import :: mld_d_dec_aggregator_type, psb_desc_type, psb_dspmat_type, psb_dpk_, & + & psb_ipk_, psb_lpk_, psb_ldspmat_type, mld_dml_parms + implicit none + class(mld_d_dec_aggregator_type), target, intent(inout) :: ag + type(mld_dml_parms), intent(inout) :: parms + type(psb_dspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_ldspmat_type), intent(inout) :: op_prol,ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + end subroutine mld_d_dec_aggregator_mat_asb + end interface contains diff --git a/mlprec/mld_d_inner_mod.f90 b/mlprec/mld_d_inner_mod.f90 index 83977dd9..dc654bfd 100644 --- a/mlprec/mld_d_inner_mod.f90 +++ b/mlprec/mld_d_inner_mod.f90 @@ -108,37 +108,9 @@ module mld_d_inner_mod integer(psb_ipk_), intent(out) :: info end subroutine mld_d_map_to_tprol end interface mld_map_to_tprol - - interface mld_lev_mat_asb - subroutine mld_d_lev_aggrmat_asb(p,a,desc_a,ilaggr,nlaggr,op_prol,info) - import :: psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_ldspmat_type - import :: mld_d_onelev_type - implicit none - type(mld_d_onelev_type), intent(inout), target :: p - type(psb_dspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:),nlaggr(:) - type(psb_ldspmat_type), intent(inout) :: op_prol - integer(psb_ipk_), intent(out) :: info - end subroutine mld_d_lev_aggrmat_asb - end interface mld_lev_mat_asb - - interface mld_aggrmat_asb - subroutine mld_daggrmat_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) - import :: psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_ldspmat_type - import :: mld_dml_parms - implicit none - type(psb_dspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) - type(mld_dml_parms), intent(inout) :: parms - type(psb_ldspmat_type), intent(out) :: ac,op_prol,op_restr - integer(psb_ipk_), intent(out) :: info - end subroutine mld_daggrmat_asb - end interface mld_aggrmat_asb abstract interface - subroutine mld_daggrmat_var_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) + subroutine mld_daggrmat_var_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) import :: psb_dspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_ldspmat_type import :: mld_d_onelev_type, mld_dml_parms implicit none @@ -148,11 +120,11 @@ module mld_d_inner_mod type(mld_dml_parms), intent(inout) :: parms type(psb_ldspmat_type), intent(out) :: ac,op_prol,op_restr integer(psb_ipk_), intent(out) :: info - end subroutine mld_daggrmat_var_asb + end subroutine mld_daggrmat_var_bld end interface - procedure(mld_daggrmat_var_asb) :: mld_daggrmat_nosmth_asb, & - & mld_daggrmat_smth_asb, mld_daggrmat_minnrg_asb, & - & mld_daggrmat_biz_asb + procedure(mld_daggrmat_var_bld) :: mld_daggrmat_nosmth_bld, & + & mld_daggrmat_smth_bld, mld_daggrmat_minnrg_bld, & + & mld_daggrmat_biz_bld end module mld_d_inner_mod diff --git a/mlprec/mld_d_onelev_mod.f90 b/mlprec/mld_d_onelev_mod.f90 index 5c8c1c82..a12609ca 100644 --- a/mlprec/mld_d_onelev_mod.f90 +++ b/mlprec/mld_d_onelev_mod.f90 @@ -488,8 +488,8 @@ contains & ilaggr,nlaggr,op_prol,ag_data,info) implicit none class(mld_d_onelev_type), intent(inout), target :: lv - type(psb_dspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_dspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_ldspmat_type), intent(out) :: op_prol type(mld_daggr_data), intent(in) :: ag_data diff --git a/mlprec/mld_d_symdec_aggregator_mod.f90 b/mlprec/mld_d_symdec_aggregator_mod.f90 index be247126..1a471821 100644 --- a/mlprec/mld_d_symdec_aggregator_mod.f90 +++ b/mlprec/mld_d_symdec_aggregator_mod.f90 @@ -70,8 +70,8 @@ module mld_d_symdec_aggregator_mod class(mld_d_symdec_aggregator_type), target, intent(inout) :: ag type(mld_dml_parms), intent(inout) :: parms type(mld_daggr_data), intent(in) :: ag_data - type(psb_dspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_dspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_ldspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info diff --git a/mlprec/mld_s_base_aggregator_mod.f90 b/mlprec/mld_s_base_aggregator_mod.f90 index d79311ec..ebc69a26 100644 --- a/mlprec/mld_s_base_aggregator_mod.f90 +++ b/mlprec/mld_s_base_aggregator_mod.f90 @@ -64,8 +64,10 @@ module mld_s_base_aggregator_mod !! !! bld_tprol - Build a tentative prolongator !! - !! mat_asb - Build the final prolongator/restrictor and the - !! coarse matrix ac + !! mat_bld - Build prolongator/restrictor and coarse matrix ac + !! + !! mat_asb - Convert prolongator/restrictor/coarse matrix + !! and fix their descriptor(s) !! !! update_next - Transfer information to the next level; default is !! to do nothing, i.e. aggregators at different @@ -83,6 +85,7 @@ module mld_s_base_aggregator_mod logical :: do_clean_zeros contains procedure, pass(ag) :: bld_tprol => mld_s_base_aggregator_build_tprol + procedure, pass(ag) :: mat_bld => mld_s_base_aggregator_mat_bld procedure, pass(ag) :: mat_asb => mld_s_base_aggregator_mat_asb procedure, pass(ag) :: bld_map => mld_s_base_aggregator_bld_map procedure, pass(ag) :: update_next => mld_s_base_aggregator_update_next @@ -273,8 +276,8 @@ contains class(mld_s_base_aggregator_type), target, intent(inout) :: ag type(mld_sml_parms), intent(inout) :: parms type(mld_saggr_data), intent(in) :: ag_data - type(psb_sspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_sspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lsspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info @@ -298,7 +301,7 @@ contains end subroutine mld_s_base_aggregator_build_tprol ! - !> Function mat_asb + !> Function mat_bld !! \memberof mld_s_base_aggregator_type !! \brief Build prolongator/restrictor/coarse matrix. !! @@ -316,7 +319,7 @@ contains !! in many cases it is the transpose of the prolongator. !! \param info Return code !! - subroutine mld_s_base_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& + subroutine mld_s_base_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& & op_prol,op_restr,info) use psb_base_mod implicit none @@ -329,6 +332,54 @@ contains type(psb_lsspmat_type), intent(out) :: ac,op_restr integer(psb_ipk_), intent(out) :: info integer(psb_ipk_) :: err_act + character(len=20) :: name='s_base_aggregator_mat_bld' + + call psb_erractionsave(err_act) + + info = psb_err_missing_override_method_ + call psb_errpush(info,name) + goto 9999 + + call psb_erractionrestore(err_act) + return + +9999 call psb_error_handler(err_act) + + return + end subroutine mld_s_base_aggregator_mat_bld + + ! + !> Function mat_asb + !! \memberof mld_s_base_aggregator_type + !! \brief Build prolongator/restrictor/coarse matrix. + !! + !! + !! \param ag The input aggregator object + !! \param parms The auxiliary parameters object + !! \param a The local matrix part + !! \param desc_a The descriptor + !! \param ilaggr Aggregation map + !! \param nlaggr Sizes of ilaggr on all processes + !! \param ac On output the coarse matrix + !! \param op_prol On input, the tentative prolongator operator, on output + !! the final prolongator + !! \param op_restr On output, the restrictor operator; + !! in many cases it is the transpose of the prolongator. + !! \param info Return code + !! + subroutine mld_s_base_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac, op_prol,op_restr,info) + use psb_base_mod + implicit none + class(mld_s_base_aggregator_type), target, intent(inout) :: ag + type(mld_sml_parms), intent(inout) :: parms + type(psb_sspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lsspmat_type), intent(inout) :: op_prol,ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_) :: err_act character(len=20) :: name='s_base_aggregator_mat_asb' call psb_erractionsave(err_act) diff --git a/mlprec/mld_s_dec_aggregator_mod.f90 b/mlprec/mld_s_dec_aggregator_mod.f90 index 0e8380f1..8dd2fef9 100644 --- a/mlprec/mld_s_dec_aggregator_mod.f90 +++ b/mlprec/mld_s_dec_aggregator_mod.f90 @@ -80,6 +80,7 @@ module mld_s_dec_aggregator_mod contains procedure, pass(ag) :: bld_tprol => mld_s_dec_aggregator_build_tprol + procedure, pass(ag) :: mat_bld => mld_s_dec_aggregator_mat_bld procedure, pass(ag) :: mat_asb => mld_s_dec_aggregator_mat_asb procedure, pass(ag) :: default => mld_s_dec_aggregator_default procedure, pass(ag) :: set_aggr_type => mld_s_dec_aggregator_set_aggr_type @@ -99,8 +100,8 @@ module mld_s_dec_aggregator_mod class(mld_s_dec_aggregator_type), target, intent(inout) :: ag type(mld_sml_parms), intent(inout) :: parms type(mld_saggr_data), intent(in) :: ag_data - type(psb_sspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_sspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lsspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info @@ -108,7 +109,7 @@ module mld_s_dec_aggregator_mod end interface interface - subroutine mld_s_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& + subroutine mld_s_dec_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& & op_prol,op_restr,info) import :: mld_s_dec_aggregator_type, psb_desc_type, psb_sspmat_type, psb_spk_, & & psb_ipk_, psb_lpk_, psb_lsspmat_type, mld_sml_parms @@ -121,9 +122,25 @@ module mld_s_dec_aggregator_mod type(psb_lsspmat_type), intent(inout) :: op_prol type(psb_lsspmat_type), intent(out) :: ac,op_restr integer(psb_ipk_), intent(out) :: info - end subroutine mld_s_dec_aggregator_mat_asb + end subroutine mld_s_dec_aggregator_mat_bld end interface + interface + subroutine mld_s_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac,op_prol,op_restr,info) + import :: mld_s_dec_aggregator_type, psb_desc_type, psb_sspmat_type, psb_spk_, & + & psb_ipk_, psb_lpk_, psb_lsspmat_type, mld_sml_parms + implicit none + class(mld_s_dec_aggregator_type), target, intent(inout) :: ag + type(mld_sml_parms), intent(inout) :: parms + type(psb_sspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lsspmat_type), intent(inout) :: op_prol,ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + end subroutine mld_s_dec_aggregator_mat_asb + end interface contains diff --git a/mlprec/mld_s_inner_mod.f90 b/mlprec/mld_s_inner_mod.f90 index 29e01457..91d16768 100644 --- a/mlprec/mld_s_inner_mod.f90 +++ b/mlprec/mld_s_inner_mod.f90 @@ -108,37 +108,9 @@ module mld_s_inner_mod integer(psb_ipk_), intent(out) :: info end subroutine mld_s_map_to_tprol end interface mld_map_to_tprol - - interface mld_lev_mat_asb - subroutine mld_s_lev_aggrmat_asb(p,a,desc_a,ilaggr,nlaggr,op_prol,info) - import :: psb_sspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lsspmat_type - import :: mld_s_onelev_type - implicit none - type(mld_s_onelev_type), intent(inout), target :: p - type(psb_sspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:),nlaggr(:) - type(psb_lsspmat_type), intent(inout) :: op_prol - integer(psb_ipk_), intent(out) :: info - end subroutine mld_s_lev_aggrmat_asb - end interface mld_lev_mat_asb - - interface mld_aggrmat_asb - subroutine mld_saggrmat_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) - import :: psb_sspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lsspmat_type - import :: mld_sml_parms - implicit none - type(psb_sspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) - type(mld_sml_parms), intent(inout) :: parms - type(psb_lsspmat_type), intent(out) :: ac,op_prol,op_restr - integer(psb_ipk_), intent(out) :: info - end subroutine mld_saggrmat_asb - end interface mld_aggrmat_asb abstract interface - subroutine mld_saggrmat_var_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) + subroutine mld_saggrmat_var_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) import :: psb_sspmat_type, psb_desc_type, psb_spk_, psb_ipk_, psb_lpk_, psb_lsspmat_type import :: mld_s_onelev_type, mld_sml_parms implicit none @@ -148,11 +120,11 @@ module mld_s_inner_mod type(mld_sml_parms), intent(inout) :: parms type(psb_lsspmat_type), intent(out) :: ac,op_prol,op_restr integer(psb_ipk_), intent(out) :: info - end subroutine mld_saggrmat_var_asb + end subroutine mld_saggrmat_var_bld end interface - procedure(mld_saggrmat_var_asb) :: mld_saggrmat_nosmth_asb, & - & mld_saggrmat_smth_asb, mld_saggrmat_minnrg_asb, & - & mld_saggrmat_biz_asb + procedure(mld_saggrmat_var_bld) :: mld_saggrmat_nosmth_bld, & + & mld_saggrmat_smth_bld, mld_saggrmat_minnrg_bld, & + & mld_saggrmat_biz_bld end module mld_s_inner_mod diff --git a/mlprec/mld_s_onelev_mod.f90 b/mlprec/mld_s_onelev_mod.f90 index d743b366..52d5376c 100644 --- a/mlprec/mld_s_onelev_mod.f90 +++ b/mlprec/mld_s_onelev_mod.f90 @@ -488,8 +488,8 @@ contains & ilaggr,nlaggr,op_prol,ag_data,info) implicit none class(mld_s_onelev_type), intent(inout), target :: lv - type(psb_sspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_sspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lsspmat_type), intent(out) :: op_prol type(mld_saggr_data), intent(in) :: ag_data diff --git a/mlprec/mld_s_symdec_aggregator_mod.f90 b/mlprec/mld_s_symdec_aggregator_mod.f90 index 026abbe6..e37c5ab2 100644 --- a/mlprec/mld_s_symdec_aggregator_mod.f90 +++ b/mlprec/mld_s_symdec_aggregator_mod.f90 @@ -70,8 +70,8 @@ module mld_s_symdec_aggregator_mod class(mld_s_symdec_aggregator_type), target, intent(inout) :: ag type(mld_sml_parms), intent(inout) :: parms type(mld_saggr_data), intent(in) :: ag_data - type(psb_sspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_sspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lsspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info diff --git a/mlprec/mld_z_base_aggregator_mod.f90 b/mlprec/mld_z_base_aggregator_mod.f90 index 00937d5e..2b8d1aca 100644 --- a/mlprec/mld_z_base_aggregator_mod.f90 +++ b/mlprec/mld_z_base_aggregator_mod.f90 @@ -64,8 +64,10 @@ module mld_z_base_aggregator_mod !! !! bld_tprol - Build a tentative prolongator !! - !! mat_asb - Build the final prolongator/restrictor and the - !! coarse matrix ac + !! mat_bld - Build prolongator/restrictor and coarse matrix ac + !! + !! mat_asb - Convert prolongator/restrictor/coarse matrix + !! and fix their descriptor(s) !! !! update_next - Transfer information to the next level; default is !! to do nothing, i.e. aggregators at different @@ -83,6 +85,7 @@ module mld_z_base_aggregator_mod logical :: do_clean_zeros contains procedure, pass(ag) :: bld_tprol => mld_z_base_aggregator_build_tprol + procedure, pass(ag) :: mat_bld => mld_z_base_aggregator_mat_bld procedure, pass(ag) :: mat_asb => mld_z_base_aggregator_mat_asb procedure, pass(ag) :: bld_map => mld_z_base_aggregator_bld_map procedure, pass(ag) :: update_next => mld_z_base_aggregator_update_next @@ -273,8 +276,8 @@ contains class(mld_z_base_aggregator_type), target, intent(inout) :: ag type(mld_dml_parms), intent(inout) :: parms type(mld_daggr_data), intent(in) :: ag_data - type(psb_zspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_zspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lzspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info @@ -298,7 +301,7 @@ contains end subroutine mld_z_base_aggregator_build_tprol ! - !> Function mat_asb + !> Function mat_bld !! \memberof mld_z_base_aggregator_type !! \brief Build prolongator/restrictor/coarse matrix. !! @@ -316,7 +319,7 @@ contains !! in many cases it is the transpose of the prolongator. !! \param info Return code !! - subroutine mld_z_base_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& + subroutine mld_z_base_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& & op_prol,op_restr,info) use psb_base_mod implicit none @@ -329,6 +332,54 @@ contains type(psb_lzspmat_type), intent(out) :: ac,op_restr integer(psb_ipk_), intent(out) :: info integer(psb_ipk_) :: err_act + character(len=20) :: name='z_base_aggregator_mat_bld' + + call psb_erractionsave(err_act) + + info = psb_err_missing_override_method_ + call psb_errpush(info,name) + goto 9999 + + call psb_erractionrestore(err_act) + return + +9999 call psb_error_handler(err_act) + + return + end subroutine mld_z_base_aggregator_mat_bld + + ! + !> Function mat_asb + !! \memberof mld_z_base_aggregator_type + !! \brief Build prolongator/restrictor/coarse matrix. + !! + !! + !! \param ag The input aggregator object + !! \param parms The auxiliary parameters object + !! \param a The local matrix part + !! \param desc_a The descriptor + !! \param ilaggr Aggregation map + !! \param nlaggr Sizes of ilaggr on all processes + !! \param ac On output the coarse matrix + !! \param op_prol On input, the tentative prolongator operator, on output + !! the final prolongator + !! \param op_restr On output, the restrictor operator; + !! in many cases it is the transpose of the prolongator. + !! \param info Return code + !! + subroutine mld_z_base_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac, op_prol,op_restr,info) + use psb_base_mod + implicit none + class(mld_z_base_aggregator_type), target, intent(inout) :: ag + type(mld_dml_parms), intent(inout) :: parms + type(psb_zspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lzspmat_type), intent(inout) :: op_prol,ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_) :: err_act character(len=20) :: name='z_base_aggregator_mat_asb' call psb_erractionsave(err_act) diff --git a/mlprec/mld_z_dec_aggregator_mod.f90 b/mlprec/mld_z_dec_aggregator_mod.f90 index 307716b0..18db135a 100644 --- a/mlprec/mld_z_dec_aggregator_mod.f90 +++ b/mlprec/mld_z_dec_aggregator_mod.f90 @@ -80,6 +80,7 @@ module mld_z_dec_aggregator_mod contains procedure, pass(ag) :: bld_tprol => mld_z_dec_aggregator_build_tprol + procedure, pass(ag) :: mat_bld => mld_z_dec_aggregator_mat_bld procedure, pass(ag) :: mat_asb => mld_z_dec_aggregator_mat_asb procedure, pass(ag) :: default => mld_z_dec_aggregator_default procedure, pass(ag) :: set_aggr_type => mld_z_dec_aggregator_set_aggr_type @@ -99,8 +100,8 @@ module mld_z_dec_aggregator_mod class(mld_z_dec_aggregator_type), target, intent(inout) :: ag type(mld_dml_parms), intent(inout) :: parms type(mld_daggr_data), intent(in) :: ag_data - type(psb_zspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_zspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lzspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info @@ -108,7 +109,7 @@ module mld_z_dec_aggregator_mod end interface interface - subroutine mld_z_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& + subroutine mld_z_dec_aggregator_mat_bld(ag,parms,a,desc_a,ilaggr,nlaggr,ac,& & op_prol,op_restr,info) import :: mld_z_dec_aggregator_type, psb_desc_type, psb_zspmat_type, psb_dpk_, & & psb_ipk_, psb_lpk_, psb_lzspmat_type, mld_dml_parms @@ -121,9 +122,25 @@ module mld_z_dec_aggregator_mod type(psb_lzspmat_type), intent(inout) :: op_prol type(psb_lzspmat_type), intent(out) :: ac,op_restr integer(psb_ipk_), intent(out) :: info - end subroutine mld_z_dec_aggregator_mat_asb + end subroutine mld_z_dec_aggregator_mat_bld end interface + interface + subroutine mld_z_dec_aggregator_mat_asb(ag,parms,a,desc_a,ilaggr,nlaggr,& + & ac,desc_ac,op_prol,op_restr,info) + import :: mld_z_dec_aggregator_type, psb_desc_type, psb_zspmat_type, psb_dpk_, & + & psb_ipk_, psb_lpk_, psb_lzspmat_type, mld_dml_parms + implicit none + class(mld_z_dec_aggregator_type), target, intent(inout) :: ag + type(mld_dml_parms), intent(inout) :: parms + type(psb_zspmat_type), intent(in) :: a + type(psb_desc_type), intent(in) :: desc_a + integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) + type(psb_lzspmat_type), intent(inout) :: op_prol,ac,op_restr + type(psb_desc_type), intent(inout) :: desc_ac + integer(psb_ipk_), intent(out) :: info + end subroutine mld_z_dec_aggregator_mat_asb + end interface contains diff --git a/mlprec/mld_z_inner_mod.f90 b/mlprec/mld_z_inner_mod.f90 index 3e46b6f0..b890049f 100644 --- a/mlprec/mld_z_inner_mod.f90 +++ b/mlprec/mld_z_inner_mod.f90 @@ -108,37 +108,9 @@ module mld_z_inner_mod integer(psb_ipk_), intent(out) :: info end subroutine mld_z_map_to_tprol end interface mld_map_to_tprol - - interface mld_lev_mat_asb - subroutine mld_z_lev_aggrmat_asb(p,a,desc_a,ilaggr,nlaggr,op_prol,info) - import :: psb_zspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_lzspmat_type - import :: mld_z_onelev_type - implicit none - type(mld_z_onelev_type), intent(inout), target :: p - type(psb_zspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:),nlaggr(:) - type(psb_lzspmat_type), intent(inout) :: op_prol - integer(psb_ipk_), intent(out) :: info - end subroutine mld_z_lev_aggrmat_asb - end interface mld_lev_mat_asb - - interface mld_aggrmat_asb - subroutine mld_zaggrmat_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) - import :: psb_zspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_lzspmat_type - import :: mld_dml_parms - implicit none - type(psb_zspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a - integer(psb_lpk_), intent(inout) :: ilaggr(:), nlaggr(:) - type(mld_dml_parms), intent(inout) :: parms - type(psb_lzspmat_type), intent(out) :: ac,op_prol,op_restr - integer(psb_ipk_), intent(out) :: info - end subroutine mld_zaggrmat_asb - end interface mld_aggrmat_asb abstract interface - subroutine mld_zaggrmat_var_asb(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) + subroutine mld_zaggrmat_var_bld(a,desc_a,ilaggr,nlaggr,parms,ac,op_prol,op_restr,info) import :: psb_zspmat_type, psb_desc_type, psb_dpk_, psb_ipk_, psb_lpk_, psb_lzspmat_type import :: mld_z_onelev_type, mld_dml_parms implicit none @@ -148,11 +120,11 @@ module mld_z_inner_mod type(mld_dml_parms), intent(inout) :: parms type(psb_lzspmat_type), intent(out) :: ac,op_prol,op_restr integer(psb_ipk_), intent(out) :: info - end subroutine mld_zaggrmat_var_asb + end subroutine mld_zaggrmat_var_bld end interface - procedure(mld_zaggrmat_var_asb) :: mld_zaggrmat_nosmth_asb, & - & mld_zaggrmat_smth_asb, mld_zaggrmat_minnrg_asb, & - & mld_zaggrmat_biz_asb + procedure(mld_zaggrmat_var_bld) :: mld_zaggrmat_nosmth_bld, & + & mld_zaggrmat_smth_bld, mld_zaggrmat_minnrg_bld, & + & mld_zaggrmat_biz_bld end module mld_z_inner_mod diff --git a/mlprec/mld_z_onelev_mod.f90 b/mlprec/mld_z_onelev_mod.f90 index a560ed5f..4b8ec326 100644 --- a/mlprec/mld_z_onelev_mod.f90 +++ b/mlprec/mld_z_onelev_mod.f90 @@ -488,8 +488,8 @@ contains & ilaggr,nlaggr,op_prol,ag_data,info) implicit none class(mld_z_onelev_type), intent(inout), target :: lv - type(psb_zspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_zspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lzspmat_type), intent(out) :: op_prol type(mld_daggr_data), intent(in) :: ag_data diff --git a/mlprec/mld_z_symdec_aggregator_mod.f90 b/mlprec/mld_z_symdec_aggregator_mod.f90 index 9fa7be0b..495ad0e0 100644 --- a/mlprec/mld_z_symdec_aggregator_mod.f90 +++ b/mlprec/mld_z_symdec_aggregator_mod.f90 @@ -70,8 +70,8 @@ module mld_z_symdec_aggregator_mod class(mld_z_symdec_aggregator_type), target, intent(inout) :: ag type(mld_dml_parms), intent(inout) :: parms type(mld_daggr_data), intent(in) :: ag_data - type(psb_zspmat_type), intent(in) :: a - type(psb_desc_type), intent(in) :: desc_a + type(psb_zspmat_type), intent(inout) :: a + type(psb_desc_type), intent(inout) :: desc_a integer(psb_lpk_), allocatable, intent(out) :: ilaggr(:),nlaggr(:) type(psb_lzspmat_type), intent(out) :: op_prol integer(psb_ipk_), intent(out) :: info