<!DOCTYPE html PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN" "http://www.w3.org/TR/html4/loose.dtd"> <html > <head><title>Examples</title> <meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1"> <meta name="generator" content="TeX4ht (https://tug.org/tex4ht/)"> <meta name="originator" content="TeX4ht (https://tug.org/tex4ht/)"> <!-- html,3 --> <meta name="src" content="userhtml.tex"> <link rel="stylesheet" type="text/css" href="userhtml.css"> </head><body > <!--l. 114--><div class="crosslinks"><p class="noindent"><span class="cmr-12">[</span><a href="userhtmlsu7.html" ><span class="cmr-12">next</span></a><span class="cmr-12">] [</span><a href="#tailuserhtmlsu6.html"><span class="cmr-12">tail</span></a><span class="cmr-12">] [</span><a href="userhtmlse4.html#userhtmlsu6.html" ><span class="cmr-12">up</span></a><span class="cmr-12">] </span></p></div> <h4 class="subsectionHead"><span class="titlemark"><span class="cmr-12">4.1 </span></span> <a id="x15-140004.1"></a><span class="cmr-12">Examples</span></h4> <!--l. 116--><p class="noindent" ><span class="cmr-12">The code reported in Figure</span><span class="cmr-12"> </span><a href="#x15-14001r1"><span class="cmr-12">1</span><!--tex4ht:ref: fig:ex1 --></a> <span class="cmr-12">shows how to set and apply the default multilevel</span> <span class="cmr-12">preconditioner available in the real double precision version of AMG4PSBLAS</span> <span class="cmr-12">(see Table</span><span class="cmr-12"> </span><a href="userhtmlse4.html#x14-13015r1"><span class="cmr-12">1</span><!--tex4ht:ref: tab:precinit --></a><span class="cmr-12">). This preconditioner is chosen by simply specifying </span><span class="lstinline"></span><span class="cmtt-12">’</span><span class="cmtt-12">ML</span><span class="cmtt-12">’</span> <span class="cmr-12">as the</span> <span class="cmr-12">second argument of </span><span class="lstinline"></span><span class="cmtt-12">P</span><span class="cmtt-12">%</span><span class="cmtt-12">init</span> <span class="cmr-12">(a call to </span><span class="lstinline"></span><span class="cmtt-12">P</span><span class="cmtt-12">%</span><span class="cmtt-12">set</span> <span class="cmr-12">is not needed) and is applied</span> <span class="cmr-12">with the CG solver provided by PSBLAS (the matrix of the system to be</span> <span class="cmr-12">solved is assumed to be positive definite). As previously observed, the modules</span> <span class="lstinline"></span><span class="cmtt-12">psb_base_mod</span><span class="cmr-12">, </span><span class="lstinline"></span><span class="cmtt-12">amg_prec_mod</span> <span class="cmr-12">and </span><span class="lstinline"></span><span class="cmtt-12">psb_krylov_mod</span> <span class="cmr-12">must be used by the example</span> <span class="cmr-12">program.</span> <!--l. 126--><p class="indent" > <span class="cmr-12">The part of the code dealing with reading and assembling the sparse matrix and the</span> <span class="cmr-12">right-hand side vector and the deallocation of the relevant data structures, performed</span> <span class="cmr-12">through the PSBLAS routines for sparse matrix and vector management,</span> <span class="cmr-12">is not reported here for the sake of conciseness. The complete code can be</span> <span class="cmr-12">found in the example program file </span><span class="obeylines-h"><span class="verb"><span class="cmtt-12">amg_dexample_ml.f90</span></span></span><span class="cmr-12">, in the directory</span> <span class="obeylines-h"><span class="verb"><span class="cmtt-12">samples/simple/file</span></span></span><span class="obeylines-h"><span class="verb"><span class="cmtt-12">read</span></span></span> <span class="cmr-12">of the AMG4PSBLAS implementation (see Section</span><span class="cmr-12"> </span><a href="userhtmlsu5.html#x13-120003.5"><span class="cmr-12">3.5</span><!--tex4ht:ref: sec:ex_and_test --></a><span class="cmr-12">). A</span> <span class="cmr-12">sample test problem along with the relevant input data is available in</span> <span class="obeylines-h"><span class="verb"><span class="cmtt-12">samples/simple/fileread/runs</span></span></span><span class="cmr-12">. For details on the use of the PSBLAS routines, see</span> <span class="cmr-12">the PSBLAS User’s Guide</span><span class="cmr-12"> </span><span class="cite"><span class="cmr-12">[</span><a href="userhtmlli5.html#XPSBLASGUIDE"><span class="cmr-12">20</span></a><span class="cmr-12">]</span></span><span class="cmr-12">.</span> <!--l. 138--><p class="indent" > <span class="cmr-12">The setup and application of the default multilevel preconditioner for the real single</span> <span class="cmr-12">precision and the complex, single and double precision, versions are obtained</span> <span class="cmr-12">with straightforward modifications of the previous example (see Section</span><span class="cmr-12"> </span><a href="userhtmlse5.html#x17-160005"><span class="cmr-12">5</span><!--tex4ht:ref: sec:userinterface --></a> <span class="cmr-12">for</span> <span class="cmr-12">details). If these versions are installed, the corresponding codes are available in</span> <span class="obeylines-h"><span class="verb"><span class="cmtt-12">samples/simple/file</span></span></span><span class="obeylines-h"><span class="verb"><span class="cmtt-12">read</span></span></span><span class="cmr-12">.</span> <!--l. 144--><p class="indent" > <a id="x15-14001r1"></a><hr class="float"><div class="float" > <div class="center" > <!--l. 145--><p class="noindent" > <div class="minipage"><pre class="verbatim" id="verbatim-6">   use psb_base_mod   use amg_prec_mod   use psb_krylov_mod ... ... ! ! sparse matrix   type(psb_dspmat_type) :: A ! sparse matrix descriptor   type(psb_desc_type)   :: desc_A ! preconditioner   type(amg_dprec_type)  :: P ! right-hand side and solution vectors   type(psb_d_vect_type) :: b, x ... ... ! ! initialize the parallel environment   call psb_init(ctxt)   call psb_info(ctxt,iam,np) ... ... ! ! read and assemble the spd matrix A and the right-hand side b ! using PSBLAS routines for sparse matrix / vector management ... ... ! ! initialize the default multilevel preconditioner, i.e. V-cycle ! with basic smoothed aggregation, 1 hybrid forward/backward ! GS sweep as pre/post-smoother and UMFPACK as coarsest-level ! solver   call P%init(ctxt,’ML’,info) ! ! build the preconditioner   call P%hierarchy_build(A,desc_A,info)   call P%smoothers_build(A,desc_A,info) ! ! set the solver parameters and the initial guess   ... ... ! ! solve Ax=b with preconditioned FCG   call psb_krylov(’FCG’,A,P,b,x,tol,desc_A,info)   ... ... ! ! deallocate the preconditioner   call P%free(info) ! ! deallocate other data structures   ... ... ! ! exit the parallel environment   call psb_exit(ctxt)   stop </pre> <!--l. 255--><p class="nopar" > </div> <br /> <div class="caption" ><span class="id">Listing 1: </span><span class="content">setup and application of the default multilevel preconditioner (example 1). </span></div><!--tex4ht:label?: x15-14001r1 --> </div> </div><hr class="endfloat" /> <!--l. 264--><p class="indent" > <span class="cmr-12">Different versions of the multilevel preconditioner can be obtained by changing the</span> <span class="cmr-12">default values of the preconditioner parameters. The code reported in Figure</span><span class="cmr-12"> </span><a href="#x15-14002r2"><span class="cmr-12">2</span><!--tex4ht:ref: fig:ex2 --></a> <span class="cmr-12">shows</span> <span class="cmr-12">how to set a V-cycle preconditioner which applies 1 block-Jacobi sweep as pre-</span> <span class="cmr-12">and post-smoother, and solves the coarsest-level system with 8 block-Jacobi</span> <span class="cmr-12">sweeps. Note that the ILU(0) factorization (plus triangular solve) is used as</span> <span class="cmr-12">local solver for the block-Jacobi sweeps, since this is the default associated</span> <span class="cmr-12">with block-Jacobi and set by</span><span class="cmr-12"> </span><span class="lstinline"></span><span class="cmtt-12">P</span><span class="cmtt-12">%</span><span class="cmtt-12">init</span><span class="cmr-12">. Furthermore, specifying block-Jacobi as</span> <span class="cmr-12">coarsest-level solver implies that the coarsest-level matrix is distributed among</span> <span class="cmr-12">the processes. Figure</span><span class="cmr-12"> </span><a href="#x15-14003r3"><span class="cmr-12">3</span><!--tex4ht:ref: fig:ex3 --></a> <span class="cmr-12">shows how to set a W-cycle preconditioner using the</span> <span class="cmr-12">Coarsening based on Compatible Weighted Matching, aggregates of size at</span> <span class="cmr-12">most 8 and smoothed prolongators. It applies 2 hybrid Gauss-Seidel sweeps as</span> <span class="cmr-12">pre- and post-smoother, and solves the coarsest-level system with the parallel</span> <span class="cmr-12">flexible Conjugate Gradient method (KRM) coupled with the block-Jacobi</span> <span class="cmr-12">preconditioner having ILU(0) on the blocks. Default parameters are used for stopping</span> <span class="cmr-12">criterion of the coarsest solver. Note that, also in this case, specifying KRM as</span> <span class="cmr-12">coarsest-level solver implies that the coarsest-level matrix is distributed among the</span> <span class="cmr-12">processes.</span> <!--l. 291--><p class="indent" > <span class="cmr-12">The code fragments shown in Figures</span><span class="cmr-12"> </span><a href="#x15-14002r2"><span class="cmr-12">2</span><!--tex4ht:ref: fig:ex2 --></a> <span class="cmr-12">and </span><a href="#x15-14003r3"><span class="cmr-12">3</span><!--tex4ht:ref: fig:ex3 --></a> <span class="cmr-12">are included in the example program</span> <span class="cmr-12">file </span><span class="obeylines-h"><span class="verb"><span class="cmtt-12">amg_dexample_ml.f90</span></span></span> <span class="cmr-12">too.</span> <!--l. 294--><p class="indent" > <span class="cmr-12">Finally, Figure</span><span class="cmr-12"> </span><a href="#x15-14004r4"><span class="cmr-12">4</span><!--tex4ht:ref: fig:ex4 --></a> <span class="cmr-12">shows the setup of a one-level additive Schwarz preconditioner,</span> <span class="cmr-12">i.e., RAS with overlap 2. Note also that a Krylov method different from CG</span> <span class="cmr-12">must be used to solve the preconditioned system, since the preconditione in</span> <span class="cmr-12">nonsymmetric. The corresponding example program is available in the file</span> <span class="obeylines-h"><span class="verb"><span class="cmtt-12">amg_dexample_1lev.f90</span></span></span><span class="cmr-12">.</span> <!--l. 301--><p class="indent" > <span class="cmr-12">For all the previous preconditioners, example programs where the sparse matrix</span> <span class="cmr-12">and the right-hand side are generated by discretizing a PDE with Dirichlet</span> <span class="cmr-12">boundary conditions are also available in the directory </span><span class="obeylines-h"><span class="verb"><span class="cmtt-12">samples/simple/pdegen</span></span></span><span class="cmr-12">.</span> <!--l. 304--><p class="indent" > <a id="x15-14002r2"></a><hr class="float"><div class="float" > <div class="center" > <!--l. 318--><p class="noindent" > <div class="minipage"><pre class="verbatim" id="verbatim-7"> ... ... ! build a V-cycle preconditioner with 1 block-Jacobi sweep (with ! ILU(0) on the blocks) as pre- and post-smoother, and 8  block-Jacobi ! sweeps (with ILU(0) on the blocks) as coarsest-level solver   call P%init(ctxt,’ML’,info)   call P%set(’SMOOTHER_TYPE’,’BJAC’,info)   call P%set(’COARSE_SOLVE’,’BJAC’,info)   call P%set(’COARSE_SWEEPS’,8,info)   call P%hierarchy_build(A,desc_A,info)   call P%smoothers_build(A,desc_A,info) ... ... </pre> <!--l. 333--><p class="nopar" ></div></div> <br /><div class="caption" ><span class="id">Listing 2: </span><span class="content">setup of a multilevel preconditioner based on the default decoupled coarsening</span></div><!--tex4ht:label?: x15-14002r2 --> </div><hr class="endfloat" /> <!--l. 340--><p class="indent" > <a id="x15-14003r3"></a><hr class="float"><div class="float" > <div class="center" > <!--l. 362--><p class="noindent" > <div class="minipage"><pre class="verbatim" id="verbatim-8"> ... ... ! build a W-cycle preconditioner with 2 hybrid Gauss-Seidel sweeps ! as pre- and post-smoother, a distributed coarsest ! matrix, and MUMPS as coarsest-level solver   call P%init(ctxt,’ML’,info)   call P%set(’PAR_AGGR_ALG’,’COUPLED’,info)   call P%set(’AGGR_TYPE’,’MATCHBOXP’,info)   call P%set(’AGGR_SIZE’,8,info)   call P%set(’ML_CYCLE’,’WCYCLE’,info)   call P%set(’SMOOTHER_TYPE’,’FBGS’,info)   call P%set(’SMOOTHER_SWEEPS’,2,info)   call P%set(’COARSE_SOLVE’,’KRM’,info)   call P%set(’COARSE_MAT’,’DIST’,info)   call P%set(’KRM_METHOD’,’FCG’,info)   call P%hierarchy_build(A,desc_A,info)   call P%smoothers_build(A,desc_A,info) ... ... </pre> <!--l. 383--><p class="nopar" ></div></div> <br /> <div class="caption" ><span class="id">Listing 3: </span><span class="content">setup of a multilevel preconditioner based on the coupled coarsening using weighted matching</span></div><!--tex4ht:label?: x15-14003r3 --> </div><hr class="endfloat" /> <!--l. 390--><p class="indent" > <a id="x15-14004r4"></a><hr class="float"><div class="float" > <div class="center" > <!--l. 402--><p class="noindent" > <div class="minipage"><pre class="verbatim" id="verbatim-9"> ... ... ! set RAS with overlap 2 and ILU(0) on the local blocks   call P%init(ctxt,’AS’,info)   call P%set(’SUB_OVR’,2,info)   call P%bld(A,desc_A,info) ... ... ! solve Ax=b with preconditioned BiCGSTAB   call psb_krylov(’BICGSTAB’,A,P,b,x,tol,desc_A,info) </pre> <!--l. 414--><p class="nopar" ></div></div> <br /> <div class="caption" ><span class="id">Listing 4: </span><span class="content">setup of a one-level Schwarz preconditioner.</span></div><!--tex4ht:label?: x15-14004r4 --> </div><hr class="endfloat" /> <!--l. 424--><div class="crosslinks"><p class="noindent"><span class="cmr-12">[</span><a href="userhtmlsu7.html" ><span class="cmr-12">next</span></a><span class="cmr-12">] [</span><a href="userhtmlsu6.html" ><span class="cmr-12">front</span></a><span class="cmr-12">] [</span><a href="userhtmlse4.html#userhtmlsu6.html" ><span class="cmr-12">up</span></a><span class="cmr-12">] </span></p></div> <!--l. 424--><p class="indent" > <a id="tailuserhtmlsu6.html"></a> </body></html>