From 212730c62d28219385c2a909410d6ddc1458cc16 Mon Sep 17 00:00:00 2001 From: Salvatore Filippone Date: Thu, 21 Sep 2017 15:51:58 +0100 Subject: [PATCH] Fixed application and description of 1lev precs. --- mlprec/impl/Makefile | 8 +- mlprec/impl/mld_cfile_prec_descr.f90 | 193 ++++++++++++++++++ mlprec/impl/mld_cprecaply.f90 | 125 ++++++------ mlprec/impl/mld_dfile_prec_descr.f90 | 193 ++++++++++++++++++ mlprec/impl/mld_dprecaply.f90 | 125 ++++++------ mlprec/impl/mld_sfile_prec_descr.f90 | 193 ++++++++++++++++++ mlprec/impl/mld_sprecaply.f90 | 125 ++++++------ mlprec/impl/mld_zfile_prec_descr.f90 | 193 ++++++++++++++++++ mlprec/impl/mld_zprecaply.f90 | 125 ++++++------ .../smoother/mld_c_base_smoother_descr.f90 | 2 +- .../smoother/mld_c_jac_smoother_descr.f90 | 14 +- .../smoother/mld_d_base_smoother_descr.f90 | 2 +- .../smoother/mld_d_jac_smoother_descr.f90 | 14 +- .../smoother/mld_s_base_smoother_descr.f90 | 2 +- .../smoother/mld_s_jac_smoother_descr.f90 | 14 +- .../smoother/mld_z_base_smoother_descr.f90 | 2 +- .../smoother/mld_z_jac_smoother_descr.f90 | 14 +- mlprec/mld_c_as_smoother.f90 | 2 +- mlprec/mld_c_prec_type.f90 | 142 +------------ mlprec/mld_d_as_smoother.f90 | 2 +- mlprec/mld_d_prec_type.f90 | 142 +------------ mlprec/mld_s_as_smoother.f90 | 2 +- mlprec/mld_s_prec_type.f90 | 142 +------------ mlprec/mld_z_as_smoother.f90 | 2 +- mlprec/mld_z_prec_type.f90 | 142 +------------ 25 files changed, 1100 insertions(+), 820 deletions(-) create mode 100644 mlprec/impl/mld_cfile_prec_descr.f90 create mode 100644 mlprec/impl/mld_dfile_prec_descr.f90 create mode 100644 mlprec/impl/mld_sfile_prec_descr.f90 create mode 100644 mlprec/impl/mld_zfile_prec_descr.f90 diff --git a/mlprec/impl/Makefile b/mlprec/impl/Makefile index d2ae7909..3f12d9c4 100644 --- a/mlprec/impl/Makefile +++ b/mlprec/impl/Makefile @@ -21,25 +21,25 @@ MPFOBJS=$(SMPFOBJS) $(DMPFOBJS) $(CMPFOBJS) $(ZMPFOBJS) MPCOBJS=mld_dslud_interface.o mld_zslud_interface.o -DINNEROBJS= mld_dmlprec_bld.o \ +DINNEROBJS= mld_dmlprec_bld.o mld_dfile_prec_descr.o \ mld_d_smoothers_bld.o mld_d_hierarchy_bld.o \ mld_dilu0_fact.o mld_diluk_fact.o mld_dilut_fact.o mld_daggrmap_bld.o \ mld_d_dec_map_bld.o mld_dmlprec_aply.o mld_daggrmat_asb.o \ $(DMPFOBJS) mld_d_extprol_bld.o mld_d_lev_aggrmap_bld.o mld_d_lev_aggrmat_asb.o -SINNEROBJS= mld_smlprec_bld.o \ +SINNEROBJS= mld_smlprec_bld.o mld_sfile_prec_descr.o \ mld_s_smoothers_bld.o mld_s_hierarchy_bld.o \ mld_silu0_fact.o mld_siluk_fact.o mld_silut_fact.o mld_saggrmap_bld.o \ mld_s_dec_map_bld.o mld_smlprec_aply.o mld_saggrmat_asb.o \ $(SMPFOBJS) mld_s_extprol_bld.o mld_s_lev_aggrmap_bld.o mld_s_lev_aggrmat_asb.o -ZINNEROBJS= mld_zmlprec_bld.o \ +ZINNEROBJS= mld_zmlprec_bld.o mld_zfile_prec_descr.o \ mld_z_smoothers_bld.o mld_z_hierarchy_bld.o \ mld_zilu0_fact.o mld_ziluk_fact.o mld_zilut_fact.o mld_zaggrmap_bld.o \ mld_z_dec_map_bld.o mld_zmlprec_aply.o mld_zaggrmat_asb.o \ $(ZMPFOBJS) mld_z_extprol_bld.o mld_z_lev_aggrmap_bld.o mld_z_lev_aggrmat_asb.o -CINNEROBJS= mld_cmlprec_bld.o \ +CINNEROBJS= mld_cmlprec_bld.o mld_cfile_prec_descr.o \ mld_c_smoothers_bld.o mld_c_hierarchy_bld.o \ mld_cilu0_fact.o mld_ciluk_fact.o mld_cilut_fact.o mld_caggrmap_bld.o \ mld_c_dec_map_bld.o mld_cmlprec_aply.o mld_caggrmat_asb.o \ diff --git a/mlprec/impl/mld_cfile_prec_descr.f90 b/mlprec/impl/mld_cfile_prec_descr.f90 new file mode 100644 index 00000000..5e90204a --- /dev/null +++ b/mlprec/impl/mld_cfile_prec_descr.f90 @@ -0,0 +1,193 @@ +! +! +! MLD2P4 version 2.1 +! MultiLevel Domain Decomposition Parallel Preconditioners Package +! based on PSBLAS (Parallel Sparse BLAS version 3.5) +! +! (C) Copyright 2008, 2010, 2012, 2015, 2017 +! +! Salvatore Filippone Cranfield University, UK +! Pasqua D'Ambra IAC-CNR, Naples, IT +! Daniela di Serafino University of Campania "L. Vanvitelli", Caserta, IT +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the MLD2P4 group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE MLD2P4 GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: mld_dfile_prec_descr.f90 +! +! +! Subroutine: mld_file_prec_descr +! Version: complex +! +! This routine prints a description of the preconditioner to the standard +! output or to a file. It must be called after the preconditioner has been +! built by mld_precbld. +! +! Arguments: +! p - type(mld_Tprec_type), input. +! The preconditioner data structure to be printed out. +! info - integer, output. +! error code. +! iout - integer, input, optional. +! The id of the file where the preconditioner description +! will be printed. If iout is not present, then the standard +! output is condidered. +! root - integer, input, optional. +! The id of the process printing the message; -1 acts as a wildcard. +! Default is psb_root_ +! +subroutine mld_cfile_prec_descr(prec,iout,root) + use psb_base_mod + use mld_c_prec_mod, mld_protect_name => mld_cfile_prec_descr + use mld_c_inner_mod + use mld_c_gs_solver + + implicit none + ! Arguments + class(mld_cprec_type), intent(in) :: prec + integer(psb_ipk_), intent(in), optional :: iout + integer(psb_ipk_), intent(in), optional :: root + + ! Local variables + integer(psb_ipk_) :: ilev, nlev, ilmin, info, nswps + integer(psb_ipk_) :: ictxt, me, np + logical :: is_symgs + character(len=20), parameter :: name='mld_file_prec_descr' + integer(psb_ipk_) :: iout_ + integer(psb_ipk_) :: root_ + + info = psb_success_ + if (present(iout)) then + iout_ = iout + else + iout_ = psb_out_unit + end if + if (iout_ < 0) iout_ = psb_out_unit + + ictxt = prec%ictxt + + if (allocated(prec%precv)) then + + call psb_info(ictxt,me,np) + if (present(root)) then + root_ = root + else + root_ = psb_root_ + end if + if (root_ == -1) root_ = me + + ! + ! The preconditioner description is printed by processor psb_root_. + ! This agrees with the fact that all the parameters defining the + ! preconditioner have the same values on all the procs (this is + ! ensured by mld_precbld). + ! + if (me == root_) then + nlev = size(prec%precv) + do ilev = 1, nlev + if (.not.allocated(prec%precv(ilev)%sm)) then + info = 3111 + write(iout_,*) ' ',name,& + & ': error: inconsistent MLPREC part, should call MLD_PRECINIT' + return + endif + end do + + write(iout_,*) + write(iout_,'(a)') 'Preconditioner description' + + if (nlev == 1) then + ! + ! Here we have a gigantic kludge just to handle Symmetrized Gauss-Seidel. + ! Will need rethinking... + ! + if (allocated(prec%precv(1)%sm2a)) then + is_symgs = .false. + select type(sv2 => prec%precv(1)%sm2a%sv) + class is (mld_c_bwgs_solver_type) + select type(sv1 => prec%precv(1)%sm%sv) + class is (mld_c_gs_solver_type) + is_symgs = .true. + end select + end select + if (is_symgs) then + write(iout_,*) ' Forward-Backward (symmetrized) Hybrid Gauss-Seidel' + else + call prec%precv(1)%sm%descr(info,iout=iout_) + call prec%precv(1)%sm2a%descr(info,iout=iout_) + end if + + else + call prec%precv(1)%sm%descr(info,iout=iout_) + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + end if + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + if (nswps > 1) write(iout_,*) ' Number of sweeps : ',nswps + write(iout_,*) + + else if (nlev > 1) then + ! + ! Print description of base preconditioner + ! + write(iout_,*) 'Multilevel Preconditioner' + write(iout_,*) 'Outer sweeps:',prec%outer_sweeps + write(iout_,*) + write(iout_,*) 'Smoother: ' + call prec%precv(1)%sm%descr(info,iout=iout_) + if (allocated(prec%precv(1)%sm2a)) then + write(iout_,*) 'Post smoother:' + call prec%precv(1)%sm2a%descr(info,iout=iout_) + end if + ! + ! Print multilevel details + ! + write(iout_,*) + write(iout_,*) 'Multilevel hierarchy: ' + write(iout_,*) ' Number of levels : ',nlev + write(iout_,*) ' Operator complexity: ',prec%get_complexity() + ilmin = 2 + if (nlev == 2) ilmin=1 + do ilev=ilmin,nlev + call prec%precv(ilev)%descr(ilev,nlev,ilmin,info,iout=iout_) + end do + write(iout_,*) + + else + write(iout_,*) trim(name), & + & ': invalid preconditioner array size ?',nlev + info = -2 + return + + end if + end if + + else + write(iout_,*) trim(name), & + & ': Error: no base preconditioner available, something is wrong!' + info = -2 + return + endif + +end subroutine mld_cfile_prec_descr diff --git a/mlprec/impl/mld_cprecaply.f90 b/mlprec/impl/mld_cprecaply.f90 index 81c75254..804b051e 100644 --- a/mlprec/impl/mld_cprecaply.f90 +++ b/mlprec/impl/mld_cprecaply.f90 @@ -89,10 +89,10 @@ subroutine mld_cprecaply(prec,x,y,desc_data,info,trans,work) ! Local variables character :: trans_ complex(psb_spk_), pointer :: work_(:) - complex(psb_spk_), allocatable :: ww(:) + complex(psb_spk_), allocatable :: w1(:), w2(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_cprecaply' @@ -143,42 +143,44 @@ subroutine mld_cprecaply(prec,x,y,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! - ! - ! Number of levels = 1: apply the base preconditioner - ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. ! Will need some rethinking. ! - allocate(ww(size(x))) - + call psb_geasb(w1,desc_data,info,scratch=.true.) + call psb_geasb(w2,desc_data,info,scratch=.true.) + + call psb_geaxpby(cone,x,czero,w1,desc_data,info) select case(trans_) case ('N') - call prec%precv(1)%sm%apply(cone,x,czero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(cone,ww,czero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(cone,w1,czero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(cone,w2,czero,w1,desc_data,trans_,& + & ione, work_,info) + end do + case('T','C') - call prec%precv(1)%sm2a%apply(cone,x,czero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(cone,ww,czero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(cone,w1,czero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(cone,w2,czero,w1,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') goto 9999 end select - deallocate(ww) - + call psb_geaxpby(cone,w1,czero,y,desc_data,info) + call psb_gefree(w1,desc_data,info) + call psb_gefree(w2,desc_data,info) + else call prec%precv(1)%sm%apply(cone,x,czero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post), & - & work_,info) + & nswps, work_,info) end if else info = psb_err_from_subroutine_ai_ @@ -255,7 +257,7 @@ subroutine mld_cprecaply1(prec,x,desc_data,info,trans) ! Local variables integer(psb_ipk_) :: ictxt,np,me integer(psb_ipk_) :: err_act - complex(psb_spk_), pointer :: WW(:), w1(:) + complex(psb_spk_), pointer :: ww(:), w1(:) character(len=20) :: name name='mld_cprecaply1' @@ -282,7 +284,7 @@ subroutine mld_cprecaply1(prec,x,desc_data,info,trans) end if x(:) = ww(:) - deallocate(ww,W1,stat=info) + deallocate(ww,w1,stat=info) if (info /= psb_success_) then info = psb_err_alloc_dealloc_ call psb_errpush(info,name) @@ -319,7 +321,7 @@ subroutine mld_cprecaply2_vect(prec,x,y,desc_data,info,trans,work) character :: trans_ complex(psb_spk_), pointer :: work_(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_cprecaply' @@ -370,42 +372,47 @@ subroutine mld_cprecaply2_vect(prec,x,y,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. ! Will need some rethinking. ! twoside: block - type(psb_c_vect_type) :: ww - call psb_geasb(ww,desc_data,info,mold=x%v,scratch=.true.) - + type(psb_c_vect_type) :: w1,w2 + call psb_geasb(w1,desc_data,info,mold=x%v,scratch=.true.) + call psb_geasb(w2,desc_data,info,mold=x%v,scratch=.true.) + call psb_geaxpby(cone,x,czero,w1,desc_data,info) select case(trans_) case ('N') - call prec%precv(1)%sm%apply(cone,x,czero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(cone,ww,czero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(cone,w1,czero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(cone,w2,czero,w1,desc_data,trans_,& + & ione, work_,info) + end do + case('T','C') - call prec%precv(1)%sm2a%apply(cone,x,czero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(cone,ww,czero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(cone,w1,czero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(cone,w2,czero,w1,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') goto 9999 end select - call psb_gefree(ww,desc_data,info) + call psb_geaxpby(cone,w1,czero,y,desc_data,info) + call psb_gefree(w1,desc_data,info) + call psb_gefree(w2,desc_data,info) end block twoside else call prec%precv(1)%sm%apply(cone,x,czero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + & nswps,work_,info) end if else @@ -454,7 +461,7 @@ subroutine mld_cprecaply1_vect(prec,x,desc_data,info,trans,work) type(psb_c_vect_type) :: ww complex(psb_spk_), pointer :: work_(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_cprecaply' @@ -506,6 +513,7 @@ subroutine mld_cprecaply1_vect(prec,x,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. @@ -513,19 +521,19 @@ subroutine mld_cprecaply1_vect(prec,x,desc_data,info,trans,work) ! select case(trans_) case ('N') - call prec%precv(1)%sm%apply(cone,x,czero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(cone,ww,czero,x,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(cone,x,czero,ww,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(cone,ww,czero,x,desc_data,trans_,& + & ione, work_,info) + end do case('T','C') - call prec%precv(1)%sm2a%apply(cone,x,czero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(cone,ww,czero,x,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(cone,x,czero,ww,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(cone,ww,czero,x,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') @@ -534,8 +542,7 @@ subroutine mld_cprecaply1_vect(prec,x,desc_data,info,trans,work) else call prec%precv(1)%sm%apply(cone,x,czero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + & nswps, work_,info) if (info == 0) call psb_geaxpby(cone,ww,czero,x,desc_data,info) end if else diff --git a/mlprec/impl/mld_dfile_prec_descr.f90 b/mlprec/impl/mld_dfile_prec_descr.f90 new file mode 100644 index 00000000..c1cf6005 --- /dev/null +++ b/mlprec/impl/mld_dfile_prec_descr.f90 @@ -0,0 +1,193 @@ +! +! +! MLD2P4 version 2.1 +! MultiLevel Domain Decomposition Parallel Preconditioners Package +! based on PSBLAS (Parallel Sparse BLAS version 3.5) +! +! (C) Copyright 2008, 2010, 2012, 2015, 2017 +! +! Salvatore Filippone Cranfield University, UK +! Pasqua D'Ambra IAC-CNR, Naples, IT +! Daniela di Serafino University of Campania "L. Vanvitelli", Caserta, IT +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the MLD2P4 group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE MLD2P4 GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: mld_dfile_prec_descr.f90 +! +! +! Subroutine: mld_file_prec_descr +! Version: real +! +! This routine prints a description of the preconditioner to the standard +! output or to a file. It must be called after the preconditioner has been +! built by mld_precbld. +! +! Arguments: +! p - type(mld_Tprec_type), input. +! The preconditioner data structure to be printed out. +! info - integer, output. +! error code. +! iout - integer, input, optional. +! The id of the file where the preconditioner description +! will be printed. If iout is not present, then the standard +! output is condidered. +! root - integer, input, optional. +! The id of the process printing the message; -1 acts as a wildcard. +! Default is psb_root_ +! +subroutine mld_dfile_prec_descr(prec,iout,root) + use psb_base_mod + use mld_d_prec_mod, mld_protect_name => mld_dfile_prec_descr + use mld_d_inner_mod + use mld_d_gs_solver + + implicit none + ! Arguments + class(mld_dprec_type), intent(in) :: prec + integer(psb_ipk_), intent(in), optional :: iout + integer(psb_ipk_), intent(in), optional :: root + + ! Local variables + integer(psb_ipk_) :: ilev, nlev, ilmin, info, nswps + integer(psb_ipk_) :: ictxt, me, np + logical :: is_symgs + character(len=20), parameter :: name='mld_file_prec_descr' + integer(psb_ipk_) :: iout_ + integer(psb_ipk_) :: root_ + + info = psb_success_ + if (present(iout)) then + iout_ = iout + else + iout_ = psb_out_unit + end if + if (iout_ < 0) iout_ = psb_out_unit + + ictxt = prec%ictxt + + if (allocated(prec%precv)) then + + call psb_info(ictxt,me,np) + if (present(root)) then + root_ = root + else + root_ = psb_root_ + end if + if (root_ == -1) root_ = me + + ! + ! The preconditioner description is printed by processor psb_root_. + ! This agrees with the fact that all the parameters defining the + ! preconditioner have the same values on all the procs (this is + ! ensured by mld_precbld). + ! + if (me == root_) then + nlev = size(prec%precv) + do ilev = 1, nlev + if (.not.allocated(prec%precv(ilev)%sm)) then + info = 3111 + write(iout_,*) ' ',name,& + & ': error: inconsistent MLPREC part, should call MLD_PRECINIT' + return + endif + end do + + write(iout_,*) + write(iout_,'(a)') 'Preconditioner description' + + if (nlev == 1) then + ! + ! Here we have a gigantic kludge just to handle Symmetrized Gauss-Seidel. + ! Will need rethinking... + ! + if (allocated(prec%precv(1)%sm2a)) then + is_symgs = .false. + select type(sv2 => prec%precv(1)%sm2a%sv) + class is (mld_d_bwgs_solver_type) + select type(sv1 => prec%precv(1)%sm%sv) + class is (mld_d_gs_solver_type) + is_symgs = .true. + end select + end select + if (is_symgs) then + write(iout_,*) ' Forward-Backward (symmetrized) Hybrid Gauss-Seidel' + else + call prec%precv(1)%sm%descr(info,iout=iout_) + call prec%precv(1)%sm2a%descr(info,iout=iout_) + end if + + else + call prec%precv(1)%sm%descr(info,iout=iout_) + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + end if + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + if (nswps > 1) write(iout_,*) ' Number of sweeps : ',nswps + write(iout_,*) + + else if (nlev > 1) then + ! + ! Print description of base preconditioner + ! + write(iout_,*) 'Multilevel Preconditioner' + write(iout_,*) 'Outer sweeps:',prec%outer_sweeps + write(iout_,*) + write(iout_,*) 'Smoother: ' + call prec%precv(1)%sm%descr(info,iout=iout_) + if (allocated(prec%precv(1)%sm2a)) then + write(iout_,*) 'Post smoother:' + call prec%precv(1)%sm2a%descr(info,iout=iout_) + end if + ! + ! Print multilevel details + ! + write(iout_,*) + write(iout_,*) 'Multilevel hierarchy: ' + write(iout_,*) ' Number of levels : ',nlev + write(iout_,*) ' Operator complexity: ',prec%get_complexity() + ilmin = 2 + if (nlev == 2) ilmin=1 + do ilev=ilmin,nlev + call prec%precv(ilev)%descr(ilev,nlev,ilmin,info,iout=iout_) + end do + write(iout_,*) + + else + write(iout_,*) trim(name), & + & ': invalid preconditioner array size ?',nlev + info = -2 + return + + end if + end if + + else + write(iout_,*) trim(name), & + & ': Error: no base preconditioner available, something is wrong!' + info = -2 + return + endif + +end subroutine mld_dfile_prec_descr diff --git a/mlprec/impl/mld_dprecaply.f90 b/mlprec/impl/mld_dprecaply.f90 index cdec46b6..41a93449 100644 --- a/mlprec/impl/mld_dprecaply.f90 +++ b/mlprec/impl/mld_dprecaply.f90 @@ -89,10 +89,10 @@ subroutine mld_dprecaply(prec,x,y,desc_data,info,trans,work) ! Local variables character :: trans_ real(psb_dpk_), pointer :: work_(:) - real(psb_dpk_), allocatable :: ww(:) + real(psb_dpk_), allocatable :: w1(:), w2(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_dprecaply' @@ -143,42 +143,44 @@ subroutine mld_dprecaply(prec,x,y,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! - ! - ! Number of levels = 1: apply the base preconditioner - ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. ! Will need some rethinking. ! - allocate(ww(size(x))) - + call psb_geasb(w1,desc_data,info,scratch=.true.) + call psb_geasb(w2,desc_data,info,scratch=.true.) + + call psb_geaxpby(done,x,dzero,w1,desc_data,info) select case(trans_) case ('N') - call prec%precv(1)%sm%apply(done,x,dzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(done,ww,dzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(done,w1,dzero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(done,w2,dzero,w1,desc_data,trans_,& + & ione, work_,info) + end do + case('T','C') - call prec%precv(1)%sm2a%apply(done,x,dzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(done,ww,dzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(done,w1,dzero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(done,w2,dzero,w1,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') goto 9999 end select - deallocate(ww) - + call psb_geaxpby(done,w1,dzero,y,desc_data,info) + call psb_gefree(w1,desc_data,info) + call psb_gefree(w2,desc_data,info) + else call prec%precv(1)%sm%apply(done,x,dzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post), & - & work_,info) + & nswps, work_,info) end if else info = psb_err_from_subroutine_ai_ @@ -255,7 +257,7 @@ subroutine mld_dprecaply1(prec,x,desc_data,info,trans) ! Local variables integer(psb_ipk_) :: ictxt,np,me integer(psb_ipk_) :: err_act - real(psb_dpk_), pointer :: WW(:), w1(:) + real(psb_dpk_), pointer :: ww(:), w1(:) character(len=20) :: name name='mld_dprecaply1' @@ -282,7 +284,7 @@ subroutine mld_dprecaply1(prec,x,desc_data,info,trans) end if x(:) = ww(:) - deallocate(ww,W1,stat=info) + deallocate(ww,w1,stat=info) if (info /= psb_success_) then info = psb_err_alloc_dealloc_ call psb_errpush(info,name) @@ -319,7 +321,7 @@ subroutine mld_dprecaply2_vect(prec,x,y,desc_data,info,trans,work) character :: trans_ real(psb_dpk_), pointer :: work_(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_dprecaply' @@ -370,42 +372,47 @@ subroutine mld_dprecaply2_vect(prec,x,y,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. ! Will need some rethinking. ! twoside: block - type(psb_d_vect_type) :: ww - call psb_geasb(ww,desc_data,info,mold=x%v,scratch=.true.) - + type(psb_d_vect_type) :: w1,w2 + call psb_geasb(w1,desc_data,info,mold=x%v,scratch=.true.) + call psb_geasb(w2,desc_data,info,mold=x%v,scratch=.true.) + call psb_geaxpby(done,x,dzero,w1,desc_data,info) select case(trans_) case ('N') - call prec%precv(1)%sm%apply(done,x,dzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(done,ww,dzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(done,w1,dzero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(done,w2,dzero,w1,desc_data,trans_,& + & ione, work_,info) + end do + case('T','C') - call prec%precv(1)%sm2a%apply(done,x,dzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(done,ww,dzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(done,w1,dzero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(done,w2,dzero,w1,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') goto 9999 end select - call psb_gefree(ww,desc_data,info) + call psb_geaxpby(done,w1,dzero,y,desc_data,info) + call psb_gefree(w1,desc_data,info) + call psb_gefree(w2,desc_data,info) end block twoside else call prec%precv(1)%sm%apply(done,x,dzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + & nswps,work_,info) end if else @@ -454,7 +461,7 @@ subroutine mld_dprecaply1_vect(prec,x,desc_data,info,trans,work) type(psb_d_vect_type) :: ww real(psb_dpk_), pointer :: work_(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_dprecaply' @@ -506,6 +513,7 @@ subroutine mld_dprecaply1_vect(prec,x,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. @@ -513,19 +521,19 @@ subroutine mld_dprecaply1_vect(prec,x,desc_data,info,trans,work) ! select case(trans_) case ('N') - call prec%precv(1)%sm%apply(done,x,dzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(done,ww,dzero,x,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(done,x,dzero,ww,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(done,ww,dzero,x,desc_data,trans_,& + & ione, work_,info) + end do case('T','C') - call prec%precv(1)%sm2a%apply(done,x,dzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(done,ww,dzero,x,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(done,x,dzero,ww,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(done,ww,dzero,x,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') @@ -534,8 +542,7 @@ subroutine mld_dprecaply1_vect(prec,x,desc_data,info,trans,work) else call prec%precv(1)%sm%apply(done,x,dzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + & nswps, work_,info) if (info == 0) call psb_geaxpby(done,ww,dzero,x,desc_data,info) end if else diff --git a/mlprec/impl/mld_sfile_prec_descr.f90 b/mlprec/impl/mld_sfile_prec_descr.f90 new file mode 100644 index 00000000..39720497 --- /dev/null +++ b/mlprec/impl/mld_sfile_prec_descr.f90 @@ -0,0 +1,193 @@ +! +! +! MLD2P4 version 2.1 +! MultiLevel Domain Decomposition Parallel Preconditioners Package +! based on PSBLAS (Parallel Sparse BLAS version 3.5) +! +! (C) Copyright 2008, 2010, 2012, 2015, 2017 +! +! Salvatore Filippone Cranfield University, UK +! Pasqua D'Ambra IAC-CNR, Naples, IT +! Daniela di Serafino University of Campania "L. Vanvitelli", Caserta, IT +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the MLD2P4 group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE MLD2P4 GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: mld_dfile_prec_descr.f90 +! +! +! Subroutine: mld_file_prec_descr +! Version: real +! +! This routine prints a description of the preconditioner to the standard +! output or to a file. It must be called after the preconditioner has been +! built by mld_precbld. +! +! Arguments: +! p - type(mld_Tprec_type), input. +! The preconditioner data structure to be printed out. +! info - integer, output. +! error code. +! iout - integer, input, optional. +! The id of the file where the preconditioner description +! will be printed. If iout is not present, then the standard +! output is condidered. +! root - integer, input, optional. +! The id of the process printing the message; -1 acts as a wildcard. +! Default is psb_root_ +! +subroutine mld_sfile_prec_descr(prec,iout,root) + use psb_base_mod + use mld_s_prec_mod, mld_protect_name => mld_sfile_prec_descr + use mld_s_inner_mod + use mld_s_gs_solver + + implicit none + ! Arguments + class(mld_sprec_type), intent(in) :: prec + integer(psb_ipk_), intent(in), optional :: iout + integer(psb_ipk_), intent(in), optional :: root + + ! Local variables + integer(psb_ipk_) :: ilev, nlev, ilmin, info, nswps + integer(psb_ipk_) :: ictxt, me, np + logical :: is_symgs + character(len=20), parameter :: name='mld_file_prec_descr' + integer(psb_ipk_) :: iout_ + integer(psb_ipk_) :: root_ + + info = psb_success_ + if (present(iout)) then + iout_ = iout + else + iout_ = psb_out_unit + end if + if (iout_ < 0) iout_ = psb_out_unit + + ictxt = prec%ictxt + + if (allocated(prec%precv)) then + + call psb_info(ictxt,me,np) + if (present(root)) then + root_ = root + else + root_ = psb_root_ + end if + if (root_ == -1) root_ = me + + ! + ! The preconditioner description is printed by processor psb_root_. + ! This agrees with the fact that all the parameters defining the + ! preconditioner have the same values on all the procs (this is + ! ensured by mld_precbld). + ! + if (me == root_) then + nlev = size(prec%precv) + do ilev = 1, nlev + if (.not.allocated(prec%precv(ilev)%sm)) then + info = 3111 + write(iout_,*) ' ',name,& + & ': error: inconsistent MLPREC part, should call MLD_PRECINIT' + return + endif + end do + + write(iout_,*) + write(iout_,'(a)') 'Preconditioner description' + + if (nlev == 1) then + ! + ! Here we have a gigantic kludge just to handle Symmetrized Gauss-Seidel. + ! Will need rethinking... + ! + if (allocated(prec%precv(1)%sm2a)) then + is_symgs = .false. + select type(sv2 => prec%precv(1)%sm2a%sv) + class is (mld_s_bwgs_solver_type) + select type(sv1 => prec%precv(1)%sm%sv) + class is (mld_s_gs_solver_type) + is_symgs = .true. + end select + end select + if (is_symgs) then + write(iout_,*) ' Forward-Backward (symmetrized) Hybrid Gauss-Seidel' + else + call prec%precv(1)%sm%descr(info,iout=iout_) + call prec%precv(1)%sm2a%descr(info,iout=iout_) + end if + + else + call prec%precv(1)%sm%descr(info,iout=iout_) + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + end if + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + if (nswps > 1) write(iout_,*) ' Number of sweeps : ',nswps + write(iout_,*) + + else if (nlev > 1) then + ! + ! Print description of base preconditioner + ! + write(iout_,*) 'Multilevel Preconditioner' + write(iout_,*) 'Outer sweeps:',prec%outer_sweeps + write(iout_,*) + write(iout_,*) 'Smoother: ' + call prec%precv(1)%sm%descr(info,iout=iout_) + if (allocated(prec%precv(1)%sm2a)) then + write(iout_,*) 'Post smoother:' + call prec%precv(1)%sm2a%descr(info,iout=iout_) + end if + ! + ! Print multilevel details + ! + write(iout_,*) + write(iout_,*) 'Multilevel hierarchy: ' + write(iout_,*) ' Number of levels : ',nlev + write(iout_,*) ' Operator complexity: ',prec%get_complexity() + ilmin = 2 + if (nlev == 2) ilmin=1 + do ilev=ilmin,nlev + call prec%precv(ilev)%descr(ilev,nlev,ilmin,info,iout=iout_) + end do + write(iout_,*) + + else + write(iout_,*) trim(name), & + & ': invalid preconditioner array size ?',nlev + info = -2 + return + + end if + end if + + else + write(iout_,*) trim(name), & + & ': Error: no base preconditioner available, something is wrong!' + info = -2 + return + endif + +end subroutine mld_sfile_prec_descr diff --git a/mlprec/impl/mld_sprecaply.f90 b/mlprec/impl/mld_sprecaply.f90 index a9c2da27..487cb6f8 100644 --- a/mlprec/impl/mld_sprecaply.f90 +++ b/mlprec/impl/mld_sprecaply.f90 @@ -89,10 +89,10 @@ subroutine mld_sprecaply(prec,x,y,desc_data,info,trans,work) ! Local variables character :: trans_ real(psb_spk_), pointer :: work_(:) - real(psb_spk_), allocatable :: ww(:) + real(psb_spk_), allocatable :: w1(:), w2(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_sprecaply' @@ -143,42 +143,44 @@ subroutine mld_sprecaply(prec,x,y,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! - ! - ! Number of levels = 1: apply the base preconditioner - ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. ! Will need some rethinking. ! - allocate(ww(size(x))) - + call psb_geasb(w1,desc_data,info,scratch=.true.) + call psb_geasb(w2,desc_data,info,scratch=.true.) + + call psb_geaxpby(sone,x,szero,w1,desc_data,info) select case(trans_) case ('N') - call prec%precv(1)%sm%apply(sone,x,szero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(sone,ww,szero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(sone,w1,szero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(sone,w2,szero,w1,desc_data,trans_,& + & ione, work_,info) + end do + case('T','C') - call prec%precv(1)%sm2a%apply(sone,x,szero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(sone,ww,szero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(sone,w1,szero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(sone,w2,szero,w1,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') goto 9999 end select - deallocate(ww) - + call psb_geaxpby(sone,w1,szero,y,desc_data,info) + call psb_gefree(w1,desc_data,info) + call psb_gefree(w2,desc_data,info) + else call prec%precv(1)%sm%apply(sone,x,szero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post), & - & work_,info) + & nswps, work_,info) end if else info = psb_err_from_subroutine_ai_ @@ -255,7 +257,7 @@ subroutine mld_sprecaply1(prec,x,desc_data,info,trans) ! Local variables integer(psb_ipk_) :: ictxt,np,me integer(psb_ipk_) :: err_act - real(psb_spk_), pointer :: WW(:), w1(:) + real(psb_spk_), pointer :: ww(:), w1(:) character(len=20) :: name name='mld_sprecaply1' @@ -282,7 +284,7 @@ subroutine mld_sprecaply1(prec,x,desc_data,info,trans) end if x(:) = ww(:) - deallocate(ww,W1,stat=info) + deallocate(ww,w1,stat=info) if (info /= psb_success_) then info = psb_err_alloc_dealloc_ call psb_errpush(info,name) @@ -319,7 +321,7 @@ subroutine mld_sprecaply2_vect(prec,x,y,desc_data,info,trans,work) character :: trans_ real(psb_spk_), pointer :: work_(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_sprecaply' @@ -370,42 +372,47 @@ subroutine mld_sprecaply2_vect(prec,x,y,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. ! Will need some rethinking. ! twoside: block - type(psb_s_vect_type) :: ww - call psb_geasb(ww,desc_data,info,mold=x%v,scratch=.true.) - + type(psb_s_vect_type) :: w1,w2 + call psb_geasb(w1,desc_data,info,mold=x%v,scratch=.true.) + call psb_geasb(w2,desc_data,info,mold=x%v,scratch=.true.) + call psb_geaxpby(sone,x,szero,w1,desc_data,info) select case(trans_) case ('N') - call prec%precv(1)%sm%apply(sone,x,szero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(sone,ww,szero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(sone,w1,szero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(sone,w2,szero,w1,desc_data,trans_,& + & ione, work_,info) + end do + case('T','C') - call prec%precv(1)%sm2a%apply(sone,x,szero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(sone,ww,szero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(sone,w1,szero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(sone,w2,szero,w1,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') goto 9999 end select - call psb_gefree(ww,desc_data,info) + call psb_geaxpby(sone,w1,szero,y,desc_data,info) + call psb_gefree(w1,desc_data,info) + call psb_gefree(w2,desc_data,info) end block twoside else call prec%precv(1)%sm%apply(sone,x,szero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + & nswps,work_,info) end if else @@ -454,7 +461,7 @@ subroutine mld_sprecaply1_vect(prec,x,desc_data,info,trans,work) type(psb_s_vect_type) :: ww real(psb_spk_), pointer :: work_(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_sprecaply' @@ -506,6 +513,7 @@ subroutine mld_sprecaply1_vect(prec,x,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. @@ -513,19 +521,19 @@ subroutine mld_sprecaply1_vect(prec,x,desc_data,info,trans,work) ! select case(trans_) case ('N') - call prec%precv(1)%sm%apply(sone,x,szero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(sone,ww,szero,x,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(sone,x,szero,ww,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(sone,ww,szero,x,desc_data,trans_,& + & ione, work_,info) + end do case('T','C') - call prec%precv(1)%sm2a%apply(sone,x,szero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(sone,ww,szero,x,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(sone,x,szero,ww,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(sone,ww,szero,x,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') @@ -534,8 +542,7 @@ subroutine mld_sprecaply1_vect(prec,x,desc_data,info,trans,work) else call prec%precv(1)%sm%apply(sone,x,szero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + & nswps, work_,info) if (info == 0) call psb_geaxpby(sone,ww,szero,x,desc_data,info) end if else diff --git a/mlprec/impl/mld_zfile_prec_descr.f90 b/mlprec/impl/mld_zfile_prec_descr.f90 new file mode 100644 index 00000000..a468761e --- /dev/null +++ b/mlprec/impl/mld_zfile_prec_descr.f90 @@ -0,0 +1,193 @@ +! +! +! MLD2P4 version 2.1 +! MultiLevel Domain Decomposition Parallel Preconditioners Package +! based on PSBLAS (Parallel Sparse BLAS version 3.5) +! +! (C) Copyright 2008, 2010, 2012, 2015, 2017 +! +! Salvatore Filippone Cranfield University, UK +! Pasqua D'Ambra IAC-CNR, Naples, IT +! Daniela di Serafino University of Campania "L. Vanvitelli", Caserta, IT +! +! Redistribution and use in source and binary forms, with or without +! modification, are permitted provided that the following conditions +! are met: +! 1. Redistributions of source code must retain the above copyright +! notice, this list of conditions and the following disclaimer. +! 2. Redistributions in binary form must reproduce the above copyright +! notice, this list of conditions, and the following disclaimer in the +! documentation and/or other materials provided with the distribution. +! 3. The name of the MLD2P4 group or the names of its contributors may +! not be used to endorse or promote products derived from this +! software without specific written permission. +! +! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED +! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE MLD2P4 GROUP OR ITS CONTRIBUTORS +! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +! POSSIBILITY OF SUCH DAMAGE. +! +! +! File: mld_dfile_prec_descr.f90 +! +! +! Subroutine: mld_file_prec_descr +! Version: complex +! +! This routine prints a description of the preconditioner to the standard +! output or to a file. It must be called after the preconditioner has been +! built by mld_precbld. +! +! Arguments: +! p - type(mld_Tprec_type), input. +! The preconditioner data structure to be printed out. +! info - integer, output. +! error code. +! iout - integer, input, optional. +! The id of the file where the preconditioner description +! will be printed. If iout is not present, then the standard +! output is condidered. +! root - integer, input, optional. +! The id of the process printing the message; -1 acts as a wildcard. +! Default is psb_root_ +! +subroutine mld_zfile_prec_descr(prec,iout,root) + use psb_base_mod + use mld_z_prec_mod, mld_protect_name => mld_zfile_prec_descr + use mld_z_inner_mod + use mld_z_gs_solver + + implicit none + ! Arguments + class(mld_zprec_type), intent(in) :: prec + integer(psb_ipk_), intent(in), optional :: iout + integer(psb_ipk_), intent(in), optional :: root + + ! Local variables + integer(psb_ipk_) :: ilev, nlev, ilmin, info, nswps + integer(psb_ipk_) :: ictxt, me, np + logical :: is_symgs + character(len=20), parameter :: name='mld_file_prec_descr' + integer(psb_ipk_) :: iout_ + integer(psb_ipk_) :: root_ + + info = psb_success_ + if (present(iout)) then + iout_ = iout + else + iout_ = psb_out_unit + end if + if (iout_ < 0) iout_ = psb_out_unit + + ictxt = prec%ictxt + + if (allocated(prec%precv)) then + + call psb_info(ictxt,me,np) + if (present(root)) then + root_ = root + else + root_ = psb_root_ + end if + if (root_ == -1) root_ = me + + ! + ! The preconditioner description is printed by processor psb_root_. + ! This agrees with the fact that all the parameters defining the + ! preconditioner have the same values on all the procs (this is + ! ensured by mld_precbld). + ! + if (me == root_) then + nlev = size(prec%precv) + do ilev = 1, nlev + if (.not.allocated(prec%precv(ilev)%sm)) then + info = 3111 + write(iout_,*) ' ',name,& + & ': error: inconsistent MLPREC part, should call MLD_PRECINIT' + return + endif + end do + + write(iout_,*) + write(iout_,'(a)') 'Preconditioner description' + + if (nlev == 1) then + ! + ! Here we have a gigantic kludge just to handle Symmetrized Gauss-Seidel. + ! Will need rethinking... + ! + if (allocated(prec%precv(1)%sm2a)) then + is_symgs = .false. + select type(sv2 => prec%precv(1)%sm2a%sv) + class is (mld_z_bwgs_solver_type) + select type(sv1 => prec%precv(1)%sm%sv) + class is (mld_z_gs_solver_type) + is_symgs = .true. + end select + end select + if (is_symgs) then + write(iout_,*) ' Forward-Backward (symmetrized) Hybrid Gauss-Seidel' + else + call prec%precv(1)%sm%descr(info,iout=iout_) + call prec%precv(1)%sm2a%descr(info,iout=iout_) + end if + + else + call prec%precv(1)%sm%descr(info,iout=iout_) + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + end if + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + if (nswps > 1) write(iout_,*) ' Number of sweeps : ',nswps + write(iout_,*) + + else if (nlev > 1) then + ! + ! Print description of base preconditioner + ! + write(iout_,*) 'Multilevel Preconditioner' + write(iout_,*) 'Outer sweeps:',prec%outer_sweeps + write(iout_,*) + write(iout_,*) 'Smoother: ' + call prec%precv(1)%sm%descr(info,iout=iout_) + if (allocated(prec%precv(1)%sm2a)) then + write(iout_,*) 'Post smoother:' + call prec%precv(1)%sm2a%descr(info,iout=iout_) + end if + ! + ! Print multilevel details + ! + write(iout_,*) + write(iout_,*) 'Multilevel hierarchy: ' + write(iout_,*) ' Number of levels : ',nlev + write(iout_,*) ' Operator complexity: ',prec%get_complexity() + ilmin = 2 + if (nlev == 2) ilmin=1 + do ilev=ilmin,nlev + call prec%precv(ilev)%descr(ilev,nlev,ilmin,info,iout=iout_) + end do + write(iout_,*) + + else + write(iout_,*) trim(name), & + & ': invalid preconditioner array size ?',nlev + info = -2 + return + + end if + end if + + else + write(iout_,*) trim(name), & + & ': Error: no base preconditioner available, something is wrong!' + info = -2 + return + endif + +end subroutine mld_zfile_prec_descr diff --git a/mlprec/impl/mld_zprecaply.f90 b/mlprec/impl/mld_zprecaply.f90 index ee417635..a87d6239 100644 --- a/mlprec/impl/mld_zprecaply.f90 +++ b/mlprec/impl/mld_zprecaply.f90 @@ -89,10 +89,10 @@ subroutine mld_zprecaply(prec,x,y,desc_data,info,trans,work) ! Local variables character :: trans_ complex(psb_dpk_), pointer :: work_(:) - complex(psb_dpk_), allocatable :: ww(:) + complex(psb_dpk_), allocatable :: w1(:), w2(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_zprecaply' @@ -143,42 +143,44 @@ subroutine mld_zprecaply(prec,x,y,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! - ! - ! Number of levels = 1: apply the base preconditioner - ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. ! Will need some rethinking. ! - allocate(ww(size(x))) - + call psb_geasb(w1,desc_data,info,scratch=.true.) + call psb_geasb(w2,desc_data,info,scratch=.true.) + + call psb_geaxpby(zone,x,zzero,w1,desc_data,info) select case(trans_) case ('N') - call prec%precv(1)%sm%apply(zone,x,zzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(zone,ww,zzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(zone,w1,zzero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(zone,w2,zzero,w1,desc_data,trans_,& + & ione, work_,info) + end do + case('T','C') - call prec%precv(1)%sm2a%apply(zone,x,zzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(zone,ww,zzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(zone,w1,zzero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(zone,w2,zzero,w1,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') goto 9999 end select - deallocate(ww) - + call psb_geaxpby(zone,w1,zzero,y,desc_data,info) + call psb_gefree(w1,desc_data,info) + call psb_gefree(w2,desc_data,info) + else call prec%precv(1)%sm%apply(zone,x,zzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post), & - & work_,info) + & nswps, work_,info) end if else info = psb_err_from_subroutine_ai_ @@ -255,7 +257,7 @@ subroutine mld_zprecaply1(prec,x,desc_data,info,trans) ! Local variables integer(psb_ipk_) :: ictxt,np,me integer(psb_ipk_) :: err_act - complex(psb_dpk_), pointer :: WW(:), w1(:) + complex(psb_dpk_), pointer :: ww(:), w1(:) character(len=20) :: name name='mld_zprecaply1' @@ -282,7 +284,7 @@ subroutine mld_zprecaply1(prec,x,desc_data,info,trans) end if x(:) = ww(:) - deallocate(ww,W1,stat=info) + deallocate(ww,w1,stat=info) if (info /= psb_success_) then info = psb_err_alloc_dealloc_ call psb_errpush(info,name) @@ -319,7 +321,7 @@ subroutine mld_zprecaply2_vect(prec,x,y,desc_data,info,trans,work) character :: trans_ complex(psb_dpk_), pointer :: work_(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_zprecaply' @@ -370,42 +372,47 @@ subroutine mld_zprecaply2_vect(prec,x,y,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) + if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. ! Will need some rethinking. ! twoside: block - type(psb_z_vect_type) :: ww - call psb_geasb(ww,desc_data,info,mold=x%v,scratch=.true.) - + type(psb_z_vect_type) :: w1,w2 + call psb_geasb(w1,desc_data,info,mold=x%v,scratch=.true.) + call psb_geasb(w2,desc_data,info,mold=x%v,scratch=.true.) + call psb_geaxpby(zone,x,zzero,w1,desc_data,info) select case(trans_) case ('N') - call prec%precv(1)%sm%apply(zone,x,zzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(zone,ww,zzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(zone,w1,zzero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(zone,w2,zzero,w1,desc_data,trans_,& + & ione, work_,info) + end do + case('T','C') - call prec%precv(1)%sm2a%apply(zone,x,zzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(zone,ww,zzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(zone,w1,zzero,w2,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(zone,w2,zzero,w1,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') goto 9999 end select - call psb_gefree(ww,desc_data,info) + call psb_geaxpby(zone,w1,zzero,y,desc_data,info) + call psb_gefree(w1,desc_data,info) + call psb_gefree(w2,desc_data,info) end block twoside else call prec%precv(1)%sm%apply(zone,x,zzero,y,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + & nswps,work_,info) end if else @@ -454,7 +461,7 @@ subroutine mld_zprecaply1_vect(prec,x,desc_data,info,trans,work) type(psb_z_vect_type) :: ww complex(psb_dpk_), pointer :: work_(:) integer(psb_ipk_) :: ictxt,np,me - integer(psb_ipk_) :: err_act,iwsz + integer(psb_ipk_) :: err_act,iwsz, k, nswps character(len=20) :: name name='mld_zprecaply' @@ -506,6 +513,7 @@ subroutine mld_zprecaply1_vect(prec,x,desc_data,info,trans,work) ! ! Number of levels = 1: apply the base preconditioner ! + nswps = max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post) if (allocated(prec%precv(1)%sm2a)) then ! ! This is a kludge for handling the symmetrized GS case. @@ -513,19 +521,19 @@ subroutine mld_zprecaply1_vect(prec,x,desc_data,info,trans,work) ! select case(trans_) case ('N') - call prec%precv(1)%sm%apply(zone,x,zzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm2a%apply(zone,ww,zzero,x,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm%apply(zone,x,zzero,ww,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm2a%apply(zone,ww,zzero,x,desc_data,trans_,& + & ione, work_,info) + end do case('T','C') - call prec%precv(1)%sm2a%apply(zone,x,zzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) - call prec%precv(1)%sm%apply(zone,ww,zzero,x,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + do k=1, nswps + call prec%precv(1)%sm2a%apply(zone,x,zzero,ww,desc_data,trans_,& + & ione, work_,info) + call prec%precv(1)%sm%apply(zone,ww,zzero,x,desc_data,trans_,& + & ione, work_,info) + end do case default info = psb_err_from_subroutine_ call psb_errpush(info,name,a_err='Invalid trans') @@ -534,8 +542,7 @@ subroutine mld_zprecaply1_vect(prec,x,desc_data,info,trans,work) else call prec%precv(1)%sm%apply(zone,x,zzero,ww,desc_data,trans_,& - & max(prec%precv(1)%parms%sweeps_pre,prec%precv(1)%parms%sweeps_post),& - & work_,info) + & nswps, work_,info) if (info == 0) call psb_geaxpby(zone,ww,zzero,x,desc_data,info) end if else diff --git a/mlprec/impl/smoother/mld_c_base_smoother_descr.f90 b/mlprec/impl/smoother/mld_c_base_smoother_descr.f90 index ffe6eff2..f59ab7f4 100644 --- a/mlprec/impl/smoother/mld_c_base_smoother_descr.f90 +++ b/mlprec/impl/smoother/mld_c_base_smoother_descr.f90 @@ -70,7 +70,7 @@ subroutine mld_c_base_smoother_descr(sm,info,iout,coarse) end if if (.not.coarse_) & - & write(iout_,*) 'Base smoother with local solver' + & write(iout_,*) 'Decoupled preconditioner/smoother with local solver' if (allocated(sm%sv)) then call sm%sv%descr(info,iout,coarse) if (info /= psb_success_) then diff --git a/mlprec/impl/smoother/mld_c_jac_smoother_descr.f90 b/mlprec/impl/smoother/mld_c_jac_smoother_descr.f90 index c54d11b0..9c91621a 100644 --- a/mlprec/impl/smoother/mld_c_jac_smoother_descr.f90 +++ b/mlprec/impl/smoother/mld_c_jac_smoother_descr.f90 @@ -74,23 +74,19 @@ subroutine mld_c_jac_smoother_descr(sm,info,iout,coarse) if (allocated(sm%sv)) then select type(smv=>sm%sv) class is (mld_c_diag_solver_type) - write(iout_,*) ' Point Jacobi smoother ' + write(iout_,*) ' Point Jacobi ' class is (mld_c_bwgs_solver_type) - write(iout_,*) ' Hybrid Backward Gauss-Seidel smoother ' -!!$ write(iout_,*) ' Local Gauss-Seidel solver details:' -!!$ call smv%descr(info,iout_,coarse=coarse) + write(iout_,*) ' Hybrid Backward Gauss-Seidel ' class is (mld_c_gs_solver_type) - write(iout_,*) ' Hybrid Forward Gauss-Seidel smoother ' -!!$ write(iout_,*) ' Local Gauss-Seidel solver details:' -!!$ call smv%descr(info,iout_,coarse=coarse) + write(iout_,*) ' Hybrid Forward Gauss-Seidel ' class default - write(iout_,*) ' Block Jacobi smoother ' + write(iout_,*) ' Block Jacobi ' write(iout_,*) ' Local solver details:' call smv%descr(info,iout_,coarse=coarse) end select else - write(iout_,*) ' Block Jacobi smoother ' + write(iout_,*) ' Block Jacobi ' end if else if (allocated(sm%sv)) then diff --git a/mlprec/impl/smoother/mld_d_base_smoother_descr.f90 b/mlprec/impl/smoother/mld_d_base_smoother_descr.f90 index ea9b282c..a55c410d 100644 --- a/mlprec/impl/smoother/mld_d_base_smoother_descr.f90 +++ b/mlprec/impl/smoother/mld_d_base_smoother_descr.f90 @@ -70,7 +70,7 @@ subroutine mld_d_base_smoother_descr(sm,info,iout,coarse) end if if (.not.coarse_) & - & write(iout_,*) 'Base smoother with local solver' + & write(iout_,*) 'Decoupled preconditioner/smoother with local solver' if (allocated(sm%sv)) then call sm%sv%descr(info,iout,coarse) if (info /= psb_success_) then diff --git a/mlprec/impl/smoother/mld_d_jac_smoother_descr.f90 b/mlprec/impl/smoother/mld_d_jac_smoother_descr.f90 index d39e5cee..757c67b5 100644 --- a/mlprec/impl/smoother/mld_d_jac_smoother_descr.f90 +++ b/mlprec/impl/smoother/mld_d_jac_smoother_descr.f90 @@ -74,23 +74,19 @@ subroutine mld_d_jac_smoother_descr(sm,info,iout,coarse) if (allocated(sm%sv)) then select type(smv=>sm%sv) class is (mld_d_diag_solver_type) - write(iout_,*) ' Point Jacobi smoother ' + write(iout_,*) ' Point Jacobi ' class is (mld_d_bwgs_solver_type) - write(iout_,*) ' Hybrid Backward Gauss-Seidel smoother ' -!!$ write(iout_,*) ' Local Gauss-Seidel solver details:' -!!$ call smv%descr(info,iout_,coarse=coarse) + write(iout_,*) ' Hybrid Backward Gauss-Seidel ' class is (mld_d_gs_solver_type) - write(iout_,*) ' Hybrid Forward Gauss-Seidel smoother ' -!!$ write(iout_,*) ' Local Gauss-Seidel solver details:' -!!$ call smv%descr(info,iout_,coarse=coarse) + write(iout_,*) ' Hybrid Forward Gauss-Seidel ' class default - write(iout_,*) ' Block Jacobi smoother ' + write(iout_,*) ' Block Jacobi ' write(iout_,*) ' Local solver details:' call smv%descr(info,iout_,coarse=coarse) end select else - write(iout_,*) ' Block Jacobi smoother ' + write(iout_,*) ' Block Jacobi ' end if else if (allocated(sm%sv)) then diff --git a/mlprec/impl/smoother/mld_s_base_smoother_descr.f90 b/mlprec/impl/smoother/mld_s_base_smoother_descr.f90 index 4d12c286..7da08114 100644 --- a/mlprec/impl/smoother/mld_s_base_smoother_descr.f90 +++ b/mlprec/impl/smoother/mld_s_base_smoother_descr.f90 @@ -70,7 +70,7 @@ subroutine mld_s_base_smoother_descr(sm,info,iout,coarse) end if if (.not.coarse_) & - & write(iout_,*) 'Base smoother with local solver' + & write(iout_,*) 'Decoupled preconditioner/smoother with local solver' if (allocated(sm%sv)) then call sm%sv%descr(info,iout,coarse) if (info /= psb_success_) then diff --git a/mlprec/impl/smoother/mld_s_jac_smoother_descr.f90 b/mlprec/impl/smoother/mld_s_jac_smoother_descr.f90 index 2a02c6d8..d1b191df 100644 --- a/mlprec/impl/smoother/mld_s_jac_smoother_descr.f90 +++ b/mlprec/impl/smoother/mld_s_jac_smoother_descr.f90 @@ -74,23 +74,19 @@ subroutine mld_s_jac_smoother_descr(sm,info,iout,coarse) if (allocated(sm%sv)) then select type(smv=>sm%sv) class is (mld_s_diag_solver_type) - write(iout_,*) ' Point Jacobi smoother ' + write(iout_,*) ' Point Jacobi ' class is (mld_s_bwgs_solver_type) - write(iout_,*) ' Hybrid Backward Gauss-Seidel smoother ' -!!$ write(iout_,*) ' Local Gauss-Seidel solver details:' -!!$ call smv%descr(info,iout_,coarse=coarse) + write(iout_,*) ' Hybrid Backward Gauss-Seidel ' class is (mld_s_gs_solver_type) - write(iout_,*) ' Hybrid Forward Gauss-Seidel smoother ' -!!$ write(iout_,*) ' Local Gauss-Seidel solver details:' -!!$ call smv%descr(info,iout_,coarse=coarse) + write(iout_,*) ' Hybrid Forward Gauss-Seidel ' class default - write(iout_,*) ' Block Jacobi smoother ' + write(iout_,*) ' Block Jacobi ' write(iout_,*) ' Local solver details:' call smv%descr(info,iout_,coarse=coarse) end select else - write(iout_,*) ' Block Jacobi smoother ' + write(iout_,*) ' Block Jacobi ' end if else if (allocated(sm%sv)) then diff --git a/mlprec/impl/smoother/mld_z_base_smoother_descr.f90 b/mlprec/impl/smoother/mld_z_base_smoother_descr.f90 index 1359e7f1..64777267 100644 --- a/mlprec/impl/smoother/mld_z_base_smoother_descr.f90 +++ b/mlprec/impl/smoother/mld_z_base_smoother_descr.f90 @@ -70,7 +70,7 @@ subroutine mld_z_base_smoother_descr(sm,info,iout,coarse) end if if (.not.coarse_) & - & write(iout_,*) 'Base smoother with local solver' + & write(iout_,*) 'Decoupled preconditioner/smoother with local solver' if (allocated(sm%sv)) then call sm%sv%descr(info,iout,coarse) if (info /= psb_success_) then diff --git a/mlprec/impl/smoother/mld_z_jac_smoother_descr.f90 b/mlprec/impl/smoother/mld_z_jac_smoother_descr.f90 index 91011454..644397c6 100644 --- a/mlprec/impl/smoother/mld_z_jac_smoother_descr.f90 +++ b/mlprec/impl/smoother/mld_z_jac_smoother_descr.f90 @@ -74,23 +74,19 @@ subroutine mld_z_jac_smoother_descr(sm,info,iout,coarse) if (allocated(sm%sv)) then select type(smv=>sm%sv) class is (mld_z_diag_solver_type) - write(iout_,*) ' Point Jacobi smoother ' + write(iout_,*) ' Point Jacobi ' class is (mld_z_bwgs_solver_type) - write(iout_,*) ' Hybrid Backward Gauss-Seidel smoother ' -!!$ write(iout_,*) ' Local Gauss-Seidel solver details:' -!!$ call smv%descr(info,iout_,coarse=coarse) + write(iout_,*) ' Hybrid Backward Gauss-Seidel ' class is (mld_z_gs_solver_type) - write(iout_,*) ' Hybrid Forward Gauss-Seidel smoother ' -!!$ write(iout_,*) ' Local Gauss-Seidel solver details:' -!!$ call smv%descr(info,iout_,coarse=coarse) + write(iout_,*) ' Hybrid Forward Gauss-Seidel ' class default - write(iout_,*) ' Block Jacobi smoother ' + write(iout_,*) ' Block Jacobi ' write(iout_,*) ' Local solver details:' call smv%descr(info,iout_,coarse=coarse) end select else - write(iout_,*) ' Block Jacobi smoother ' + write(iout_,*) ' Block Jacobi ' end if else if (allocated(sm%sv)) then diff --git a/mlprec/mld_c_as_smoother.f90 b/mlprec/mld_c_as_smoother.f90 index 0562eb73..2c067aad 100644 --- a/mlprec/mld_c_as_smoother.f90 +++ b/mlprec/mld_c_as_smoother.f90 @@ -461,7 +461,7 @@ contains implicit none character(len=32) :: val - val = "Schwarz smoother" + val = "Additive Schwarz" end function c_as_smoother_get_fmt function c_as_smoother_get_id() result(val) diff --git a/mlprec/mld_c_prec_type.f90 b/mlprec/mld_c_prec_type.f90 index cb0dbb8b..f83c096a 100644 --- a/mlprec/mld_c_prec_type.f90 +++ b/mlprec/mld_c_prec_type.f90 @@ -156,7 +156,14 @@ module mld_c_prec_type interface mld_precdescr - module procedure mld_cfile_prec_descr + subroutine mld_cfile_prec_descr(prec,iout,root) + import :: mld_cprec_type, psb_ipk_ + implicit none + ! Arguments + class(mld_cprec_type), intent(in) :: prec + integer(psb_ipk_), intent(in), optional :: iout + integer(psb_ipk_), intent(in), optional :: root + end subroutine mld_cfile_prec_descr end interface interface mld_sizeof @@ -489,139 +496,6 @@ contains prec%op_complexity = num/den end subroutine mld_c_cmp_compl - ! - ! Subroutine: mld_file_prec_descr - ! Version: complex - ! - ! This routine prints a description of the preconditioner to the standard - ! output or to a file. It must be called after the preconditioner has been - ! built by mld_precbld. - ! - ! Arguments: - ! p - type(mld_Tprec_type), input. - ! The preconditioner data structure to be printed out. - ! info - integer, output. - ! error code. - ! iout - integer, input, optional. - ! The id of the file where the preconditioner description - ! will be printed. If iout is not present, then the standard - ! output is condidered. - ! root - integer, input, optional. - ! The id of the process printing the message; -1 acts as a wildcard. - ! Default is psb_root_ - ! - subroutine mld_cfile_prec_descr(prec,iout,root) - implicit none - ! Arguments - class(mld_cprec_type), intent(in) :: prec - integer(psb_ipk_), intent(in), optional :: iout - integer(psb_ipk_), intent(in), optional :: root - - ! Local variables - integer(psb_ipk_) :: ilev, nlev, ilmin, info - integer(psb_ipk_) :: ictxt, me, np - character(len=20), parameter :: name='mld_file_prec_descr' - integer(psb_ipk_) :: iout_ - integer(psb_ipk_) :: root_ - - info = psb_success_ - if (present(iout)) then - iout_ = iout - else - iout_ = psb_out_unit - end if - if (iout_ < 0) iout_ = psb_out_unit - - ictxt = prec%ictxt - - if (allocated(prec%precv)) then - - call psb_info(ictxt,me,np) - if (present(root)) then - root_ = root - else - root_ = psb_root_ - end if - if (root_ == -1) root_ = me - - ! - ! The preconditioner description is printed by processor psb_root_. - ! This agrees with the fact that all the parameters defining the - ! preconditioner have the same values on all the procs (this is - ! ensured by mld_precbld). - ! - if (me == root_) then - nlev = size(prec%precv) - do ilev = 1, nlev - if (.not.allocated(prec%precv(ilev)%sm)) then - info = 3111 - write(iout_,*) ' ',name,& - & ': error: inconsistent MLPREC part, should call MLD_PRECINIT' - return - endif - end do - - write(iout_,*) - write(iout_,'(a)') 'Preconditioner description' - if (nlev >= 1) then - ! - ! Print description of base preconditioner - ! - if (nlev > 1) then - write(iout_,*) 'Multilevel Preconditioner' - write(iout_,*) 'Outer sweeps:',prec%outer_sweeps - write(iout_,*) - write(iout_,*) 'Base preconditioner (smoother) details' - endif - call prec%precv(1)%sm%descr(info,iout=iout_) - if (nlev == 1) then - if (prec%precv(1)%parms%sweeps_pre > 1) then - write(iout_,*) ' Number of smoother sweeps_pre : ',& - & prec%precv(1)%parms%sweeps_pre - end if - write(iout_,*) - end if - if (allocated(prec%precv(1)%sm2a)) then - write(iout_,*) 'Post smoother details' - call prec%precv(1)%sm2a%descr(info,iout=iout_) - if (nlev == 1) then - if (prec%precv(1)%parms%sweeps_post > 1) then - write(iout_,*) ' Number of smoother sweeps_post : ',& - & prec%precv(1)%parms%sweeps_post - end if - write(iout_,*) - end if - end if - end if - - if (nlev > 1) then - ! - ! Print multilevel details - ! - write(iout_,*) - write(iout_,*) 'Multilevel details' - write(iout_,*) ' Number of levels : ',nlev - write(iout_,*) ' Operator complexity: ',prec%get_complexity() - ilmin = 2 - if (nlev == 2) ilmin=1 - do ilev=ilmin,nlev - call prec%precv(ilev)%descr(ilev,nlev,ilmin,info,iout=iout_) - end do - write(iout_,*) - end if - end if - - else - write(iout_,*) trim(name), & - & ': Error: no base preconditioner available, something is wrong!' - info = -2 - return - endif - - - end subroutine mld_cfile_prec_descr - - ! ! Subroutines: mld_Tprec_free ! Version: complex diff --git a/mlprec/mld_d_as_smoother.f90 b/mlprec/mld_d_as_smoother.f90 index 0ed795d3..096fdaa9 100644 --- a/mlprec/mld_d_as_smoother.f90 +++ b/mlprec/mld_d_as_smoother.f90 @@ -461,7 +461,7 @@ contains implicit none character(len=32) :: val - val = "Schwarz smoother" + val = "Additive Schwarz" end function d_as_smoother_get_fmt function d_as_smoother_get_id() result(val) diff --git a/mlprec/mld_d_prec_type.f90 b/mlprec/mld_d_prec_type.f90 index bc489f7c..b3419c30 100644 --- a/mlprec/mld_d_prec_type.f90 +++ b/mlprec/mld_d_prec_type.f90 @@ -156,7 +156,14 @@ module mld_d_prec_type interface mld_precdescr - module procedure mld_dfile_prec_descr + subroutine mld_dfile_prec_descr(prec,iout,root) + import :: mld_dprec_type, psb_ipk_ + implicit none + ! Arguments + class(mld_dprec_type), intent(in) :: prec + integer(psb_ipk_), intent(in), optional :: iout + integer(psb_ipk_), intent(in), optional :: root + end subroutine mld_dfile_prec_descr end interface interface mld_sizeof @@ -489,139 +496,6 @@ contains prec%op_complexity = num/den end subroutine mld_d_cmp_compl - ! - ! Subroutine: mld_file_prec_descr - ! Version: real - ! - ! This routine prints a description of the preconditioner to the standard - ! output or to a file. It must be called after the preconditioner has been - ! built by mld_precbld. - ! - ! Arguments: - ! p - type(mld_Tprec_type), input. - ! The preconditioner data structure to be printed out. - ! info - integer, output. - ! error code. - ! iout - integer, input, optional. - ! The id of the file where the preconditioner description - ! will be printed. If iout is not present, then the standard - ! output is condidered. - ! root - integer, input, optional. - ! The id of the process printing the message; -1 acts as a wildcard. - ! Default is psb_root_ - ! - subroutine mld_dfile_prec_descr(prec,iout,root) - implicit none - ! Arguments - class(mld_dprec_type), intent(in) :: prec - integer(psb_ipk_), intent(in), optional :: iout - integer(psb_ipk_), intent(in), optional :: root - - ! Local variables - integer(psb_ipk_) :: ilev, nlev, ilmin, info - integer(psb_ipk_) :: ictxt, me, np - character(len=20), parameter :: name='mld_file_prec_descr' - integer(psb_ipk_) :: iout_ - integer(psb_ipk_) :: root_ - - info = psb_success_ - if (present(iout)) then - iout_ = iout - else - iout_ = psb_out_unit - end if - if (iout_ < 0) iout_ = psb_out_unit - - ictxt = prec%ictxt - - if (allocated(prec%precv)) then - - call psb_info(ictxt,me,np) - if (present(root)) then - root_ = root - else - root_ = psb_root_ - end if - if (root_ == -1) root_ = me - - ! - ! The preconditioner description is printed by processor psb_root_. - ! This agrees with the fact that all the parameters defining the - ! preconditioner have the same values on all the procs (this is - ! ensured by mld_precbld). - ! - if (me == root_) then - nlev = size(prec%precv) - do ilev = 1, nlev - if (.not.allocated(prec%precv(ilev)%sm)) then - info = 3111 - write(iout_,*) ' ',name,& - & ': error: inconsistent MLPREC part, should call MLD_PRECINIT' - return - endif - end do - - write(iout_,*) - write(iout_,'(a)') 'Preconditioner description' - if (nlev >= 1) then - ! - ! Print description of base preconditioner - ! - if (nlev > 1) then - write(iout_,*) 'Multilevel Preconditioner' - write(iout_,*) 'Outer sweeps:',prec%outer_sweeps - write(iout_,*) - write(iout_,*) 'Base preconditioner (smoother) details' - endif - call prec%precv(1)%sm%descr(info,iout=iout_) - if (nlev == 1) then - if (prec%precv(1)%parms%sweeps_pre > 1) then - write(iout_,*) ' Number of smoother sweeps_pre : ',& - & prec%precv(1)%parms%sweeps_pre - end if - write(iout_,*) - end if - if (allocated(prec%precv(1)%sm2a)) then - write(iout_,*) 'Post smoother details' - call prec%precv(1)%sm2a%descr(info,iout=iout_) - if (nlev == 1) then - if (prec%precv(1)%parms%sweeps_post > 1) then - write(iout_,*) ' Number of smoother sweeps_post : ',& - & prec%precv(1)%parms%sweeps_post - end if - write(iout_,*) - end if - end if - end if - - if (nlev > 1) then - ! - ! Print multilevel details - ! - write(iout_,*) - write(iout_,*) 'Multilevel details' - write(iout_,*) ' Number of levels : ',nlev - write(iout_,*) ' Operator complexity: ',prec%get_complexity() - ilmin = 2 - if (nlev == 2) ilmin=1 - do ilev=ilmin,nlev - call prec%precv(ilev)%descr(ilev,nlev,ilmin,info,iout=iout_) - end do - write(iout_,*) - end if - end if - - else - write(iout_,*) trim(name), & - & ': Error: no base preconditioner available, something is wrong!' - info = -2 - return - endif - - - end subroutine mld_dfile_prec_descr - - ! ! Subroutines: mld_Tprec_free ! Version: real diff --git a/mlprec/mld_s_as_smoother.f90 b/mlprec/mld_s_as_smoother.f90 index c1e8669c..79cd9c94 100644 --- a/mlprec/mld_s_as_smoother.f90 +++ b/mlprec/mld_s_as_smoother.f90 @@ -461,7 +461,7 @@ contains implicit none character(len=32) :: val - val = "Schwarz smoother" + val = "Additive Schwarz" end function s_as_smoother_get_fmt function s_as_smoother_get_id() result(val) diff --git a/mlprec/mld_s_prec_type.f90 b/mlprec/mld_s_prec_type.f90 index cc1160a5..b3d74b94 100644 --- a/mlprec/mld_s_prec_type.f90 +++ b/mlprec/mld_s_prec_type.f90 @@ -156,7 +156,14 @@ module mld_s_prec_type interface mld_precdescr - module procedure mld_sfile_prec_descr + subroutine mld_sfile_prec_descr(prec,iout,root) + import :: mld_sprec_type, psb_ipk_ + implicit none + ! Arguments + class(mld_sprec_type), intent(in) :: prec + integer(psb_ipk_), intent(in), optional :: iout + integer(psb_ipk_), intent(in), optional :: root + end subroutine mld_sfile_prec_descr end interface interface mld_sizeof @@ -489,139 +496,6 @@ contains prec%op_complexity = num/den end subroutine mld_s_cmp_compl - ! - ! Subroutine: mld_file_prec_descr - ! Version: real - ! - ! This routine prints a description of the preconditioner to the standard - ! output or to a file. It must be called after the preconditioner has been - ! built by mld_precbld. - ! - ! Arguments: - ! p - type(mld_Tprec_type), input. - ! The preconditioner data structure to be printed out. - ! info - integer, output. - ! error code. - ! iout - integer, input, optional. - ! The id of the file where the preconditioner description - ! will be printed. If iout is not present, then the standard - ! output is condidered. - ! root - integer, input, optional. - ! The id of the process printing the message; -1 acts as a wildcard. - ! Default is psb_root_ - ! - subroutine mld_sfile_prec_descr(prec,iout,root) - implicit none - ! Arguments - class(mld_sprec_type), intent(in) :: prec - integer(psb_ipk_), intent(in), optional :: iout - integer(psb_ipk_), intent(in), optional :: root - - ! Local variables - integer(psb_ipk_) :: ilev, nlev, ilmin, info - integer(psb_ipk_) :: ictxt, me, np - character(len=20), parameter :: name='mld_file_prec_descr' - integer(psb_ipk_) :: iout_ - integer(psb_ipk_) :: root_ - - info = psb_success_ - if (present(iout)) then - iout_ = iout - else - iout_ = psb_out_unit - end if - if (iout_ < 0) iout_ = psb_out_unit - - ictxt = prec%ictxt - - if (allocated(prec%precv)) then - - call psb_info(ictxt,me,np) - if (present(root)) then - root_ = root - else - root_ = psb_root_ - end if - if (root_ == -1) root_ = me - - ! - ! The preconditioner description is printed by processor psb_root_. - ! This agrees with the fact that all the parameters defining the - ! preconditioner have the same values on all the procs (this is - ! ensured by mld_precbld). - ! - if (me == root_) then - nlev = size(prec%precv) - do ilev = 1, nlev - if (.not.allocated(prec%precv(ilev)%sm)) then - info = 3111 - write(iout_,*) ' ',name,& - & ': error: inconsistent MLPREC part, should call MLD_PRECINIT' - return - endif - end do - - write(iout_,*) - write(iout_,'(a)') 'Preconditioner description' - if (nlev >= 1) then - ! - ! Print description of base preconditioner - ! - if (nlev > 1) then - write(iout_,*) 'Multilevel Preconditioner' - write(iout_,*) 'Outer sweeps:',prec%outer_sweeps - write(iout_,*) - write(iout_,*) 'Base preconditioner (smoother) details' - endif - call prec%precv(1)%sm%descr(info,iout=iout_) - if (nlev == 1) then - if (prec%precv(1)%parms%sweeps_pre > 1) then - write(iout_,*) ' Number of smoother sweeps_pre : ',& - & prec%precv(1)%parms%sweeps_pre - end if - write(iout_,*) - end if - if (allocated(prec%precv(1)%sm2a)) then - write(iout_,*) 'Post smoother details' - call prec%precv(1)%sm2a%descr(info,iout=iout_) - if (nlev == 1) then - if (prec%precv(1)%parms%sweeps_post > 1) then - write(iout_,*) ' Number of smoother sweeps_post : ',& - & prec%precv(1)%parms%sweeps_post - end if - write(iout_,*) - end if - end if - end if - - if (nlev > 1) then - ! - ! Print multilevel details - ! - write(iout_,*) - write(iout_,*) 'Multilevel details' - write(iout_,*) ' Number of levels : ',nlev - write(iout_,*) ' Operator complexity: ',prec%get_complexity() - ilmin = 2 - if (nlev == 2) ilmin=1 - do ilev=ilmin,nlev - call prec%precv(ilev)%descr(ilev,nlev,ilmin,info,iout=iout_) - end do - write(iout_,*) - end if - end if - - else - write(iout_,*) trim(name), & - & ': Error: no base preconditioner available, something is wrong!' - info = -2 - return - endif - - - end subroutine mld_sfile_prec_descr - - ! ! Subroutines: mld_Tprec_free ! Version: real diff --git a/mlprec/mld_z_as_smoother.f90 b/mlprec/mld_z_as_smoother.f90 index 7de01a94..24452e2f 100644 --- a/mlprec/mld_z_as_smoother.f90 +++ b/mlprec/mld_z_as_smoother.f90 @@ -461,7 +461,7 @@ contains implicit none character(len=32) :: val - val = "Schwarz smoother" + val = "Additive Schwarz" end function z_as_smoother_get_fmt function z_as_smoother_get_id() result(val) diff --git a/mlprec/mld_z_prec_type.f90 b/mlprec/mld_z_prec_type.f90 index 1f8e9442..df130c6f 100644 --- a/mlprec/mld_z_prec_type.f90 +++ b/mlprec/mld_z_prec_type.f90 @@ -156,7 +156,14 @@ module mld_z_prec_type interface mld_precdescr - module procedure mld_zfile_prec_descr + subroutine mld_zfile_prec_descr(prec,iout,root) + import :: mld_zprec_type, psb_ipk_ + implicit none + ! Arguments + class(mld_zprec_type), intent(in) :: prec + integer(psb_ipk_), intent(in), optional :: iout + integer(psb_ipk_), intent(in), optional :: root + end subroutine mld_zfile_prec_descr end interface interface mld_sizeof @@ -489,139 +496,6 @@ contains prec%op_complexity = num/den end subroutine mld_z_cmp_compl - ! - ! Subroutine: mld_file_prec_descr - ! Version: complex - ! - ! This routine prints a description of the preconditioner to the standard - ! output or to a file. It must be called after the preconditioner has been - ! built by mld_precbld. - ! - ! Arguments: - ! p - type(mld_Tprec_type), input. - ! The preconditioner data structure to be printed out. - ! info - integer, output. - ! error code. - ! iout - integer, input, optional. - ! The id of the file where the preconditioner description - ! will be printed. If iout is not present, then the standard - ! output is condidered. - ! root - integer, input, optional. - ! The id of the process printing the message; -1 acts as a wildcard. - ! Default is psb_root_ - ! - subroutine mld_zfile_prec_descr(prec,iout,root) - implicit none - ! Arguments - class(mld_zprec_type), intent(in) :: prec - integer(psb_ipk_), intent(in), optional :: iout - integer(psb_ipk_), intent(in), optional :: root - - ! Local variables - integer(psb_ipk_) :: ilev, nlev, ilmin, info - integer(psb_ipk_) :: ictxt, me, np - character(len=20), parameter :: name='mld_file_prec_descr' - integer(psb_ipk_) :: iout_ - integer(psb_ipk_) :: root_ - - info = psb_success_ - if (present(iout)) then - iout_ = iout - else - iout_ = psb_out_unit - end if - if (iout_ < 0) iout_ = psb_out_unit - - ictxt = prec%ictxt - - if (allocated(prec%precv)) then - - call psb_info(ictxt,me,np) - if (present(root)) then - root_ = root - else - root_ = psb_root_ - end if - if (root_ == -1) root_ = me - - ! - ! The preconditioner description is printed by processor psb_root_. - ! This agrees with the fact that all the parameters defining the - ! preconditioner have the same values on all the procs (this is - ! ensured by mld_precbld). - ! - if (me == root_) then - nlev = size(prec%precv) - do ilev = 1, nlev - if (.not.allocated(prec%precv(ilev)%sm)) then - info = 3111 - write(iout_,*) ' ',name,& - & ': error: inconsistent MLPREC part, should call MLD_PRECINIT' - return - endif - end do - - write(iout_,*) - write(iout_,'(a)') 'Preconditioner description' - if (nlev >= 1) then - ! - ! Print description of base preconditioner - ! - if (nlev > 1) then - write(iout_,*) 'Multilevel Preconditioner' - write(iout_,*) 'Outer sweeps:',prec%outer_sweeps - write(iout_,*) - write(iout_,*) 'Base preconditioner (smoother) details' - endif - call prec%precv(1)%sm%descr(info,iout=iout_) - if (nlev == 1) then - if (prec%precv(1)%parms%sweeps_pre > 1) then - write(iout_,*) ' Number of smoother sweeps_pre : ',& - & prec%precv(1)%parms%sweeps_pre - end if - write(iout_,*) - end if - if (allocated(prec%precv(1)%sm2a)) then - write(iout_,*) 'Post smoother details' - call prec%precv(1)%sm2a%descr(info,iout=iout_) - if (nlev == 1) then - if (prec%precv(1)%parms%sweeps_post > 1) then - write(iout_,*) ' Number of smoother sweeps_post : ',& - & prec%precv(1)%parms%sweeps_post - end if - write(iout_,*) - end if - end if - end if - - if (nlev > 1) then - ! - ! Print multilevel details - ! - write(iout_,*) - write(iout_,*) 'Multilevel details' - write(iout_,*) ' Number of levels : ',nlev - write(iout_,*) ' Operator complexity: ',prec%get_complexity() - ilmin = 2 - if (nlev == 2) ilmin=1 - do ilev=ilmin,nlev - call prec%precv(ilev)%descr(ilev,nlev,ilmin,info,iout=iout_) - end do - write(iout_,*) - end if - end if - - else - write(iout_,*) trim(name), & - & ': Error: no base preconditioner available, something is wrong!' - info = -2 - return - endif - - - end subroutine mld_zfile_prec_descr - - ! ! Subroutines: mld_Tprec_free ! Version: complex