diff --git a/krylov/Makefile b/krylov/Makefile index bebb82b3..5f0cbf7c 100644 --- a/krylov/Makefile +++ b/krylov/Makefile @@ -14,7 +14,8 @@ PSBKRYLDIR=$(PSBLASDIR)/krylov METHDOBJS=psb_dcgstab.o psb_dcg.o psb_dcgs.o \ psb_dbicg.o psb_dcgstabl.o psb_drgmres.o\ - psb_zcgstab.o psb_zcgs.o psb_zrgmres.o + psb_zcgstab.o psb_zcg.o psb_zcgs.o \ + psb_zbicg.o psb_zcgstabl.o psb_zrgmres.o LIBMOD=psb_krylov_mod$(.mod) psb_prec_mod$(.mod) MODOBJS=$(LIBMOD:$(.mod)=.o) @@ -29,7 +30,7 @@ lib: $(OBJS) /bin/cp -p $(HERE)/$(LIBNAME) $(LIBDIR) /bin/cp -p $(LIBMOD) $(LIBDIR) -$(METHDOBJS): psb_prec_mod.o +$(METHDOBJS): psb_prec_mod.o psb_krylov_mod.o symlink: (/bin/ln -fs $(PSBKRYLDIR)/*.f90 . ) diff --git a/mlprec/mld_das_bld.f90 b/mlprec/mld_das_bld.f90 index 8e274a42..63fc8db2 100644 --- a/mlprec/mld_das_bld.f90 +++ b/mlprec/mld_das_bld.f90 @@ -41,7 +41,9 @@ ! ! This routine builds the Additive Schwarz (AS) preconditioner. ! If the preconditioner is the block-Jacobi one, the routine makes only a copy of -! the descriptor of the original matrix. +! the descriptor of the original matrix and then proceeds to call mld_fact_bld +! for LU or incomplete LU factorization of the diagonal blocks of the +! distributed matrix. ! ! ! Arguments: diff --git a/mlprec/mld_dbaseprec_bld.f90 b/mlprec/mld_dbaseprec_bld.f90 index 9e19076f..6ad6a6c5 100644 --- a/mlprec/mld_dbaseprec_bld.f90 +++ b/mlprec/mld_dbaseprec_bld.f90 @@ -154,32 +154,7 @@ subroutine mld_dbaseprc_bld(a,desc_a,p,info,upd) goto 9999 end if - case(mld_bjac_) - ! Block Jacobi preconditioner/smoother - - call mld_check_def(p%iprcparm(mld_sub_ren_),'renumbering',& - & mld_renum_none_,is_legal_renum) - call mld_check_def(p%iprcparm(mld_sub_solve_),'fact',& - & mld_ilu_n_,is_legal_ml_fact) - - call psb_cdcpy(desc_a,p%desc_data,info) - if(info /= 0) then - info=4010 - ch_err='psb_cdcpy' - call psb_errpush(info,name,a_err=ch_err) - goto 9999 - end if - - ! Build the local part of the base preconditioner/smoother - call mld_fact_bld(a,p,iupd,info) - if(info /= 0) then - info=4010 - call psb_errpush(info,name,a_err='mld_fact_bld') - goto 9999 - end if - - - case(mld_as_) + case(mld_bjac_,mld_as_) ! Additive Schwarz preconditioners/smoothers call mld_check_def(p%iprcparm(mld_n_ovr_),'overlap',& diff --git a/mlprec/mld_dfact_bld.f90 b/mlprec/mld_dfact_bld.f90 index 26a71151..6063a794 100644 --- a/mlprec/mld_dfact_bld.f90 +++ b/mlprec/mld_dfact_bld.f90 @@ -39,34 +39,35 @@ ! Subroutine: mld_dfact_bld ! Version: real ! -! This routine computes an LU or incomplete LU factorization of the input -! matrix, according to the value of p%iprcparm(iprcparm(sub_solve_), +! This routine computes an LU or incomplete LU factorization of the diagonal blocks +! of a distributed matrix, according to the value of p%iprcparm(iprcparm(sub_solve_), ! set by the user through mld_dprecinit or mld_dprecset. -! It may also split the local matrix into its block-diagonal and +! It may also split the matrix into its block-diagonal and ! off block-diagonal parts, for the future application of multiple ! block-Jacobi sweeps. ! -! This routine is used by mld_dbaseprec_bld, to build a 'base' block-Jacobi or +! This routine is used by mld_as_bld, to build a 'base' block-Jacobi or ! Additive Schwarz (AS) preconditioner at any level of a multilevel preconditioner, ! or a block-Jacobi or LU or ILU solver at the coarsest level of a multilevel ! preconditioner. For the Additive Schwarz, it is called from mld_as_bld, ! which prepares the overlap descriptor and retrieves the remote rows into blck. ! -! More precisely, the routine performs one of the following tasks: +! More precisely, the routine performs one of the following tasks: ! -! 1. construction of a block-Jacobi preconditioner associated -! to a matrix A distributed among the processes (allowed at any level); +! 1. LU or ILU factorization of the diagonal blocks of the distributed matrix +! for the construction of a block-Jacobi or AS preconditioners +! (allowed at any level); ! ! 2. setup of block-Jacobi sweeps to compute an approximate solution of a ! linear system ! A*Y = X, ! distributed among the processes (allowed only at the coarsest level); ! -! 3. LU factorization of a linear system +! 3. LU factorization of the matrix of a linear system ! A*Y = X, ! distributed among the processes (allowed only at the coarsest level); ! -! 4. LU or incomplete LU factorization of a linear system +! 4. LU or incomplete LU factorization of the matrix of a linear system ! A*Y = X, ! replicated on the processes (allowed only at the coarsest level). ! @@ -84,7 +85,7 @@ ! Arguments: ! a - type(psb_dspmat_type), input. ! The sparse matrix structure containing the local part of the -! matrix to be preconditioned or factorized. +! distributed matrix. ! p - type(mld_dbaseprec_type), input/output. ! The 'base preconditioner' data structure containing the local ! part of the preconditioner or solver at the current level. @@ -98,7 +99,7 @@ ! the new preconditioner. ! blck - type(psb_dspmat_type), input, optional. ! The sparse matrix structure containing the remote rows of the -! matrix to be factorized, that have been retrieved by mld_as_bld +! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 blck is empty. ! @@ -290,8 +291,6 @@ subroutine mld_dfact_bld(a,p,upd,info,blck) ! No reordering of the local matrix is required ! case(0) - - ! ! In case of multiple block-Jacobi sweeps, clip into p%av(ap_nd_) ! the off block-diagonal part of the local extended matrix. The @@ -333,7 +332,6 @@ subroutine mld_dfact_bld(a,p,upd,info,blck) goto 9999 end if end if - ! ! Compute a factorization of the diagonal block of the local matrix, ! according to the choice made by the user by setting p%iprcparm(sub_solve_) diff --git a/mlprec/mld_dilu0_fact.f90 b/mlprec/mld_dilu0_fact.f90 index faaf9fe3..38fd94f4 100644 --- a/mlprec/mld_dilu0_fact.f90 +++ b/mlprec/mld_dilu0_fact.f90 @@ -41,20 +41,20 @@ ! Contains: mld_dilu0_factint, ilu_copyin ! ! This routine computes either the ILU(0) or the MILU(0) factorization of the -! local part of the matrix stored into a. These factorizations are used to -! build the 'base preconditioner' (block-Jacobi preconditioner/solver, Additive -! Schwarz preconditioner) corresponding to a certain level of a multilevel +! diagonal blocks of a distributed matrix. These factorizations +! are used to build the 'base preconditioner' (block-Jacobi preconditioner/solver, +! Additive Schwarz preconditioner) corresponding to a given level of a multilevel ! preconditioner. ! ! Details on the above factorizations can be found in ! Y. Saad, Iterative Methods for Sparse Linear Systems, Second Edition, ! SIAM, 2003, Chapter 10. ! -! The local matrix to be factorized is stored into a and blck, as specified -! in the description of the arguments below. The storage format for both the -! L and U factors is CSR. The diagonal of the U factor is stored separately -! (actually, the inverse of the diagonal entries is stored; this is then -! managed in the solve stage associated to the ILU(0)/MILU(0) factorization). +! The local matrix is stored into a and blck, as specified in the description +! of the arguments below. The storage format for both the L and U factors is CSR. +! The diagonal of the U factor is stored separately (actually, the inverse of the +! diagonal entries is stored; this is then managed in the solve stage associated +! to the ILU(0)/MILU(0) factorization). ! ! The routine copies and factors "on the fly" from a and blck into l (L factor), ! u (U factor, except its diagonal) and d (diagonal of U). @@ -69,11 +69,11 @@ ! The MILU(0) factorization is computed if ialg = 2 (= mld_milu_n_); ! the ILU(0) factorization otherwise. ! a - type(psb_dspmat_type), input. -! The sparse matrix structure containing the local matrix to be -! factorized. Note that if the 'base' Additive Schwarz preconditioner +! The sparse matrix structure containing the local matrix. +! Note that if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered -! (see mld_bjac_bld), then a contains only the 'original' local part -! of the matrix to be factorized, i.e. the rows of the matrix held +! (see mld_as_bld), then a contains only the 'original' local part +! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! l - type(psb_dspmat_type), input/output. ! The L factor in the incomplete factorization. @@ -92,10 +92,10 @@ ! Error code. ! blck - type(psb_dspmat_type), input, optional, target. ! The sparse matrix structure containing the remote rows of the -! matrix to be factorized, that have been retrieved by mld_asmat_bld +! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered -! (see mld_bjac_bld), then blck is empty. +! (see mld_fact_bld), then blck is empty. ! subroutine mld_dilu0_fact(ialg,a,l,u,d,info,blck) @@ -205,12 +205,13 @@ contains ! Version: real ! Note: internal subroutine of mld_dilu0_fact. ! - ! This routine computes either the ILU(0) or the MILU(0) factorization of the - ! local part of the matrix stored into a. These factorizations are used to build - ! the 'base preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz - ! preconditioner) corresponding to a certain level of a multilevel preconditioner. + ! This routine computes either the ILU(0) or the MILU(0) factorization of the + ! diagonal blocks of a distributed matrix. + ! These factorizations are used to build the 'base preconditioner' + ! (block-Jacobi preconditioner/solver, Additive Schwarz + ! preconditioner) corresponding to a given level of a multilevel preconditioner. ! - ! The local matrix to be factorized is stored into a and b, as specified in the + ! The local matrix is stored into a and b, as specified in the ! description of the arguments below. The storage format for both the L and U ! factors is CSR. The diagonal of the U factor is stored separately (actually, ! the inverse of the diagonal entries is stored; this is then managed in the @@ -232,20 +233,20 @@ contains ! ma - integer, input ! The number of rows of the local submatrix stored into a. ! a - type(psb_dspmat_type), input. - ! The sparse matrix structure containing the local matrix to be - ! factorized. Note that, if the 'base' Additive Schwarz preconditioner + ! The sparse matrix structure containing the local matrix. + ! Note that, if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered - ! (see mld_bjac_bld), then a contains only the 'original' local part - ! of the matrix to be factorized, i.e. the rows of the matrix held + ! (see mld_fact_bld), then a contains only the 'original' local part + ! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! mb - integer, input. ! The number of rows of the local submatrix stored into b. ! b - type(psb_dspmat_type), input. ! The sparse matrix structure containing the remote rows of the - ! matrix to be factorized, that have been retrieved by mld_asmat_bld + ! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been - ! reordered (see mld_bjac_bld), then b does not contain any row. + ! reordered (see mld_fact_bld), then b does not contain any row. ! d - real(kind(1.d0)), dimension(:), output. ! The inverse of the diagonal entries of the U factor in the ! incomplete factorization. @@ -486,7 +487,7 @@ contains ! copied into laspk, dia, uaspk row by row, through successive calls to ! ilu_copyin. ! - ! The routine is used by mld_dilu0_factin in the computation of the ILU(0)/MILU(0) + ! The routine is used by mld_dilu0_factint in the computation of the ILU(0)/MILU(0) ! factorization of a local sparse matrix. ! ! TODO: modify the routine to allow copying into output L and U that are diff --git a/mlprec/mld_dilu_bld.f90 b/mlprec/mld_dilu_bld.f90 index d3ae2696..f1c0aebc 100644 --- a/mlprec/mld_dilu_bld.f90 +++ b/mlprec/mld_dilu_bld.f90 @@ -39,9 +39,9 @@ ! Subroutine: mld_dilu_bld ! Version: real ! -! This routine computes an incomplete LU (ILU) factorization of the local part -! of the matrix stored into a. This factorization is used to build the 'base -! preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz +! This routine computes an incomplete LU (ILU) factorization of the diagonal blocks +! of a distributed matrix. This factorization is used to build +! the 'base preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz ! preconditioner) corresponding to a certain level of a multilevel preconditioner. ! The following factorizations are available: ! - ILU(k), i.e. ILU factorization with fill-in level k, @@ -62,12 +62,12 @@ ! ! Arguments: ! a - type(psb_dspmat_type), input. -! The sparse matrix structure containing the local matrix to be -! factorized. Note that if p%iprcparm(mld_n_ovr_) > 0, i.e. the +! The sparse matrix structure containing the local matrix. +! Note that if p%iprcparm(mld_n_ovr_) > 0, i.e. the ! 'base' Additive Schwarz preconditioner has overlap greater than -! 0, and p%iprcparm(mld_sub_ren_) = 0, i.e. a reordering of the -! matrix has not been performed (see mld_bjac_bld), then a contains -! only the 'original' local part of the matrix to be factorized, +! 0, and p%iprcparm(mld_sub_ren_) = 0, i.e. a reordering of the +! matrix has not been performed (see mld_fact_bld), then a contains +! only the 'original' local part of the distributed matrix, ! i.e. the rows of the matrix held by the calling process according ! to the initial data distribution. ! p - type(mld_dbaseprc_type), input/output. @@ -81,10 +81,10 @@ ! Error code. ! blck - type(psb_dspmat_type), input, optional. ! The sparse matrix structure containing the remote rows of the -! matrix to be factorized, that have been retrieved by mld_asmat_bld +! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered -! (see mld_bjac_bld), then blck does not contain any row. +! (see mld_fact_bld), then blck does not contain any row. ! subroutine mld_dilu_bld(a,p,upd,info,blck) diff --git a/mlprec/mld_diluk_fact.f90 b/mlprec/mld_diluk_fact.f90 index adf24ce2..c915e64c 100644 --- a/mlprec/mld_diluk_fact.f90 +++ b/mlprec/mld_diluk_fact.f90 @@ -41,16 +41,16 @@ ! Contains: mld_diluk_factint, iluk_copyin, iluk_fact, iluk_copyout. ! ! This routine computes either the ILU(k) or the MILU(k) factorization of the -! local part of the matrix stored into a. These factorizations are used to -! build the 'base preconditioner' (block-Jacobi preconditioner/solver, Additive -! Schwarz preconditioner) corresponding to a certain level of a multilevel +! diagonal blocks of a distributed matrix. These factorizations are used to build +! the 'base preconditioner' (block-Jacobi preconditioner/solver, +! Additive Schwarz preconditioner) corresponding to a certain level of a multilevel ! preconditioner. ! ! Details on the above factorizations can be found in ! Y. Saad, Iterative Methods for Sparse Linear Systems, Second Edition, ! SIAM, 2003, Chapter 10. ! -! The local matrix to be factorized is stored into a and blck, as specified in +! The local matrix is stored into a and blck, as specified in ! the description of the arguments below. The storage format for both the L and ! U factors is CSR. The diagonal of the U factor is stored separately (actually, ! the inverse of the diagonal entries is stored; this is then managed in the solve @@ -66,11 +66,11 @@ ! the MILU(k) one if ialg = 2 (= mld_milu_n_); other values are ! not allowed. ! a - type(psb_dspmat_type), input. -! The sparse matrix structure containing the local matrix to be -! factorized. Note that if the 'base' Additive Schwarz preconditioner +! The sparse matrix structure containing the local matrix. +! Note that if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered -! (see mld_bjac_bld), then a contains only the 'original' local part -! of the matrix to be factorized, i.e. the rows of the matrix held +! (see mld_fact_bld), then a contains only the 'original' local part +! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! l - type(psb_dspmat_type), input/output. ! The L factor in the incomplete factorization. @@ -89,10 +89,10 @@ ! Error code. ! blck - type(psb_dspmat_type), input, optional, target. ! The sparse matrix structure containing the remote rows of the -! matrix to be factorized, that have been retrieved by mld_asmat_bld +! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered -! (see mld_bjac_bld), then blck does not contain any row. +! (see mld_fact_bld), then blck does not contain any row. ! subroutine mld_diluk_fact(fill_in,ialg,a,l,u,d,info,blck) @@ -201,11 +201,11 @@ contains ! Note: internal subroutine of mld_diluk_fact ! ! This routine computes either the ILU(k) or the MILU(k) factorization of the - ! local part of the matrix stored into a. These factorizations are used to build + ! diagonal blocks of a distributed matrix. These factorizations are used to build ! the 'base preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz ! preconditioner) corresponding to a certain level of a multilevel preconditioner. ! - ! The local matrix to be factorized is stored into a and b, as specified in the + ! The local matrix is stored into a and b, as specified in the ! description of the arguments below. The storage format for both the L and U ! factors is CSR. The diagonal of the U factor is stored separately (actually, ! the inverse of the diagonal entries is stored; this is then managed in the @@ -223,18 +223,18 @@ contains ! The total number of rows of the local matrix to be factorized, ! i.e. ma+mb. ! a - type(psb_dspmat_type), input. - ! The sparse matrix structure containing the local matrix to be - ! factorized. Note that, if the 'base' Additive Schwarz preconditioner + ! The sparse matrix structure containing the local matrix. + ! Note that, if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered - ! (see mld_bjac_bld), then a contains only the 'original' local part - ! of the matrix to be factorized, i.e. the rows of the matrix held + ! (see mld_fact_bld), then a contains only the 'original' local part + ! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! b - type(psb_dspmat_type), input. ! The sparse matrix structure containing the remote rows of the - ! matrix to be factorized, that have been retrieved by mld_asmat_bld + ! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered - ! (see mld_bjac_bld), then b does not contain any row. + ! (see mld_fact_bld), then b does not contain any row. ! d - real(kind(1.d0)), dimension(:), output. ! The inverse of the diagonal entries of the U factor in the incomplete ! factorization. diff --git a/mlprec/mld_dilut_fact.f90 b/mlprec/mld_dilut_fact.f90 index 7b2f48c6..0b9e8a43 100644 --- a/mlprec/mld_dilut_fact.f90 +++ b/mlprec/mld_dilut_fact.f90 @@ -40,8 +40,8 @@ ! Version: real ! Contains: mld_dilut_factint, ilut_copyin, ilut_fact, ilut_copyout ! -! This routine computes the ILU(k,t) factorization of the local part of the -! matrix stored into a. This factorization is used to build the 'base +! This routine computes the ILU(k,t) factorization of the diagonal blocks of a +! distributed matrix. This factorization is used to build the 'base ! preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz ! preconditioner) corresponding to a certain level of a multilevel preconditioner. ! @@ -49,7 +49,7 @@ ! Y. Saad, Iterative Methods for Sparse Linear Systems, Second Edition, ! SIAM, 2003, Chapter 10. ! -! The local matrix to be factorized is stored into a and blck, as specified in +! The local matrix is stored into a and blck, as specified in ! the description of the arguments below. The storage format for both the L and ! U factors is CSR. The diagonal of the U factor is stored separately (actually, ! the inverse of the diagonal entries is stored; this is then managed in the solve @@ -62,11 +62,11 @@ ! thres - integer, input. ! The threshold t, i.e. the drop tolerance, in ILU(k,t). ! a - type(psb_dspmat_type), input. -! The sparse matrix structure containing the local matrix to be -! factorized. Note that if the 'base' Additive Schwarz preconditioner +! The sparse matrix structure containing the local matrix. +! Note that if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered -! (see mld_bjac_bld), then a contains only the 'original' local part -! of the matrix to be factorized, i.e. the rows of the matrix held +! (see mld_fact_bld), then a contains only the 'original' local part +! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! l - type(psb_dspmat_type), input/output. ! The L factor in the incomplete factorization. @@ -85,10 +85,10 @@ ! Error code. ! blck - type(psb_dspmat_type), input, optional, target. ! The sparse matrix structure containing the remote rows of the -! matrix to be factorized, that have been retrieved by mld_asmat_bld +! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered -! (see mld_bjac_bld), then blck does not contain any row. +! (see mld_fact_bld), then blck does not contain any row. ! subroutine mld_dilut_fact(fill_in,thres,a,l,u,d,info,blck) @@ -203,8 +203,8 @@ contains ! Version: real ! Note: internal subroutine of mld_dilut_fact ! - ! This routine computes the ILU(k,t) factorization of the local part of the - ! matrix stored into a. These factorization is used to build the 'base + ! This routine computes the ILU(k,t) factorization of the diagonal blocks of a + ! distributed matrix. This factorization is used to build the 'base ! preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz ! preconditioner) corresponding to a certain level of a multilevel preconditioner. ! @@ -224,18 +224,18 @@ contains ! The total number of rows of the local matrix to be factorized, ! i.e. ma+mb. ! a - type(psb_dspmat_type), input. - ! The sparse matrix structure containing the local matrix to be - ! factorized. Note that, if the 'base' Additive Schwarz preconditioner + ! The sparse matrix structure containing the local matrix. + ! Note that, if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered - ! (see mld_bjac_bld), then a contains only the 'original' local part - ! of the matrix to be factorized, i.e. the rows of the matrix held + ! (see mld_fact_bld), then a contains only the 'original' local part + ! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! b - type(psb_dspmat_type), input. ! The sparse matrix structure containing the remote rows of the - ! matrix to be factorized, that have been retrieved by mld_asmat_bld + ! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered - ! (see mld_bjac_bld), then b does not contain any row. + ! (see mld_fact_bld), then b does not contain any row. ! d - real(kind(1.d0)), dimension(:), output. ! The inverse of the diagonal entries of the U factor in the incomplete ! factorization. @@ -666,8 +666,7 @@ contains ! Note: this argument is intent(inout) and not only intent(out) ! to retain its allocation, done by this routine. ! - subroutine ilut_fact(thres,i,nrmi,row,heap,& - & d,uia1,uia2,uaspk,nidx,idxs,info) + subroutine ilut_fact(thres,i,nrmi,row,heap,d,uia1,uia2,uaspk,nidx,idxs,info) use psb_base_mod @@ -787,7 +786,7 @@ contains ! - the array row is re-initialized for future use in mld_ilut_fact(see also ! ilut_copyin and ilut_fact). ! - ! This routine is used by mld_dilut_factint in the computation of the ILU(k,t) + ! This routine is used by mld_dilut_factint in the computation of the ILU(k,t) ! factorization of a local sparse matrix. ! ! diff --git a/mlprec/mld_dmlprec_aply.f90 b/mlprec/mld_dmlprec_aply.f90 index 0614e710..466b7d5c 100644 --- a/mlprec/mld_dmlprec_aply.f90 +++ b/mlprec/mld_dmlprec_aply.f90 @@ -49,15 +49,15 @@ ! - X and Y are vectors, ! - alpha and beta are scalars. ! -! For each level we have as many subdomains as processes (except for the coarsest +! For each level we have as many submatrices as processes (except for the coarsest ! level where we might have a replicated index space) and each process takes care -! of one subdomain. +! of one submatrix. ! ! The multilevel preconditioner M is regarded as an array of 'base preconditioners', ! each representing the part of the preconditioner associated to a certain level. ! For each level ilev, the base preconditioner K(ilev) is stored in baseprecv(ilev) ! and is associated to a matrix A(ilev), obtained by 'tranferring' the original -! matrix A (i.e. the matrix to be preconditioned) to level ilev, through smoothed +! matrix A (i.e. the matrix to be preconditioned) to the level ilev, through smoothed ! aggregation. ! ! The levels are numbered in increasing order starting from the finest one, i.e. @@ -157,8 +157,8 @@ ! Error code. ! ! Note that when the LU factorization of the matrix A(ilev) is computed instead of -! the ILU one, by using UMFPACK or SuperLU_dist, the corresponding L and U factors -! are stored in data structures provided by UMFPACK or SuperLU_dist and pointed by +! the ILU one, by using UMFPACK or SuperLU, the corresponding L and U factors +! are stored in data structures provided by UMFPACK or SuperLU and pointed by ! baseprecv(ilev)%iprcparm(mld_umf_ptr) or baseprecv(ilev)%iprcparm(mld_slu_ptr), ! respectively. ! @@ -285,8 +285,45 @@ subroutine mld_dmlprec_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info) return contains - + ! + ! Subroutine: add_ml_aply + ! Version: real + ! Note: internal subroutine of mld_dmlprec_aply. + ! This routine computes + ! + ! Y = beta*Y + alpha*op(M^(-1))*X, + ! where + ! - M is a multilevel domain decomposition (Schwarz) preconditioner associated + ! to a certain matrix A and stored in the array baseprecv, + ! - op(M^(-1)) is M^(-1) or its transpose, according to the value of trans, + ! - X and Y are vectors, + ! - alpha and beta are scalars. + ! + ! For each level we have as many submatrices as processes (except for the coarsest + ! level where we might have a replicated index space) and each process takes care + ! of one submatrix. + ! + ! The multilevel preconditioner M is regarded as an array of 'base preconditioners', + ! each representing the part of the preconditioner associated to a certain level. + ! For each level ilev, the base preconditioner K(ilev) is stored in baseprecv(ilev) + ! and is associated to a matrix A(ilev), obtained by 'tranferring' the original + ! matrix A (i.e. the matrix to be preconditioned) to the level ilev, through smoothed + ! aggregation. + ! + ! The levels are numbered in increasing order starting from the finest one, i.e. + ! level 1 is the finest level and A(1) is the matrix A. + ! This routine applies the multilevel preconditioner in an additive + ! way (additive through the levels and additive on the same level). + ! For details on the additive multilevel Schwarz preconditioner see + ! the Algorithm 3.1.1 in the book: + ! - B.F. Smith, P.E. Bjorstad & W.D. Gropp, + ! Domain decomposition: parallel multilevel methods for elliptic partial + ! differential equations, Cambridge University Press, 1996. + ! + ! For a detailed description of the arguments, see mld_dmlprec_aply. + ! subroutine add_ml_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info) + ! implicit none ! Arguments type(psb_desc_type),intent(in) :: desc_data @@ -329,39 +366,6 @@ contains call psb_errpush(4010,name,a_err='Allocate') goto 9999 end if - - ! - ! Additive multilevel - ! - ! 1. ! Apply the base preconditioner at level 1. - ! ! The sum over the subdomains is carried out in the - ! ! application of K(1). - ! X(1) = Xest - ! Y(1) = (K(1)^(-1))*X(1) - ! - ! 2. DO ilev=2,nlev - ! - ! ! Transfer X(ilev-1) to the next coarser level. - ! X(ilev) = AV(ilev; sm_pr_t_)*X(ilev-1) - ! - ! ! Apply the base preconditioner at the current level. - ! ! The sum over the subdomains is carried out in the - ! ! application of K(ilev). - ! Y(ilev) = (K(ilev)^(-1))*X(ilev) - ! - ! ENDDO - ! - ! 3. DO ilev=nlev-1,1,-1 - ! - ! ! Transfer Y(ilev+1) to the next finer level. - ! Y(ilev) = AV(ilev+1; sm_pr_)*Y(ilev+1) - ! - ! ENDDO - ! - ! 4. Yext = beta*Yext + alpha*Y(1) - ! - - ! ! STEP 1 ! @@ -385,7 +389,6 @@ contains call psb_errpush(4010,name,a_err='baseprec_aply') goto 9999 end if - ! ! STEP 2 ! @@ -524,9 +527,46 @@ contains return end subroutine add_ml_aply - - + ! + ! Subroutine: mlt_pre_ml_aply + ! Version: real + ! Note: internal subroutine of mld_dmlprec_aply. + ! This routine computes + ! + ! Y = beta*Y + alpha*op(M^(-1))*X, + ! where + ! - M is a multilevel domain decomposition (Schwarz) preconditioner associated + ! to a certain matrix A and stored in the array baseprecv, + ! - op(M^(-1)) is M^(-1) or its transpose, according to the value of trans, + ! - X and Y are vectors, + ! - alpha and beta are scalars. + ! + ! For each level we have as many submatrices as processes (except for the coarsest + ! level where we might have a replicated index space) and each process takes care + ! of one submatrix. + ! + ! The multilevel preconditioner M is regarded as an array of 'base preconditioners', + ! each representing the part of the preconditioner associated to a certain level. + ! For each level ilev, the base preconditioner K(ilev) is stored in baseprecv(ilev) + ! and is associated to a matrix A(ilev), obtained by 'tranferring' the original + ! matrix A (i.e. the matrix to be preconditioned) to the level ilev, through smoothed + ! aggregation. + ! + ! The levels are numbered in increasing order starting from the finest one, i.e. + ! level 1 is the finest level and A(1) is the matrix A. + ! + ! This routine applies the multilevel preconditioner in a hybrid way + ! (multiplicative through the levels and additive on the same level) + ! and pre-smoothing. + ! For details on pre-smoothed hybrid multiplicative multilevel Schwarz preconditioner, + ! see the Algorithm 3.2.1 in the book: + ! - B.F. Smith, P.E. Bjorstad & W.D. Gropp, + ! Domain decomposition: parallel multilevel methods for elliptic partial + ! differential equations, Cambridge University Press, 1996. + ! + ! For a detailed description of the arguments, see mld_dmlprec_aply. subroutine mlt_pre_ml_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info) + ! implicit none ! Arguments type(psb_desc_type),intent(in) :: desc_data @@ -570,43 +610,6 @@ contains goto 9999 end if - ! - ! Pre-smoothing - ! - ! 1. X(1) = Xext - ! - ! 2. ! Apply the base preconditioner at the finest level. - ! Y(1) = (K(1)^(-1))*X(1) - ! - ! 3. ! Compute the residual at the finest level. - ! TX(1) = X(1) - A(1)*Y(1) - ! - ! 4. DO ilev=2, nlev - ! - ! ! Transfer the residual to the current (coarser) level. - ! X(ilev) = AV(ilev; sm_pr_t_)*TX(ilev-1) - ! - ! ! Apply the base preconditioner at the current level. - ! ! The sum over the subdomains is carried out in the - ! ! application of K(ilev). - ! Y(ilev) = (K(ilev)^(-1))*X(ilev) - ! - ! ! Compute the residual at the current level (except at - ! ! the coarsest level). - ! IF (ilev < nlev) - ! TX(ilev) = (X(ilev)-A(ilev)*Y(ilev)) - ! - ! ENDDO - ! - ! 5. DO ilev=nlev-1,1,-1 - ! - ! ! Transfer Y(ilev+1) to the next finer level - ! Y(ilev) = Y(ilev) + AV(ilev+1; sm_pr_)*Y(ilev+1) - ! - ! ENDDO - ! - ! 6. Yext = beta*Yext + alpha*Y(1) - ! ! ! STEP 1 @@ -800,9 +803,45 @@ contains end if return end subroutine mlt_pre_ml_aply - - + ! + ! Subroutine: mlt_post_ml_aply + ! Version: real + ! Note: internal subroutine of mld_dmlprec_aply. + ! This routine computes + ! + ! Y = beta*Y + alpha*op(M^(-1))*X, + ! where + ! - M is a multilevel domain decomposition (Schwarz) preconditioner associated + ! to a certain matrix A and stored in the array baseprecv, + ! - op(M^(-1)) is M^(-1) or its transpose, according to the value of trans, + ! - X and Y are vectors, + ! - alpha and beta are scalars. + ! + ! For each level we have as many submatrices as processes (except for the coarsest + ! level where we might have a replicated index space) and each process takes care + ! of one submatrix. + ! + ! The multilevel preconditioner M is regarded as an array of 'base preconditioners', + ! each representing the part of the preconditioner associated to a certain level. + ! For each level ilev, the base preconditioner K(ilev) is stored in baseprecv(ilev) + ! and is associated to a matrix A(ilev), obtained by 'tranferring' the original + ! matrix A (i.e. the matrix to be preconditioned) to the level ilev, through smoothed + ! aggregation. + ! + ! The levels are numbered in increasing order starting from the finest one, i.e. + ! level 1 is the finest level and A(1) is the matrix A. + ! + ! This routine applies the multilevel preconditioner in a hybrid way + ! (multiplicative through the levels and additive on the same level) + ! and post-smoothing. + ! For details on hybrid multiplicative multilevel Schwarz preconditioners, see + ! - B.F. Smith, P.E. Bjorstad & W.D. Gropp, + ! Domain decomposition: parallel multilevel methods for elliptic partial + ! differential equations, Cambridge University Press, 1996. + ! + ! For a detailed description of the arguments, see mld_dmlprec_aply. subroutine mlt_post_ml_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info) + ! implicit none ! Arguments type(psb_desc_type),intent(in) :: desc_data @@ -845,37 +884,6 @@ contains call psb_errpush(4010,name,a_err='Allocate') goto 9999 end if - - ! - ! Post-smoothing - ! - ! 1. X(1) = Xext - ! - ! 2. DO ilev=2, nlev - ! - ! ! Transfer X(ilev-1) to the next coarser level. - ! X(ilev) = AV(ilev; sm_pr_t_)*X(ilev-1) - ! - ! ENDDO - ! - ! 3.! Apply the preconditioner at the coarsest level. - ! Y(nlev) = (K(nlev)^(-1))*X(nlev) - ! - ! 4. DO ilev=nlev-1,1,-1 - ! - ! ! Transfer Y(ilev+1) to the next finer level. - ! Y(ilev) = AV(ilev+1; sm_pr_)*Y(ilev+1) - ! - ! ! Compute the residual at the current level and apply to it the - ! ! base preconditioner. The sum over the subdomains is carried out - ! ! in the application of K(ilev). - ! Y(ilev) = Y(ilev) + (K(ilev)^(-1))*(X(ilev)-A(ilev)*Y(ilev)) - ! - ! ENDDO - ! - ! 5. Yext = beta*Yext + alpha*Y(1) - ! - ! ! STEP 1 ! @@ -1091,9 +1099,47 @@ contains end if return end subroutine mlt_post_ml_aply - - + ! + ! Subroutine: mlt_twoside_ml_aply + ! Version: real + ! Note: internal subroutine of mld_dmlprec_aply. + ! This routine computes + ! + ! Y = beta*Y + alpha*op(M^(-1))*X, + ! where + ! - M is a multilevel domain decomposition (Schwarz) preconditioner associated + ! to a certain matrix A and stored in the array baseprecv, + ! - op(M^(-1)) is M^(-1) or its transpose, according to the value of trans, + ! - X and Y are vectors, + ! - alpha and beta are scalars. + ! + ! For each level we have as many submatrices as processes (except for the coarsest + ! level where we might have a replicated index space) and each process takes care + ! of one submatrix. + ! + ! The multilevel preconditioner M is regarded as an array of 'base preconditioners', + ! each representing the part of the preconditioner associated to a certain level. + ! For each level ilev, the base preconditioner K(ilev) is stored in baseprecv(ilev) + ! and is associated to a matrix A(ilev), obtained by 'tranferring' the original + ! matrix A (i.e. the matrix to be preconditioned) to the level ilev, through smoothed + ! aggregation. + ! + ! The levels are numbered in increasing order starting from the finest one, i.e. + ! level 1 is the finest level and A(1) is the matrix A. + ! + ! This routine applies the multilevel preconditioner in a symmetrized hybrid way + ! (multiplicative through the levels and additive on the same level, with pre and + ! post-smoothing). + ! For details on the symmetrized hybrid multiplicative multilevel Schwarz + ! preconditioners, see the Algorithm 3.2.2 of the book: + ! - B.F. Smith, P.E. Bjorstad & W.D. Gropp, + ! Domain decomposition: parallel multilevel methods for elliptic partial + ! differential equations, Cambridge University Press, 1996. + ! + ! For a detailed description of the arguments, see mld_dmlprec_aply. + ! subroutine mlt_twoside_ml_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info) + ! implicit none ! Arguments type(psb_desc_type),intent(in) :: desc_data @@ -1136,49 +1182,6 @@ contains call psb_errpush(4010,name,a_err='Allocate') goto 9999 end if - - ! - ! Pre- and post-smoothing (symmetrized) - ! - ! 1. X(1) = Xext - ! - ! 2. ! Apply the base peconditioner at the finest level - ! Y(1) = (K(1)^(-1))*X(1) - ! - ! 3. ! Compute the residual at the finest level - ! TX(1) = X(1) - A(1)*Y(1) - ! - ! 4. DO ilev=2, nlev - ! - ! ! Transfer the residual to the current (coarser) level - ! X(ilev) = AV(ilev; sm_pr_t)*TX(ilev-1) - ! - ! ! Apply the base preconditioner at the current level. - ! ! The sum over the subdomains is carried out in the - ! ! application of K(ilev) - ! Y(ilev) = (K(ilev)^(-1))*X(ilev) - ! - ! ! Compute the residual at the current level - ! TX(ilev) = (X(ilev)-A(ilev)*Y(ilev)) - ! - ! ENDDO - ! - ! 5. DO ilev=NLEV-1,1,-1 - ! - ! ! Transfer Y(ilev+1) to the next finer level - ! Y(ilev) = Y(ilev) + AV(ilev+1; sm_pr_)*Y(ilev+1) - ! - ! ! Compute the residual at the current level and apply to it the - ! ! base preconditioner. The sum over the subdomains is carried out - ! ! in the application of K(ilev) - ! Y(ilev) = Y(ilev) + (K(ilev)**(-1))*(X(ilev)-A(ilev)*Y(ilev)) - ! - ! ENDDO - ! - ! 6. Yext = beta*Yext + alpha*Y(1) - ! - - ! ! STEP 1 ! ! Copy the input vector X diff --git a/mlprec/mld_dprec_aply.f90 b/mlprec/mld_dprec_aply.f90 index d882bd72..ab27fb9b 100644 --- a/mlprec/mld_dprec_aply.f90 +++ b/mlprec/mld_dprec_aply.f90 @@ -141,7 +141,7 @@ subroutine mld_dprec_aply(prec,x,y,desc_data,info,trans,work) endif ! If the original distribution has an overlap we should fix that. - call psb_ovrl(y,desc_data,info,update=psb_avg_) + call psb_halo(y,desc_data,info,data=psb_comm_mov_) if (present(work)) then diff --git a/mlprec/mld_dprecinit.f90 b/mlprec/mld_dprecinit.f90 index 1d46808d..a1efe40a 100644 --- a/mlprec/mld_dprecinit.f90 +++ b/mlprec/mld_dprecinit.f90 @@ -210,7 +210,7 @@ subroutine mld_dprecinit(p,ptype,info,nlev) p%baseprecv(ilev_)%iprcparm(mld_n_ovr_) = 0 p%baseprecv(ilev_)%iprcparm(mld_ml_type_) = mld_mult_ml_ p%baseprecv(ilev_)%iprcparm(mld_aggr_alg_) = mld_dec_aggr_ - p%baseprecv(ilev_)%iprcparm(mld_aggr_kind_) = mld_smooth_prol_ + p%baseprecv(ilev_)%iprcparm(mld_aggr_kind_) = mld_smooth_prol_ p%baseprecv(ilev_)%iprcparm(mld_coarse_mat_) = mld_distr_mat_ p%baseprecv(ilev_)%iprcparm(mld_smooth_pos_) = mld_post_smooth_ p%baseprecv(ilev_)%iprcparm(mld_aggr_eig_) = mld_max_norm_ diff --git a/mlprec/mld_dprecset.f90 b/mlprec/mld_dprecset.f90 index b45615d8..3bc4cce5 100644 --- a/mlprec/mld_dprecset.f90 +++ b/mlprec/mld_dprecset.f90 @@ -58,7 +58,7 @@ ! numbers, as reported in MLD2P4 user's guide. ! val - integer, input. ! The value of the parameter to be set. The list of allowed -! values is reported in MLD2P4 user's guide. +! values is reported in MLD2P4 user's guide. ! info - integer, output. ! Error code. ! ilev - integer, optional, input. @@ -276,7 +276,7 @@ subroutine mld_dprecsetc(p,what,string,info,ilev) implicit none -! Arguments + ! Arguments type(mld_dprec_type), intent(inout) :: p integer, intent(in) :: what character(len=*), intent(in) :: string @@ -527,7 +527,7 @@ subroutine mld_dprecsetd(p,what,val,info,ilev) implicit none -! Arguments + ! Arguments type(mld_dprec_type), intent(inout) :: p integer, intent(in) :: what real(kind(1.d0)), intent(in) :: val diff --git a/mlprec/mld_dsp_renum.f90 b/mlprec/mld_dsp_renum.f90 index 3282e0d1..0307cd82 100644 --- a/mlprec/mld_dsp_renum.f90 +++ b/mlprec/mld_dsp_renum.f90 @@ -50,7 +50,7 @@ ! description of the arguments below. ! ! If required by the user (p%iprcparm(mld_sub_ren_) /= 0), the routine is -! used by mld_bjac_bld in building the block-Jacobi and Additive Schwarz +! used by mld_fact_bld in building the block-Jacobi and Additive Schwarz ! 'base preconditioners' corresponding to any level of a multilevel ! preconditioner. ! @@ -63,7 +63,7 @@ ! distribution. ! blck - type(psb_dspmat_type), input. ! The sparse matrix structure containing the remote rows of the -! matrix to be reordered, that have been retrieved by mld_asmat_bld +! matrix to be reordered, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0.If the overlap is 0, then blck does not contain ! any row. diff --git a/mlprec/mld_dsub_aply.f90 b/mlprec/mld_dsub_aply.f90 index f384289b..b96ecb04 100644 --- a/mlprec/mld_dsub_aply.f90 +++ b/mlprec/mld_dsub_aply.f90 @@ -81,29 +81,28 @@ ! K^(-1), alpha = 1 and beta = 0. ! ! The block-Jacobi preconditioner or solver and the L and U factors of the LU -! or ILU factorizations have been built by the routine mld_dbjac_bld and stored -! into the 'base preconditioner' data structure prec. See mld_dbjac_bld for more +! or ILU factorizations have been built by the routine mld_fact_bld and stored +! into the 'base preconditioner' data structure prec. See mld_fact_bld for more ! details. ! -! This routine is used by mld_dbaseprec_aply, to apply a 'base' block-Jacobi or +! This routine is used by mld_as_aply, to apply a 'base' block-Jacobi or ! Additive Schwarz (AS) preconditioner at any level of a multilevel preconditioner, ! or a block-Jacobi or LU or ILU solver at the coarsest level of a multilevel ! preconditioner. ! -! Inside mld_dbaseprec_aply, tasks 1, 3 and 4 may be selected if -! prec%iprcparm(smooth_sweeps_) = 1, while task 2 if prec%iprcparm(smooth_sweeps_) -! > 1. Furthermore, tasks 1, 2 and 3 may be performed if the matrix A is +! Tasks 1, 3 and 4 are selected when prec%iprcparm(smooth_sweeps_) = 1, +! while task 2 is selected when prec%iprcparm(smooth_sweeps_) > 1. Furthermore +! Tasks 1, 2 and 3 may be performed when the matrix A is ! distributed among the processes (prec%iprcparm(mld_coarse_mat_) = mld_distr_mat_), -! while task 4 may be performed if A is replicated on the processes +! while task 4 may be performed when A is replicated on the processes ! (prec%iprcparm(mld_coarse_mat_) = mld_repl_mat_). Note that the matrix A is ! distributed among the processes at each level of the multilevel preconditioner, ! except the coarsest one, where it may be either distributed or replicated on -! the processes. Furthermore, the tasks 2, 3 and 4 are performed only at the -! coarsest level. Note also that this routine manages implicitly the fact that +! the processes. Tasks 2, 3 and 4 are performed only at the coarsest level. +! Note also that this routine manages implicitly the fact that ! the matrix is distributed or replicated, i.e. it does not make any explicit ! reference to the value of prec%iprcparm(mld_coarse_mat_). ! -! ! Arguments: ! ! alpha - real(kind(0.d0)), input. @@ -207,9 +206,7 @@ subroutine mld_dsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info) endif else if (prec%iprcparm(mld_smooth_sweeps_) > 1) then - ! - ! TASK 2 ! ! Apply prec%iprcparm(smooth_sweeps_) sweeps of a block-Jacobi solver ! to compute an approximate solution of a linear system. diff --git a/mlprec/mld_dsub_solve.f90 b/mlprec/mld_dsub_solve.f90 index 6c75c1ac..2090c481 100644 --- a/mlprec/mld_dsub_solve.f90 +++ b/mlprec/mld_dsub_solve.f90 @@ -44,11 +44,14 @@ ! Y = beta*Y + alpha*op(K^(-1))*X, ! ! where -! - K is a factored matrix, as specified below, +! - K is a factored matrix, as specified below, ! - op(K^(-1)) is K^(-1) or its transpose, according to the value of trans, ! - X and Y are vectors, ! - alpha and beta are scalars. ! +! Depending on K, alpha, beta (and on the communication descriptor desc_data +! - see the arguments below), the above computation may correspond to one of +! the following tasks: ! ! 1. Solution of a linear system with sparse factors LU generated by means ! of an incomplete factorization approximating @@ -58,7 +61,7 @@ ! they are also block-diagonal) or replicated. ! ! 2. Solution of a linear system with sparse factors LU generated by means -! of an complete factorization +! of a complete factorization ! ! A*Y = X, ! @@ -70,6 +73,11 @@ ! and block diagonal or replicated; in case c. the matrix A and its ! factors are distributed. ! +! This routine is used by mld_dsub_aply, to apply a 'base' block-Jacobi or +! Additive Schwarz (AS) preconditioner at any level of a multilevel preconditioner, +! or a block-Jacobi or LU or ILU solver at the coarsest level of a multilevel +! preconditioner. +! ! ! Arguments: ! diff --git a/mlprec/mld_zas_bld.f90 b/mlprec/mld_zas_bld.f90 index 2c41ca37..dc9a6251 100644 --- a/mlprec/mld_zas_bld.f90 +++ b/mlprec/mld_zas_bld.f90 @@ -41,7 +41,9 @@ ! ! This routine builds the Additive Schwarz (AS) preconditioner. ! If the preconditioner is the block-Jacobi one, the routine makes only a copy of -! the descriptor of the original matrix. +! the descriptor of the original matrix and then proceeds to call mld_fact_bld +! for LU or incomplete LU factorization of the diagonal blocks of the +! distributed matrix. ! ! ! Arguments: @@ -69,7 +71,6 @@ subroutine mld_zas_bld(a,desc_a,p,upd,info) Implicit None - ! Arguments ! Arguments type(psb_zspmat_type), intent(in), target :: a Type(psb_desc_type), Intent(in) :: desc_a diff --git a/mlprec/mld_zbaseprec_bld.f90 b/mlprec/mld_zbaseprec_bld.f90 index 1de7a9c1..e69dbbb5 100644 --- a/mlprec/mld_zbaseprec_bld.f90 +++ b/mlprec/mld_zbaseprec_bld.f90 @@ -154,32 +154,8 @@ subroutine mld_zbaseprc_bld(a,desc_a,p,info,upd) goto 9999 end if - case(mld_bjac_) - - call mld_check_def(p%iprcparm(mld_sub_ren_),'renumbering',& - & mld_renum_none_,is_legal_renum) - call mld_check_def(p%iprcparm(mld_sub_solve_),'fact',& - & mld_ilu_n_,is_legal_ml_fact) - - call psb_cdcpy(desc_a,p%desc_data,info) - if(info /= 0) then - info=4010 - ch_err='psb_cdcpy' - call psb_errpush(info,name,a_err=ch_err) - goto 9999 - end if - - ! Build the local part of the base preconditioner - call mld_fact_bld(a,p,iupd,info) - if(info /= 0) then - info=4010 - call psb_errpush(info,name,a_err='mld_fact_bld') - goto 9999 - end if - - - case(mld_as_) - ! Block Jacobi and additive Schwarz preconditioners/smoothers + case(mld_bjac_,mld_as_) + ! Additive Schwarz preconditioners/smoothers call mld_check_def(p%iprcparm(mld_n_ovr_),'overlap',& & 0,is_legal_n_ovr) @@ -202,7 +178,7 @@ subroutine mld_zbaseprc_bld(a,desc_a,p,info,upd) & write(debug_unit,*) me,' ',trim(name),& & ': Calling mld_as_bld' - ! Build the local part of the base preconditioner + ! Build the local part of the base preconditioner/smoother call mld_as_bld(a,desc_a,p,iupd,info) if(info /= 0) then info=4010 diff --git a/mlprec/mld_zfact_bld.f90 b/mlprec/mld_zfact_bld.f90 index d788fd0d..6c913c86 100644 --- a/mlprec/mld_zfact_bld.f90 +++ b/mlprec/mld_zfact_bld.f90 @@ -39,34 +39,35 @@ ! Subroutine: mld_zfact_bld ! Version: complex ! -! This routine computes an LU or incomplete LU factorization of the input -! matrix, according to the value of p%iprcparm(iprcparm(sub_solve_), +! This routine computes an LU or incomplete LU factorization of the diagonal blocks +! of a distributed matrix, according to the value of p%iprcparm(iprcparm(sub_solve_), ! set by the user through mld_dprecinit or mld_dprecset. -! It may also split the local matrix into its block-diagonal and +! It may also split the matrix into its block-diagonal and ! off block-diagonal parts, for the future application of multiple ! block-Jacobi sweeps. ! -! This routine is used by mld_dbaseprec_bld, to build a 'base' block-Jacobi or +! This routine is used by mld_as_bld, to build a 'base' block-Jacobi or ! Additive Schwarz (AS) preconditioner at any level of a multilevel preconditioner, ! or a block-Jacobi or LU or ILU solver at the coarsest level of a multilevel ! preconditioner. For the Additive Schwarz, it is called from mld_as_bld, ! which prepares the overlap descriptor and retrieves the remote rows into blck. ! -! More precisely, the routine performs one of the following tasks: +! More precisely, the routine performs one of the following tasks: ! -! 1. construction of a block-Jacobi preconditioner associated -! to a matrix A distributed among the processes (allowed at any level); +! 1. LU or ILU factorization of the diagonal blocks of the distributed matrix +! for the construction of a block-Jacobi or AS preconditioners +! (allowed at any level); ! ! 2. setup of block-Jacobi sweeps to compute an approximate solution of a ! linear system ! A*Y = X, ! distributed among the processes (allowed only at the coarsest level); ! -! 3. LU factorization of a linear system +! 3. LU factorization of the matrix of a linear system ! A*Y = X, ! distributed among the processes (allowed only at the coarsest level); ! -! 4. LU or incomplete LU factorization of a linear system +! 4. LU or incomplete LU factorization of the matrix of a linear system ! A*Y = X, ! replicated on the processes (allowed only at the coarsest level). ! @@ -84,7 +85,7 @@ ! Arguments: ! a - type(psb_zspmat_type), input. ! The sparse matrix structure containing the local part of the -! matrix to be preconditioned or factorized. +! distributed matrix. ! p - type(mld_zbaseprec_type), input/output. ! The 'base preconditioner' data structure containing the local ! part of the preconditioner or solver at the current level. @@ -99,9 +100,9 @@ ! the new preconditioner. ! blck - type(psb_zspmat_type), input, optional. ! The sparse matrix structure containing the remote rows of the -! matrix to be factorized, that have been retrieved by mld_as_bld +! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap -! greater than 0. If the overlap is 0 blck is empty. +! greater than 0. If the overlap is 0 blck is empty. ! subroutine mld_zfact_bld(a,p,upd,info,blck) @@ -291,8 +292,6 @@ subroutine mld_zfact_bld(a,p,upd,info,blck) ! No reordering of the local matrix is required ! case(0) - - ! ! In case of multiple block-Jacobi sweeps, clip into p%av(ap_nd_) ! the off block-diagonal part of the local extended matrix. The @@ -334,7 +333,6 @@ subroutine mld_zfact_bld(a,p,upd,info,blck) goto 9999 end if end if - ! ! Compute a factorization of the diagonal block of the local matrix, ! according to the choice made by the user by setting p%iprcparm(sub_solve_) diff --git a/mlprec/mld_zilu0_fact.f90 b/mlprec/mld_zilu0_fact.f90 index 62345172..021c317a 100644 --- a/mlprec/mld_zilu0_fact.f90 +++ b/mlprec/mld_zilu0_fact.f90 @@ -41,20 +41,20 @@ ! Contains: mld_zilu0_factint, ilu_copyin ! ! This routine computes either the ILU(0) or the MILU(0) factorization of the -! local part of the matrix stored into a. These factorizations are used to -! build the 'base preconditioner' (block-Jacobi preconditioner/solver, Additive -! Schwarz preconditioner) corresponding to a certain level of a multilevel +! diagonal blocks of a distributed matrix. These factorizations +! are used to build the 'base preconditioner' (block-Jacobi preconditioner/solver, +! Additive Schwarz preconditioner) corresponding to a given level of a multilevel ! preconditioner. ! ! Details on the above factorizations can be found in ! Y. Saad, Iterative Methods for Sparse Linear Systems, Second Edition, ! SIAM, 2003, Chapter 10. ! -! The local matrix to be factorized is stored into a and blck, as specified -! in the description of the arguments below. The storage format for both the -! L and U factors is CSR. The diagonal of the U factor is stored separately -! (actually, the inverse of the diagonal entries is stored; this is then -! managed in the solve stage associated to the ILU(0)/MILU(0) factorization). +! The local matrix is stored into a and blck, as specified in the description +! of the arguments below. The storage format for both the L and U factors is CSR. +! The diagonal of the U factor is stored separately (actually, the inverse of the +! diagonal entries is stored; this is then managed in the solve stage associated +! to the ILU(0)/MILU(0) factorization). ! ! The routine copies and factors "on the fly" from a and blck into l (L factor), ! u (U factor, except its diagonal) and d (diagonal of U). @@ -69,11 +69,11 @@ ! The MILU(0) factorization is computed if ialg = 2 (= mld_milu_n_); ! the ILU(0) factorization otherwise. ! a - type(psb_zspmat_type), input. -! The sparse matrix structure containing the local matrix to be -! factorized. Note that if the 'base' Additive Schwarz preconditioner +! The sparse matrix structure containing the local matrix. +! Note that if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered -! (see mld_bjac_bld), then a contains only the 'original' local part -! of the matrix to be factorized, i.e. the rows of the matrix held +! (see mld_as_bld), then a contains only the 'original' local part +! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! l - type(psb_zspmat_type), input/output. ! The L factor in the incomplete factorization. @@ -92,15 +92,16 @@ ! Error code. ! blck - type(psb_zspmat_type), input, optional, target. ! The sparse matrix structure containing the remote rows of the -! matrix to be factorized, that have been retrieved by mld_asmat_bld +! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered -! (see mld_bjac_bld), then blck is empty. +! (see mld_fact_bld), then blck is empty. ! subroutine mld_zilu0_fact(ialg,a,l,u,d,info,blck) use psb_base_mod use mld_inner_mod, mld_protect_name => mld_zilu0_fact + implicit none ! Arguments @@ -202,14 +203,15 @@ contains ! ! Subroutine: mld_zilu0_factint ! Version: complex - ! Note: internal subroutine of mld_zilu0_fact + ! Note: internal subroutine of mld_zilu0_fact. ! - ! This routine computes either the ILU(0) or the MILU(0) factorization of the - ! local part of the matrix stored into a. These factorizations are used to build - ! the 'base preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz - ! preconditioner) corresponding to a certain level of a multilevel preconditioner. + ! This routine computes either the ILU(0) or the MILU(0) factorization of the + ! diagonal blocks of a distributed matrix. + ! These factorizations are used to build the 'base preconditioner' + ! (block-Jacobi preconditioner/solver, Additive Schwarz + ! preconditioner) corresponding to a given level of a multilevel preconditioner. ! - ! The local matrix to be factorized is stored into a and b, as specified in the + ! The local matrix is stored into a and b, as specified in the ! description of the arguments below. The storage format for both the L and U ! factors is CSR. The diagonal of the U factor is stored separately (actually, ! the inverse of the diagonal entries is stored; this is then managed in the @@ -222,28 +224,29 @@ contains ! Arguments: ! ialg - integer, input. ! The type of incomplete factorization to be performed. - ! The MILU(0) factorization is computed if ialg = 2 (= mld_milu_n_); - ! the ILU(0) factorization otherwise. + ! The ILU(0) factorization is computed if ialg = 1 (= mld_ilu_n_), + ! the MILU(0) one if ialg = 2 (= mld_milu_n_); other values + ! are not allowed. ! m - integer, output. ! The total number of rows of the local matrix to be factorized, ! i.e. ma+mb. ! ma - integer, input ! The number of rows of the local submatrix stored into a. ! a - type(psb_zspmat_type), input. - ! The sparse matrix structure containing the local matrix to be - ! factorized. Note that, if the 'base' Additive Schwarz preconditioner + ! The sparse matrix structure containing the local matrix. + ! Note that, if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered - ! (see mld_bjac_bld), then a contains only the 'original' local part - ! of the matrix to be factorized, i.e. the rows of the matrix held + ! (see mld_fact_bld), then a contains only the 'original' local part + ! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! mb - integer, input. ! The number of rows of the local submatrix stored into b. ! b - type(psb_zspmat_type), input. ! The sparse matrix structure containing the remote rows of the - ! matrix to be factorized, that have been retrieved by mld_asmat_bld + ! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been - ! reordered (see mld_bjac_bld), then b does not contain any row. + ! reordered (see mld_fact_bld), then b does not contain any row. ! d - complex(kind(1.d0)), dimension(:), output. ! The inverse of the diagonal entries of the U factor in the ! incomplete factorization. @@ -484,7 +487,7 @@ contains ! copied into laspk, dia, uaspk row by row, through successive calls to ! ilu_copyin. ! - ! The routine is used by mld_zilu0_factin in the computation of the ILU(0)/MILU(0) + ! The routine is used by mld_zilu0_factint in the computation of the ILU(0)/MILU(0) ! factorization of a local sparse matrix. ! ! TODO: modify the routine to allow copying into output L and U that are diff --git a/mlprec/mld_zilu_bld.f90 b/mlprec/mld_zilu_bld.f90 index 47e32063..707c7f0d 100644 --- a/mlprec/mld_zilu_bld.f90 +++ b/mlprec/mld_zilu_bld.f90 @@ -39,9 +39,9 @@ ! Subroutine: mld_zilu_bld ! Version: complex ! -! This routine computes an incomplete LU (ILU) factorization of the local part -! of the matrix stored into a. This factorization is used to build the 'base -! preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz +! This routine computes an incomplete LU (ILU) factorization of the diagonal blocks +! of a distributed matrix. This factorization is used to build +! the 'base preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz ! preconditioner) corresponding to a certain level of a multilevel preconditioner. ! The following factorizations are available: ! - ILU(k), i.e. ILU factorization with fill-in level k, @@ -62,12 +62,12 @@ ! ! Arguments: ! a - type(psb_zspmat_type), input. -! The sparse matrix structure containing the local matrix to be -! factorized. Note that if p%iprcparm(mld_n_ovr_) > 0, i.e. the +! The sparse matrix structure containing the local matrix. +! Note that if p%iprcparm(mld_n_ovr_) > 0, i.e. the ! 'base' Additive Schwarz preconditioner has overlap greater than -! 0, and p%iprcparm(mld_sub_ren_) = 0, i.e. a reordering of the -! matrix has not been performed (see mld_bjac_bld), then a contains -! only the 'original' local part of the matrix to be factorized, +! 0, and p%iprcparm(mld_sub_ren_) = 0, i.e. a reordering of the +! matrix has not been performed (see mld_fact_bld), then a contains +! only the 'original' local part of the distributed matrix, ! i.e. the rows of the matrix held by the calling process according ! to the initial data distribution. ! p - type(mld_zbaseprc_type), input/output. @@ -81,10 +81,10 @@ ! Error code. ! blck - type(psb_zspmat_type), input, optional. ! The sparse matrix structure containing the remote rows of the -! matrix to be factorized, that have been retrieved by mld_asmat_bld +! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered -! (see mld_bjac_bld), then blck does not contain any row. +! (see mld_fact_bld), then blck does not contain any row. ! subroutine mld_zilu_bld(a,p,upd,info,blck) @@ -272,7 +272,6 @@ subroutine mld_zilu_bld(a,p,upd,info,blck) end if return - end subroutine mld_zilu_bld diff --git a/mlprec/mld_ziluk_fact.f90 b/mlprec/mld_ziluk_fact.f90 index 37040478..b5afd602 100644 --- a/mlprec/mld_ziluk_fact.f90 +++ b/mlprec/mld_ziluk_fact.f90 @@ -41,16 +41,16 @@ ! Contains: mld_ziluk_factint, iluk_copyin, iluk_fact, iluk_copyout ! ! This routine computes either the ILU(k) or the MILU(k) factorization of the -! local part of the matrix stored into a. These factorizations are used to -! build the 'base preconditioner' (block-Jacobi preconditioner/solver, Additive -! Schwarz preconditioner) corresponding to a certain level of a multilevel +! diagonal blocks of a distributed matrix. These factorizations are used to build +! the 'base preconditioner' (block-Jacobi preconditioner/solver, +! Additive Schwarz preconditioner) corresponding to a certain level of a multilevel ! preconditioner. ! ! Details on the above factorizations can be found in ! Y. Saad, Iterative Methods for Sparse Linear Systems, Second Edition, ! SIAM, 2003, Chapter 10. ! -! The local matrix to be factorized is stored into a and blck, as specified in +! The local matrix is stored into a and blck, as specified in ! the description of the arguments below. The storage format for both the L and ! U factors is CSR. The diagonal of the U factor is stored separately (actually, ! the inverse of the diagonal entries is stored; this is then managed in the solve @@ -62,14 +62,15 @@ ! The fill-in level k in ILU(k)/MILU(k). ! ialg - integer, input. ! The type of incomplete factorization to be performed. -! The MILU(k) factorization is computed if ialg = 2 (= mld_milu_n_); -! the ILU(k) factorization otherwise. +! The ILU(k) factorization is computed if ialg = 1 (= mld_ilu_n_); +! the MILU(k) one if ialg = 2 (= mld_milu_n_); other values are +! not allowed. ! a - type(psb_zspmat_type), input. -! The sparse matrix structure containing the local matrix to be -! factorized. Note that if the 'base' Additive Schwarz preconditioner +! The sparse matrix structure containing the local matrix. +! Note that if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered -! (see mld_bjac_bld), then a contains only the 'original' local part -! of the matrix to be factorized, i.e. the rows of the matrix held +! (see mld_fact_bld), then a contains only the 'original' local part +! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! l - type(psb_zspmat_type), input/output. ! The L factor in the incomplete factorization. @@ -88,15 +89,16 @@ ! Error code. ! blck - type(psb_zspmat_type), input, optional, target. ! The sparse matrix structure containing the remote rows of the -! matrix to be factorized, that have been retrieved by mld_asmat_bld +! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered -! (see mld_bjac_bld), then blck does not contain any row. +! (see mld_fact_bld), then blck does not contain any row. ! subroutine mld_ziluk_fact(fill_in,ialg,a,l,u,d,info,blck) use psb_base_mod use mld_inner_mod, mld_protect_name => mld_ziluk_fact + implicit none ! Arguments @@ -199,11 +201,11 @@ contains ! Note: internal subroutine of mld_ziluk_fact ! ! This routine computes either the ILU(k) or the MILU(k) factorization of the - ! local part of the matrix stored into a. These factorizations are used to build + ! diagonal blocks of a distributed matrix. These factorizations are used to build ! the 'base preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz ! preconditioner) corresponding to a certain level of a multilevel preconditioner. ! - ! The local matrix to be factorized is stored into a and b, as specified in the + ! The local matrix is stored into a and b, as specified in the ! description of the arguments below. The storage format for both the L and U ! factors is CSR. The diagonal of the U factor is stored separately (actually, ! the inverse of the diagonal entries is stored; this is then managed in the @@ -221,18 +223,18 @@ contains ! The total number of rows of the local matrix to be factorized, ! i.e. ma+mb. ! a - type(psb_zspmat_type), input. - ! The sparse matrix structure containing the local matrix to be - ! factorized. Note that, if the 'base' Additive Schwarz preconditioner + ! The sparse matrix structure containing the local matrix. + ! Note that, if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered - ! (see mld_bjac_bld), then a contains only the 'original' local part - ! of the matrix to be factorized, i.e. the rows of the matrix held + ! (see mld_fact_bld), then a contains only the 'original' local part + ! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! b - type(psb_zspmat_type), input. ! The sparse matrix structure containing the remote rows of the - ! matrix to be factorized, that have been retrieved by mld_asmat_bld + ! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered - ! (see mld_bjac_bld), then b does not contain any row. + ! (see mld_fact_bld), then b does not contain any row. ! d - complex(kind(1.d0)), dimension(:), output. ! The inverse of the diagonal entries of the U factor in the incomplete ! factorization. diff --git a/mlprec/mld_zilut_fact.f90 b/mlprec/mld_zilut_fact.f90 index f5b844c3..7d9b05d9 100644 --- a/mlprec/mld_zilut_fact.f90 +++ b/mlprec/mld_zilut_fact.f90 @@ -40,8 +40,8 @@ ! Version: real ! Contains: mld_zilut_factint, ilut_copyin, ilut_fact, ilut_copyout ! -! This routine computes the ILU(k,t) factorization of the local part of the -! matrix stored into a. This factorization is used to build the 'base +! This routine computes the ILU(k,t) factorization of the diagonal blocks of a +! distributed matrix. This factorization is used to build the 'base ! preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz ! preconditioner) corresponding to a certain level of a multilevel preconditioner. ! @@ -49,7 +49,7 @@ ! Y. Saad, Iterative Methods for Sparse Linear Systems, Second Edition, ! SIAM, 2003, Chapter 10. ! -! The local matrix to be factorized is stored into a and blck, as specified in +! The local matrix is stored into a and blck, as specified in ! the description of the arguments below. The storage format for both the L and ! U factors is CSR. The diagonal of the U factor is stored separately (actually, ! the inverse of the diagonal entries is stored; this is then managed in the solve @@ -62,11 +62,11 @@ ! thres - integer, input. ! The threshold t, i.e. the drop tolerance, in ILU(k,t). ! a - type(psb_zspmat_type), input. -! The sparse matrix structure containing the local matrix to be -! factorized. Note that if the 'base' Additive Schwarz preconditioner +! The sparse matrix structure containing the local matrix. +! Note that if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered -! (see mld_bjac_bld), then a contains only the 'original' local part -! of the matrix to be factorized, i.e. the rows of the matrix held +! (see mld_fact_bld), then a contains only the 'original' local part +! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! l - type(psb_zspmat_type), input/output. ! The L factor in the incomplete factorization. @@ -85,15 +85,16 @@ ! Error code. ! blck - type(psb_zspmat_type), input, optional, target. ! The sparse matrix structure containing the remote rows of the -! matrix to be factorized, that have been retrieved by mld_asmat_bld +! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered -! (see mld_bjac_bld), then blck does not contain any row. +! (see mld_fact_bld), then blck does not contain any row. ! subroutine mld_zilut_fact(fill_in,thres,a,l,u,d,info,blck) use psb_base_mod use mld_inner_mod, mld_protect_name => mld_zilut_fact + implicit none ! Arguments @@ -202,8 +203,8 @@ contains ! Version: real ! Note: internal subroutine of mld_zilut_fact ! - ! This routine computes the ILU(k,t) factorization of the local part of the - ! matrix stored into a. These factorization is used to build the 'base + ! This routine computes the ILU(k,t) factorization of the diagonal blocks of a + ! distributed matrix. This factorization is used to build the 'base ! preconditioner' (block-Jacobi preconditioner/solver, Additive Schwarz ! preconditioner) corresponding to a certain level of a multilevel preconditioner. ! @@ -223,18 +224,18 @@ contains ! The total number of rows of the local matrix to be factorized, ! i.e. ma+mb. ! a - type(psb_zspmat_type), input. - ! The sparse matrix structure containing the local matrix to be - ! factorized. Note that, if the 'base' Additive Schwarz preconditioner + ! The sparse matrix structure containing the local matrix. + ! Note that, if the 'base' Additive Schwarz preconditioner ! has overlap greater than 0 and the matrix has not been reordered - ! (see mld_bjac_bld), then a contains only the 'original' local part - ! of the matrix to be factorized, i.e. the rows of the matrix held + ! (see mld_fact_bld), then a contains only the 'original' local part + ! of the distributed matrix, i.e. the rows of the matrix held ! by the calling process according to the initial data distribution. ! b - type(psb_zspmat_type), input. ! The sparse matrix structure containing the remote rows of the - ! matrix to be factorized, that have been retrieved by mld_asmat_bld + ! distributed matrix, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0. If the overlap is 0 or the matrix has been reordered - ! (see mld_bjac_bld), then b does not contain any row. + ! (see mld_fact_bld), then b does not contain any row. ! d - complex(kind(1.d0)), dimension(:), output. ! The inverse of the diagonal entries of the U factor in the incomplete ! factorization. @@ -286,6 +287,7 @@ contains type(psb_int_heap) :: heap type(psb_zspmat_type) :: trw character(len=20), parameter :: name='mld_zilut_factint' + character(len=20) :: ch_err if (psb_get_errstatus() /= 0) return info = 0 @@ -785,7 +787,7 @@ contains ! - the array row is re-initialized for future use in mld_ilut_fact(see also ! ilut_copyin and ilut_fact). ! - ! This routine is used by mld_zilut_factint in the computation of the ILU(k,t) + ! This routine is used by mld_zilut_factint in the computation of the ILU(k,t) ! factorization of a local sparse matrix. ! ! diff --git a/mlprec/mld_zmlprec_aply.f90 b/mlprec/mld_zmlprec_aply.f90 index 237d9235..06efb0c5 100644 --- a/mlprec/mld_zmlprec_aply.f90 +++ b/mlprec/mld_zmlprec_aply.f90 @@ -49,15 +49,15 @@ ! - X and Y are vectors, ! - alpha and beta are scalars. ! -! For each level we have as many subdomains as processes (except for the coarsest +! For each level we have as many submatrices as processes (except for the coarsest ! level where we might have a replicated index space) and each process takes care -! of one subdomain. +! of one submatrix. ! ! The multilevel preconditioner M is regarded as an array of 'base preconditioners', ! each representing the part of the preconditioner associated to a certain level. ! For each level ilev, the base preconditioner K(ilev) is stored in baseprecv(ilev) ! and is associated to a matrix A(ilev), obtained by 'tranferring' the original -! matrix A (i.e. the matrix to be preconditioned) to level ilev, through smoothed +! matrix A (i.e. the matrix to be preconditioned) to the level ilev, through smoothed ! aggregation. ! ! The levels are numbered in increasing order starting from the finest one, i.e. @@ -157,8 +157,8 @@ ! Error code. ! ! Note that when the LU factorization of the matrix A(ilev) is computed instead of -! the ILU one, by using UMFPACK or SuperLU_dist, the corresponding L and U factors -! are stored in data structures provided by UMFPACK or SuperLU_dist and pointed by +! the ILU one, by using UMFPACK or SuperLU, the corresponding L and U factors +! are stored in data structures provided by UMFPACK or SuperLU and pointed by ! baseprecv(ilev)%iprcparm(mld_umf_ptr) or baseprecv(ilev)%iprcparm(mld_slu_ptr), ! respectively. ! @@ -285,7 +285,44 @@ subroutine mld_zmlprec_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info) return contains - + ! + ! Subroutine: add_ml_aply + ! Version: complex + ! Note: internal subroutine of mld_dmlprec_aply. + ! This routine computes + ! + ! Y = beta*Y + alpha*op(M^(-1))*X, + ! where + ! - M is a multilevel domain decomposition (Schwarz) preconditioner associated + ! to a certain matrix A and stored in the array baseprecv, + ! - op(M^(-1)) is M^(-1) or its (conjugate) transpose, according to + ! the value of trans, + ! - X and Y are vectors, + ! - alpha and beta are scalars. + ! + ! For each level we have as many submatrices as processes (except for the coarsest + ! level where we might have a replicated index space) and each process takes care + ! of one submatrix. + ! + ! The multilevel preconditioner M is regarded as an array of 'base preconditioners', + ! each representing the part of the preconditioner associated to a certain level. + ! For each level ilev, the base preconditioner K(ilev) is stored in baseprecv(ilev) + ! and is associated to a matrix A(ilev), obtained by 'tranferring' the original + ! matrix A (i.e. the matrix to be preconditioned) to the level ilev, through smoothed + ! aggregation. + ! + ! The levels are numbered in increasing order starting from the finest one, i.e. + ! level 1 is the finest level and A(1) is the matrix A. + ! This routine applies the multilevel preconditioner in an additive + ! way (additive through the levels and additive on the same level). + ! For details on the additive multilevel Schwarz preconditioner see + ! the Algorithm 3.1.1 in the book: + ! - B.F. Smith, P.E. Bjorstad & W.D. Gropp, + ! Domain decomposition: parallel multilevel methods for elliptic partial + ! differential equations, Cambridge University Press, 1996. + ! + ! For a detailed description of the arguments, see mld_dmlprec_aply. + ! subroutine add_ml_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info) implicit none ! Arguments @@ -329,39 +366,6 @@ contains call psb_errpush(4010,name,a_err='Allocate') goto 9999 end if - - ! - ! Additive multilevel - ! - ! 1. ! Apply the base preconditioner at level 1. - ! ! The sum over the subdomains is carried out in the - ! ! application of K(1). - ! X(1) = Xest - ! Y(1) = (K(1)^(-1))*X(1) - ! - ! 2. DO ilev=2,nlev - ! - ! ! Transfer X(ilev-1) to the next coarser level. - ! X(ilev) = AV(ilev; sm_pr_t_)*X(ilev-1) - ! - ! ! Apply the base preconditioner at the current level. - ! ! The sum over the subdomains is carried out in the - ! ! application of K(ilev). - ! Y(ilev) = (K(ilev)^(-1))*X(ilev) - ! - ! ENDDO - ! - ! 3. DO ilev=nlev-1,1,-1 - ! - ! ! Transfer Y(ilev+1) to the next finer level. - ! Y(ilev) = AV(ilev+1; sm_pr_)*Y(ilev+1) - ! - ! ENDDO - ! - ! 4. Yext = beta*Yext + alpha*Y(1) - ! - - ! ! STEP 1 ! @@ -385,7 +389,6 @@ contains call psb_errpush(4010,name,a_err='baseprec_aply') goto 9999 end if - ! ! STEP 2 ! @@ -524,9 +527,47 @@ contains return end subroutine add_ml_aply - - + ! + ! Subroutine: mlt_pre_ml_aply + ! Version: complex + ! Note: internal subroutine of mld_dmlprec_aply. + ! This routine computes + ! + ! Y = beta*Y + alpha*op(M^(-1))*X, + ! where + ! - M is a multilevel domain decomposition (Schwarz) preconditioner associated + ! to a certain matrix A and stored in the array baseprecv, + ! - op(M^(-1)) is M^(-1) or its (conjugate) transpose, according to + ! the value of trans, + ! - X and Y are vectors, + ! - alpha and beta are scalars. + ! + ! For each level we have as many submatrices as processes (except for the coarsest + ! level where we might have a replicated index space) and each process takes care + ! of one submatrix. + ! + ! The multilevel preconditioner M is regarded as an array of 'base preconditioners', + ! each representing the part of the preconditioner associated to a certain level. + ! For each level ilev, the base preconditioner K(ilev) is stored in baseprecv(ilev) + ! and is associated to a matrix A(ilev), obtained by 'tranferring' the original + ! matrix A (i.e. the matrix to be preconditioned) to the level ilev, through smoothed + ! aggregation. + ! + ! The levels are numbered in increasing order starting from the finest one, i.e. + ! level 1 is the finest level and A(1) is the matrix A. + ! + ! This routine applies the multilevel preconditioner in a hybrid way + ! (multiplicative through the levels and additive on the same level) + ! and pre-smoothing. + ! For details on pre-smoothed hybrid multiplicative multilevel Schwarz preconditioner, + ! see the Algorithm 3.2.1 in the book: + ! - B.F. Smith, P.E. Bjorstad & W.D. Gropp, + ! Domain decomposition: parallel multilevel methods for elliptic partial + ! differential equations, Cambridge University Press, 1996. + ! + ! For a detailed description of the arguments, see mld_dmlprec_aply. subroutine mlt_pre_ml_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info) + ! implicit none ! Arguments type(psb_desc_type),intent(in) :: desc_data @@ -570,43 +611,6 @@ contains goto 9999 end if - ! - ! Pre-smoothing - ! - ! 1. X(1) = Xext - ! - ! 2. ! Apply the base preconditioner at the finest level. - ! Y(1) = (K(1)^(-1))*X(1) - ! - ! 3. ! Compute the residual at the finest level. - ! TX(1) = X(1) - A(1)*Y(1) - ! - ! 4. DO ilev=2, nlev - ! - ! ! Transfer the residual to the current (coarser) level. - ! X(ilev) = AV(ilev; sm_pr_t_)*TX(ilev-1) - ! - ! ! Apply the base preconditioner at the current level. - ! ! The sum over the subdomains is carried out in the - ! ! application of K(ilev). - ! Y(ilev) = (K(ilev)^(-1))*X(ilev) - ! - ! ! Compute the residual at the current level (except at - ! ! the coarsest level). - ! IF (ilev < nlev) - ! TX(ilev) = (X(ilev)-A(ilev)*Y(ilev)) - ! - ! ENDDO - ! - ! 5. DO ilev=nlev-1,1,-1 - ! - ! ! Transfer Y(ilev+1) to the next finer level - ! Y(ilev) = Y(ilev) + AV(ilev+1; sm_pr_)*Y(ilev+1) - ! - ! ENDDO - ! - ! 6. Yext = beta*Yext + alpha*Y(1) - ! ! ! STEP 1 @@ -800,9 +804,46 @@ contains end if return end subroutine mlt_pre_ml_aply - - + ! + ! Subroutine: mlt_post_ml_aply + ! Version: complex + ! Note: internal subroutine of mld_dmlprec_aply. + ! This routine computes + ! + ! Y = beta*Y + alpha*op(M^(-1))*X, + ! where + ! - M is a multilevel domain decomposition (Schwarz) preconditioner associated + ! to a certain matrix A and stored in the array baseprecv, + ! - op(M^(-1)) is M^(-1) or its (conjugate) transpose, according to + ! the value of trans, + ! - X and Y are vectors, + ! - alpha and beta are scalars. + ! + ! For each level we have as many submatrices as processes (except for the coarsest + ! level where we might have a replicated index space) and each process takes care + ! of one submatrix. + ! + ! The multilevel preconditioner M is regarded as an array of 'base preconditioners', + ! each representing the part of the preconditioner associated to a certain level. + ! For each level ilev, the base preconditioner K(ilev) is stored in baseprecv(ilev) + ! and is associated to a matrix A(ilev), obtained by 'tranferring' the original + ! matrix A (i.e. the matrix to be preconditioned) to the level ilev, through smoothed + ! aggregation. + ! + ! The levels are numbered in increasing order starting from the finest one, i.e. + ! level 1 is the finest level and A(1) is the matrix A. + ! + ! This routine applies the multilevel preconditioner in a hybrid way + ! (multiplicative through the levels and additive on the same level) + ! and post-smoothing. + ! For details on hybrid multiplicative multilevel Schwarz preconditioners, see + ! - B.F. Smith, P.E. Bjorstad & W.D. Gropp, + ! Domain decomposition: parallel multilevel methods for elliptic partial + ! differential equations, Cambridge University Press, 1996. + ! + ! For a detailed description of the arguments, see mld_dmlprec_aply. subroutine mlt_post_ml_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info) + ! implicit none ! Arguments type(psb_desc_type),intent(in) :: desc_data @@ -845,37 +886,6 @@ contains call psb_errpush(4010,name,a_err='Allocate') goto 9999 end if - - ! - ! Post-smoothing - ! - ! 1. X(1) = Xext - ! - ! 2. DO ilev=2, nlev - ! - ! ! Transfer X(ilev-1) to the next coarser level. - ! X(ilev) = AV(ilev; sm_pr_t_)*X(ilev-1) - ! - ! ENDDO - ! - ! 3.! Apply the preconditioner at the coarsest level. - ! Y(nlev) = (K(nlev)^(-1))*X(nlev) - ! - ! 4. DO ilev=nlev-1,1,-1 - ! - ! ! Transfer Y(ilev+1) to the next finer level. - ! Y(ilev) = AV(ilev+1; sm_pr_)*Y(ilev+1) - ! - ! ! Compute the residual at the current level and apply to it the - ! ! base preconditioner. The sum over the subdomains is carried out - ! ! in the application of K(ilev). - ! Y(ilev) = Y(ilev) + (K(ilev)^(-1))*(X(ilev)-A(ilev)*Y(ilev)) - ! - ! ENDDO - ! - ! 5. Yext = beta*Yext + alpha*Y(1) - ! - ! ! STEP 1 ! @@ -1090,9 +1100,48 @@ contains end if return end subroutine mlt_post_ml_aply - - + ! + ! Subroutine: mlt_twoside_ml_aply + ! Version: complex + ! Note: internal subroutine of mld_dmlprec_aply. + ! This routine computes + ! + ! Y = beta*Y + alpha*op(M^(-1))*X, + ! where + ! - M is a multilevel domain decomposition (Schwarz) preconditioner associated + ! to a certain matrix A and stored in the array baseprecv, + ! - op(M^(-1)) is M^(-1) or its (conjugate) transpose, according to + ! the value of trans, + ! - X and Y are vectors, + ! - alpha and beta are scalars. + ! + ! For each level we have as many submatrices as processes (except for the coarsest + ! level where we might have a replicated index space) and each process takes care + ! of one submatrix. + ! + ! The multilevel preconditioner M is regarded as an array of 'base preconditioners', + ! each representing the part of the preconditioner associated to a certain level. + ! For each level ilev, the base preconditioner K(ilev) is stored in baseprecv(ilev) + ! and is associated to a matrix A(ilev), obtained by 'tranferring' the original + ! matrix A (i.e. the matrix to be preconditioned) to the level ilev, through smoothed + ! aggregation. + ! + ! The levels are numbered in increasing order starting from the finest one, i.e. + ! level 1 is the finest level and A(1) is the matrix A. + ! + ! This routine applies the multilevel preconditioner in a symmetrized hybrid way + ! (multiplicative through the levels and additive on the same level, with pre and + ! post-smoothing). + ! For details on the symmetrized hybrid multiplicative multilevel Schwarz + ! preconditioners, see the Algorithm 3.2.2 of the book: + ! - B.F. Smith, P.E. Bjorstad & W.D. Gropp, + ! Domain decomposition: parallel multilevel methods for elliptic partial + ! differential equations, Cambridge University Press, 1996. + ! + ! For a detailed description of the arguments, see mld_dmlprec_aply. + ! subroutine mlt_twoside_ml_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info) + ! implicit none ! Arguments type(psb_desc_type),intent(in) :: desc_data @@ -1135,49 +1184,6 @@ contains call psb_errpush(4010,name,a_err='Allocate') goto 9999 end if - - ! - ! Pre- and post-smoothing (symmetrized) - ! - ! 1. X(1) = Xext - ! - ! 2. ! Apply the base peconditioner at the finest level - ! Y(1) = (K(1)^(-1))*X(1) - ! - ! 3. ! Compute the residual at the finest level - ! TX(1) = X(1) - A(1)*Y(1) - ! - ! 4. DO ilev=2, nlev - ! - ! ! Transfer the residual to the current (coarser) level - ! X(ilev) = AV(ilev; sm_pr_t)*TX(ilev-1) - ! - ! ! Apply the base preconditioner at the current level. - ! ! The sum over the subdomains is carried out in the - ! ! application of K(ilev) - ! Y(ilev) = (K(ilev)^(-1))*X(ilev) - ! - ! ! Compute the residual at the current level - ! TX(ilev) = (X(ilev)-A(ilev)*Y(ilev)) - ! - ! ENDDO - ! - ! 5. DO ilev=NLEV-1,1,-1 - ! - ! ! Transfer Y(ilev+1) to the next finer level - ! Y(ilev) = Y(ilev) + AV(ilev+1; sm_pr_)*Y(ilev+1) - ! - ! ! Compute the residual at the current level and apply to it the - ! ! base preconditioner. The sum over the subdomains is carried out - ! ! in the application of K(ilev) - ! Y(ilev) = Y(ilev) + (K(ilev)**(-1))*(X(ilev)-A(ilev)*Y(ilev)) - ! - ! ENDDO - ! - ! 6. Yext = beta*Yext + alpha*Y(1) - ! - - ! ! STEP 1 ! ! Copy the input vector X diff --git a/mlprec/mld_zprec_aply.f90 b/mlprec/mld_zprec_aply.f90 index aa91cea2..1925a1ae 100644 --- a/mlprec/mld_zprec_aply.f90 +++ b/mlprec/mld_zprec_aply.f90 @@ -141,7 +141,8 @@ subroutine mld_zprec_aply(prec,x,y,desc_data,info,trans,work) endif ! If the original distribution has an overlap we should fix that. - call psb_ovrl(y,desc_data,info,update=psb_avg_) + call psb_halo(y,desc_data,info,data=psb_comm_mov_) + if (present(work)) then else diff --git a/mlprec/mld_zprecinit.f90 b/mlprec/mld_zprecinit.f90 index 668aa2fc..1b64b660 100644 --- a/mlprec/mld_zprecinit.f90 +++ b/mlprec/mld_zprecinit.f90 @@ -210,7 +210,7 @@ subroutine mld_zprecinit(p,ptype,info,nlev) p%baseprecv(ilev_)%iprcparm(mld_n_ovr_) = 0 p%baseprecv(ilev_)%iprcparm(mld_ml_type_) = mld_mult_ml_ p%baseprecv(ilev_)%iprcparm(mld_aggr_alg_) = mld_dec_aggr_ - p%baseprecv(ilev_)%iprcparm(mld_aggr_kind_) = mld_smooth_prol_ + p%baseprecv(ilev_)%iprcparm(mld_aggr_kind_) = mld_smooth_prol_ p%baseprecv(ilev_)%iprcparm(mld_coarse_mat_) = mld_distr_mat_ p%baseprecv(ilev_)%iprcparm(mld_smooth_pos_) = mld_post_smooth_ p%baseprecv(ilev_)%iprcparm(mld_aggr_eig_) = mld_max_norm_ diff --git a/mlprec/mld_zprecset.f90 b/mlprec/mld_zprecset.f90 index 9806f284..54376786 100644 --- a/mlprec/mld_zprecset.f90 +++ b/mlprec/mld_zprecset.f90 @@ -58,7 +58,7 @@ ! numbers, as reported in MLD2P4 user's guide. ! val - integer, input. ! The value of the parameter to be set. The list of allowed -! values is reported in MLD2P4 user's guide. +! values is reported in MLD2P4 user's guide. ! info - integer, output. ! Error code. ! ilev - integer, optional, input. @@ -275,6 +275,8 @@ subroutine mld_zprecsetc(p,what,string,info,ilev) use mld_prec_mod, mld_protect_name => mld_zprecsetc implicit none + + ! Arguments type(mld_zprec_type), intent(inout) :: p integer, intent(in) :: what character(len=*), intent(in) :: string @@ -524,6 +526,8 @@ subroutine mld_zprecsetd(p,what,val,info,ilev) use mld_prec_mod, mld_protect_name => mld_zprecsetd implicit none + + ! Arguments type(mld_zprec_type), intent(inout) :: p integer, intent(in) :: what real(kind(1.d0)), intent(in) :: val diff --git a/mlprec/mld_zsp_renum.f90 b/mlprec/mld_zsp_renum.f90 index 548bf372..179ce0eb 100644 --- a/mlprec/mld_zsp_renum.f90 +++ b/mlprec/mld_zsp_renum.f90 @@ -50,7 +50,7 @@ ! description of the arguments below. ! ! If required by the user (p%iprcparm(mld_sub_ren_) /= 0), the routine is -! used by mld_bjac_bld in building the block-Jacobi and Additive Schwarz +! used by mld_fact_bld in building the block-Jacobi and Additive Schwarz ! 'base preconditioners' corresponding to any level of a multilevel ! preconditioner. ! @@ -63,7 +63,7 @@ ! distribution. ! blck - type(psb_zspmat_type), input. ! The sparse matrix structure containing the remote rows of the -! matrix to be reordered, that have been retrieved by mld_asmat_bld +! matrix to be reordered, that have been retrieved by mld_as_bld ! to build an Additive Schwarz base preconditioner with overlap ! greater than 0.If the overlap is 0, then blck does not contain ! any row. diff --git a/mlprec/mld_zsub_aply.f90 b/mlprec/mld_zsub_aply.f90 index 751c7077..70945af8 100644 --- a/mlprec/mld_zsub_aply.f90 +++ b/mlprec/mld_zsub_aply.f90 @@ -81,29 +81,28 @@ ! K^(-1), alpha = 1 and beta = 0. ! ! The block-Jacobi preconditioner or solver and the L and U factors of the LU -! or ILU factorizations have been built by the routine mld_dbjac_bld and stored -! into the 'base preconditioner' data structure prec. See mld_dbjac_bld for more +! or ILU factorizations have been built by the routine mld_fact_bld and stored +! into the 'base preconditioner' data structure prec. See mld_fact_bld for more ! details. ! -! This routine is used by mld_dbaseprec_aply, to apply a 'base' block-Jacobi or +! This routine is used by mld_as_aply, to apply a 'base' block-Jacobi or ! Additive Schwarz (AS) preconditioner at any level of a multilevel preconditioner, ! or a block-Jacobi or LU or ILU solver at the coarsest level of a multilevel ! preconditioner. ! -! Inside mld_dbaseprec_aply, tasks 1, 3 and 4 may be selected if -! prec%iprcparm(smooth_sweeps_) = 1, while task 2 if prec%iprcparm(smooth_sweeps_) -! > 1. Furthermore, tasks 1, 2 and 3 may be performed if the matrix A is +! Tasks 1, 3 and 4 are selected when prec%iprcparm(smooth_sweeps_) = 1, +! while task 2 is selected when prec%iprcparm(smooth_sweeps_) > 1. Furthermore +! Tasks 1, 2 and 3 may be performed when the matrix A is ! distributed among the processes (prec%iprcparm(mld_coarse_mat_) = mld_distr_mat_), -! while task 4 may be performed if A is replicated on the processes +! while task 4 may be performed when A is replicated on the processes ! (prec%iprcparm(mld_coarse_mat_) = mld_repl_mat_). Note that the matrix A is ! distributed among the processes at each level of the multilevel preconditioner, ! except the coarsest one, where it may be either distributed or replicated on -! the processes. Furthermore, the tasks 2, 3 and 4 are performed only at the -! coarsest level. Note also that this routine manages implicitly the fact that +! the processes. Tasks 2, 3 and 4 are performed only at the coarsest level. +! Note also that this routine manages implicitly the fact that ! the matrix is distributed or replicated, i.e. it does not make any explicit ! reference to the value of prec%iprcparm(mld_coarse_mat_). ! -! ! Arguments: ! ! alpha - complex(kind(0.d0)), input. @@ -208,9 +207,7 @@ subroutine mld_zsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info) endif else if (prec%iprcparm(mld_smooth_sweeps_) > 1) then - ! - ! TASK 2 ! ! Apply prec%iprcparm(smooth_sweeps_) sweeps of a block-Jacobi solver ! to compute an approximate solution of a linear system. diff --git a/mlprec/mld_zsub_solve.f90 b/mlprec/mld_zsub_solve.f90 index 3adca02a..1a1f9d9e 100644 --- a/mlprec/mld_zsub_solve.f90 +++ b/mlprec/mld_zsub_solve.f90 @@ -44,11 +44,14 @@ ! Y = beta*Y + alpha*op(K^(-1))*X, ! ! where -! - K is a factored matrix, as specified below, +! - K is a factored matrix, as specified below, ! - op(K^(-1)) is K^(-1) or its transpose, according to the value of trans, ! - X and Y are vectors, ! - alpha and beta are scalars. ! +! Depending on K, alpha, beta (and on the communication descriptor desc_data +! - see the arguments below), the above computation may correspond to one of +! the following tasks: ! ! 1. Solution of a linear system with sparse factors LU generated by means ! of an incomplete factorization approximating @@ -58,7 +61,7 @@ ! they are also block-diagonal) or replicated. ! ! 2. Solution of a linear system with sparse factors LU generated by means -! of an complete factorization +! of a complete factorization ! ! A*Y = X, ! @@ -70,6 +73,11 @@ ! and block diagonal or replicated; in case c. the matrix A and its ! factors are distributed. ! +! This routine is used by mld_dsub_aply, to apply a 'base' block-Jacobi or +! Additive Schwarz (AS) preconditioner at any level of a multilevel preconditioner, +! or a block-Jacobi or LU or ILU solver at the coarsest level of a multilevel +! preconditioner. +! ! ! Arguments: !