next up previous contents
Next: Subroutine mld_precbld Up: User Interface Previous: Subroutine mld_precinit   Contents


Subroutine mld_precset

mld_precset(p,what,val,info)

This routine sets the parameters defining the preconditioner. More precisely, the parameter identified by what is assigned the value contained in val.

Arguments

p type(mld_xprec_type), intent(inout).
  The preconditioner data structure. Note that x must be chosen according to the real/complex, single/double precision version of MLD2P4 under use.
what integer, intent(in).
  The number identifying the parameter to be set. A mnemonic constant has been associated to each of these numbers, as reported in Tables 2-5.
val integer or character(len=*) or real(psb_spk_) or real(psb_dpk_), intent(in).
  The value of the parameter to be set. The list of allowed values and the corresponding data types is given in Tables 2-5. When the value is of type character(len=*), it is also treated as case insensitive.
info integer, intent(out).
  Error code. If no error, 0 is returned. See Section 7 for details.


A variety of (one-level and multi-level) preconditioners can be obtained by a suitable setting of the preconditioner parameters. These parameters can be logically divided into four groups, i.e. parameters defining

  1. the type of multi-level preconditioner;
  2. the one-level preconditioner used as smoother;
  3. the aggregation algorithm;
  4. the coarse-space correction at the coarsest level.
A list of the parameters that can be set, along with their allowed and default values, is given in Tables 2-5. For a detailed description of the meaning of the parameters, please refer to Section 4.


Table 2: Parameters defining the type of multi-level preconditioner.
what DATA TYPE val DEFAULT COMMENTS
mld_ml_type_ character(len=*) 'ADD' 'MULT' 'MULT' Basic multi-level framework: additive or multiplicative among the levels (always additive inside a level).
mld_smoother_type_ character(len=*) 'DIAG' 'BJAC' 'AS' 'AS' Basic one-level preconditioner (i.e. smoother): diagonal, block Jacobi, AS.
mld_smoother_pos_ character(len=*) 'PRE' 'POST' 'TWOSIDE' 'POST' ``Position'' of the smoother: pre-smoother, post-smoother, pre- and post-smoother.



Table 3: Parameters defining the one-level preconditioner used as smoother.
what DATA TYPE val DEFAULT COMMENTS
mld_sub_ovr_ integer any int. num. $\ge 0$ 1 Number of overlap layers.
mld_sub_restr_ character(len=*) 'HALO' 'NONE' 'HALO' Type of restriction operator: 'HALO' for taking into account the overlap, 'NONE' for neglecting it.
mld_sub_prol_ character(len=*) 'SUM' 'NONE' 'NONE' Type of prolongation operator: 'SUM' for adding the contributions from the overlap, 'NONE' for neglecting them.
mld_sub_solve_ character(len=*) 'ILU' 'MILU' 'ILUT' 'UMF' 'SLU' 'ILU' Local solver: ILU($p$), MILU($p$), ILU($p,t$), LU from UMFPACK, LU from SuperLU (plus triangular solve).
mld_sub_fillin_ integer Any int. num. $\ge 0$ 0 Fill-in level $p$ of the incomplete LU factorizations.
mld_sub_iluthrs_ real(kind_parameter) Any real num. $\ge 0$ 0 Drop tolerance $t$ in the ILU($p,t$) factorization.
mld_sub_ren_ character(len=*) 'RENUM_NONE' 'RENUM_GLOBAL' 'RENUM_NONE' Row and column reordering of the local submatrices: no reordering, reordering according to the global numbering of the rows and columns of the whole matrix.



Table 4: Parameters defining the aggregation algorithm.
what DATA TYPE val DEFAULT COMMENTS
mld_aggr_alg_ character(len=*) 'DEC' 'DEC' Aggregation algorithm. Currently, only the decoupled aggregation is available.
mld_aggr_kind_ character(len=*) 'SMOOTH' 'RAW' 'SMOOTH' Type of aggregation: smoothed, raw (i.e. using the tentative prolongator).
mld_aggr_thresh_ real(kind_parameter) Any real num. $\in [0, 1]$ 0 Threshold $\theta$ in the aggregation algorithm.
mld_aggr_omega_alg_ character(len=*) 'EIG_EST' 'USER_CHOICE' 'EIG_EST' How the damping parameter $\omega$ in the smoothed aggregation should be computed: either via an estimate of the spectral radius of $D^{-1}A$, or explicily specified by the user.
mld_aggr_eig_ character(len=*) 'A_NORMI' 'A_NORMI' How to estimate the spectral radius of $D^{-1}A$. Currently only the infinity norm estimate is available.
mld_aggr_omega_val_ real(kind_parameter) Any nonnegative real num. $4/(3\rho(D^{-1}A))$ Damping parameter $\omega$ in the smoothed aggregation algorithm. It must be set by the user if USER_CHOICE was specified for mld_aggr_omega_alg_, otherwise it is computed by the library, using the selected estimate of the spectral radius $\rho(D^{-1}A)$ of $D^{-1}A$.



Table 5: Parameters defining the coarse-space correction at the coarsest level.
what DATA TYPE val DEFAULT COMMENTS
mld_coarse_mat_ character(len=*) 'DISTR' 'REPL' 'DISTR' Coarsest matrix: distributed among the processors or replicated on each of them.
mld_coarse_solve_ character(len=*) 'BJAC' 'UMF' 'SLU' 'SLUDIST' 'BJAC' Solver used at the coarsest level: block Jacobi, sequential LU from UMFPACK, sequential LU from SuperLU, distributed LU from SuperLU_Dist. 'SLUDIST' requires the coarsest matrix to be distributed, while 'UMF' and 'SLU' require it to be replicated.
mld_coarse_subsolve_ character(len=*) 'ILU' 'MILU' 'ILUT' 'UMF' 'SLU' See note Solver for the diagonal blocks of the coarse matrix, in case the block Jacobi solver is chosen as coarsest-level solver: ILU($p$), MILU($p$), ILU($p,t$), LU from UMFPACK, LU from SuperLU, plus triangular solve.
mld_coarse_sweeps_ integer Any int. num. $> 0$ 4 Number of Block-Jacobi sweeps when 'BJAC' is used as coarsest-level solver.
mld_coarse_fillin_ integer Any int. num. $\ge 0$ 0 Fill-in level $p$ of the incomplete LU factorizations.
mld_coarse_iluthrs_ real(kind_parameter) Any real. num. $\ge 0$ 0 Drop tolerance $t$ in the ILU($p,t$) factorization.
Note Defaults for mld_coarse_subsolve_ are chosen as
Single precision version: 'SLU' if installed, 'ILU' otherwise
Double precision version: 'UMF' if installed, else 'SLU' if installed, 'ILU' otherwise



next up previous contents
Next: Subroutine mld_precbld Up: User Interface Previous: Subroutine mld_precinit   Contents