From 80c58b32eb08ab5b86c99df23de2907d92acd684 Mon Sep 17 00:00:00 2001 From: Salvatore Filippone Date: Sun, 8 May 2016 19:41:25 +0000 Subject: [PATCH] mld2p4-smooth-2side: mlprec/impl/mld_dprecset.F90 mlprec/mld_d_prec_mod.f90 mlprec/mld_d_prec_type.f90 tests/pdegen/Makefile tests/pdegen/runs/ppde.inp Fix dec_map XZERO. --- mlprec/impl/mld_dprecset.F90 | 181 ++++++++++++++++++++++++++--------- mlprec/mld_d_prec_mod.f90 | 29 +++--- mlprec/mld_d_prec_type.f90 | 38 ++++++-- tests/pdegen/Makefile | 12 ++- tests/pdegen/runs/ppde.inp | 6 +- 5 files changed, 192 insertions(+), 74 deletions(-) diff --git a/mlprec/impl/mld_dprecset.F90 b/mlprec/impl/mld_dprecset.F90 index 2f328351..730afd4f 100644 --- a/mlprec/impl/mld_dprecset.F90 +++ b/mlprec/impl/mld_dprecset.F90 @@ -76,7 +76,7 @@ ! For this reason, the interface mld_precset to this routine has been built in ! such a way that ilev is not visible to the user (see mld_prec_mod.f90). ! -subroutine mld_dprecseti(p,what,val,info,ilev) +subroutine mld_dprecseti(p,what,val,info,ilev,pos) use psb_base_mod use mld_d_prec_mod, mld_protect_name => mld_dprecseti @@ -107,6 +107,7 @@ subroutine mld_dprecseti(p,what,val,info,ilev) integer(psb_ipk_), intent(in) :: val integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), optional, intent(in) :: ilev + character(len=*), optional, intent(in) :: pos ! Local variables integer(psb_ipk_) :: ilev_, nlev_ @@ -671,7 +672,7 @@ contains end subroutine mld_dprecseti -subroutine mld_dprecsetsm(p,val,info,ilev) +subroutine mld_dprecsetsm(p,val,info,ilev,pos) use psb_base_mod use mld_d_prec_mod, mld_protect_name => mld_dprecsetsm @@ -679,15 +680,16 @@ subroutine mld_dprecsetsm(p,val,info,ilev) implicit none ! Arguments - class(mld_dprec_type), intent(inout) :: p - class(mld_d_base_smoother_type), intent(in) :: val - integer(psb_ipk_), intent(out) :: info - integer(psb_ipk_), optional, intent(in) :: ilev - + class(mld_dprec_type), target, intent(inout) :: p + class(mld_d_base_smoother_type), intent(in) :: val + integer(psb_ipk_), intent(out) :: info + integer(psb_ipk_), optional, intent(in) :: ilev + character(len=*), optional, intent(in) :: pos + ! Local variables - integer(psb_ipk_) :: ilev_, nlev_, ilmin, ilmax + integer(psb_ipk_) :: ilev_, nlev_, ilmin, ilmax, ipos_ character(len=*), parameter :: name='mld_precseti' - + info = psb_success_ if (.not.allocated(p%precv)) then @@ -708,7 +710,19 @@ subroutine mld_dprecsetsm(p,val,info,ilev) ilmin = 1 ilmax = nlev_ end if - + if (present(pos)) then + select case(psb_toupper(trim(pos))) + case('PRE') + ipos_ = mld_pre_smooth_ + case('POST') + ipos_ = mld_post_smooth_ + case default + ipos_ = mld_pre_smooth_ + end select + else + ipos_ = mld_pre_smooth_ + end if + if ((ilev_<1).or.(ilev_ > nlev_)) then info = -1 write(psb_err_unit,*) name,& @@ -716,25 +730,44 @@ subroutine mld_dprecsetsm(p,val,info,ilev) return endif - - do ilev_ = ilmin, ilmax - if (allocated(p%precv(ilev_)%sm)) then - if (allocated(p%precv(ilev_)%sm%sv)) then - deallocate(p%precv(ilev_)%sm%sv) - endif - deallocate(p%precv(ilev_)%sm) - end if + select case(ipos_) + case(mld_pre_smooth_) + do ilev_ = ilmin, ilmax + if (allocated(p%precv(ilev_)%sm)) then + if (allocated(p%precv(ilev_)%sm%sv)) then + deallocate(p%precv(ilev_)%sm%sv) + endif + deallocate(p%precv(ilev_)%sm) + end if #ifdef HAVE_MOLD - allocate(p%precv(ilev_)%sm,mold=val) + allocate(p%precv(ilev_)%sm,mold=val) #else - allocate(p%precv(ilev_)%sm,source=val) + allocate(p%precv(ilev_)%sm,source=val) #endif - call p%precv(ilev_)%sm%default() - end do - + call p%precv(ilev_)%sm%default() + p%precv(ilev_)%sm2 => p%precv(ilev_)%sm + end do + case(mld_post_smooth_) + do ilev_ = ilmin, ilmax + if (allocated(p%precv(ilev_)%sm2a)) then + if (allocated(p%precv(ilev_)%sm2a%sv)) then + deallocate(p%precv(ilev_)%sm2a%sv) + endif + deallocate(p%precv(ilev_)%sm2a) + end if +#ifdef HAVE_MOLD + allocate(p%precv(ilev_)%sm2a,mold=val) +#else + allocate(p%precv(ilev_)%sm2a,source=val) +#endif + call p%precv(ilev_)%sm2a%default() + p%precv(ilev_)%sm2 => p%precv(ilev_)%sm2a + end do + end select + end subroutine mld_dprecsetsm -subroutine mld_dprecsetsv(p,val,info,ilev) +subroutine mld_dprecsetsv(p,val,info,ilev,pos) use psb_base_mod use mld_d_prec_mod, mld_protect_name => mld_dprecsetsv @@ -746,10 +779,11 @@ subroutine mld_dprecsetsv(p,val,info,ilev) class(mld_d_base_solver_type), intent(in) :: val integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), optional, intent(in) :: ilev - + character(len=*), optional, intent(in) :: pos + ! Local variables - integer(psb_ipk_) :: ilev_, nlev_, ilmin, ilmax - character(len=*), parameter :: name='mld_precseti' + integer(psb_ipk_) :: ilev_, nlev_, ilmin, ilmax, ipos_ + character(len=*), parameter :: name='mld_precseti' info = psb_success_ @@ -772,6 +806,19 @@ subroutine mld_dprecsetsv(p,val,info,ilev) ilmax = nlev_ end if + if (present(pos)) then + select case(psb_toupper(trim(pos))) + case('PRE') + ipos_ = mld_pre_smooth_ + case('POST') + ipos_ = mld_post_smooth_ + case default + ipos_ = mld_pre_smooth_ + end select + else + ipos_ = mld_pre_smooth_ + end if + if ((ilev_<1).or.(ilev_ > nlev_)) then info = -1 @@ -781,16 +828,20 @@ subroutine mld_dprecsetsv(p,val,info,ilev) endif - do ilev_ = ilmin, ilmax - if (allocated(p%precv(ilev_)%sm)) then - if (allocated(p%precv(ilev_)%sm%sv)) then - if (.not.same_type_as(p%precv(ilev_)%sm%sv,val)) then - deallocate(p%precv(ilev_)%sm%sv,stat=info) - if (info /= 0) then - info = 3111 - return + select case(ipos_) + case(mld_pre_smooth_) + do ilev_ = ilmin, ilmax + if (allocated(p%precv(ilev_)%sm)) then + if (allocated(p%precv(ilev_)%sm%sv)) then + if (.not.same_type_as(p%precv(ilev_)%sm%sv,val)) then + deallocate(p%precv(ilev_)%sm%sv,stat=info) + if (info /= 0) then + info = 3111 + return + end if end if end if + if (.not.allocated(p%precv(ilev_)%sm%sv)) then #ifdef HAVE_MOLD allocate(p%precv(ilev_)%sm%sv,mold=val,stat=info) @@ -802,19 +853,55 @@ subroutine mld_dprecsetsv(p,val,info,ilev) return end if end if + call p%precv(ilev_)%sm%sv%default() + else + info = 3111 + write(psb_err_unit,*) name,& + &': Error: uninitialized preconditioner component,',& + &' should call MLD_PRECINIT/MLD_PRECSET' + return + end if - call p%precv(ilev_)%sm%sv%default() - else - info = 3111 - write(psb_err_unit,*) name,& - &': Error: uninitialized preconditioner component,',& - &' should call MLD_PRECINIT/MLD_PRECSET' - return - - end if - - end do - + + end do + + case(mld_post_smooth_) + do ilev_ = ilmin, ilmax + if (allocated(p%precv(ilev_)%sm2a)) then + if (allocated(p%precv(ilev_)%sm2a%sv)) then + write(0,*)p%precv(ilev_)%sm2a%sv%get_fmt(),val%get_fmt() + if (.not.same_type_as(p%precv(ilev_)%sm2a%sv,val)) then + deallocate(p%precv(ilev_)%sm2a%sv,stat=info) + if (info /= 0) then + info = 3111 + return + end if + end if + end if + if (.not.allocated(p%precv(ilev_)%sm2a%sv)) then +#ifdef HAVE_MOLD + allocate(p%precv(ilev_)%sm2a%sv,mold=val,stat=info) +#else + allocate(p%precv(ilev_)%sm2a%sv,source=val,stat=info) +#endif + if (info /= 0) then + info = 3111 + return + end if + end if + call p%precv(ilev_)%sm2a%sv%default() + write(0,*)p%precv(ilev_)%sm2a%sv%get_fmt(),val%get_fmt() + else + info = 3111 + write(psb_err_unit,*) name,& + &': Error: uninitialized preconditioner component,',& + &' should call MLD_PRECINIT/MLD_PRECSET' + return + + end if + + end do + end select end subroutine mld_dprecsetsv diff --git a/mlprec/mld_d_prec_mod.f90 b/mlprec/mld_d_prec_mod.f90 index a9474b96..dd6195b1 100644 --- a/mlprec/mld_d_prec_mod.f90 +++ b/mlprec/mld_d_prec_mod.f90 @@ -91,72 +91,79 @@ module mld_d_prec_mod contains - subroutine mld_d_iprecsetsm(p,val,info) + subroutine mld_d_iprecsetsm(p,val,info,pos) type(mld_dprec_type), intent(inout) :: p class(mld_d_base_smoother_type), intent(in) :: val integer(psb_ipk_), intent(out) :: info + character(len=*), optional, intent(in) :: pos - call p%set(val,info) + call p%set(val,info,pos=pos) end subroutine mld_d_iprecsetsm - subroutine mld_d_iprecsetsv(p,val,info) + subroutine mld_d_iprecsetsv(p,val,info,pos) type(mld_dprec_type), intent(inout) :: p class(mld_d_base_solver_type), intent(in) :: val integer(psb_ipk_), intent(out) :: info - - call p%set(val,info) + character(len=*), optional, intent(in) :: pos + call p%set(val,info, pos=pos) end subroutine mld_d_iprecsetsv - subroutine mld_d_iprecseti(p,what,val,info) + subroutine mld_d_iprecseti(p,what,val,info,pos) type(mld_dprec_type), intent(inout) :: p integer(psb_ipk_), intent(in) :: what integer(psb_ipk_), intent(in) :: val integer(psb_ipk_), intent(out) :: info + character(len=*), optional, intent(in) :: pos call p%set(what,val,info) end subroutine mld_d_iprecseti - subroutine mld_d_iprecsetr(p,what,val,info) + subroutine mld_d_iprecsetr(p,what,val,info,pos) type(mld_dprec_type), intent(inout) :: p integer(psb_ipk_), intent(in) :: what real(psb_dpk_), intent(in) :: val integer(psb_ipk_), intent(out) :: info + character(len=*), optional, intent(in) :: pos call p%set(what,val,info) end subroutine mld_d_iprecsetr - subroutine mld_d_iprecsetc(p,what,val,info) + subroutine mld_d_iprecsetc(p,what,val,info,pos) type(mld_dprec_type), intent(inout) :: p integer(psb_ipk_), intent(in) :: what character(len=*), intent(in) :: val integer(psb_ipk_), intent(out) :: info + character(len=*), optional, intent(in) :: pos call p%set(what,val,info) end subroutine mld_d_iprecsetc - subroutine mld_d_cprecseti(p,what,val,info) + subroutine mld_d_cprecseti(p,what,val,info,pos) type(mld_dprec_type), intent(inout) :: p character(len=*), intent(in) :: what integer(psb_ipk_), intent(in) :: val integer(psb_ipk_), intent(out) :: info + character(len=*), optional, intent(in) :: pos call p%set(what,val,info) end subroutine mld_d_cprecseti - subroutine mld_d_cprecsetr(p,what,val,info) + subroutine mld_d_cprecsetr(p,what,val,info,pos) type(mld_dprec_type), intent(inout) :: p character(len=*), intent(in) :: what real(psb_dpk_), intent(in) :: val integer(psb_ipk_), intent(out) :: info + character(len=*), optional, intent(in) :: pos call p%set(what,val,info) end subroutine mld_d_cprecsetr - subroutine mld_d_cprecsetc(p,what,val,info) + subroutine mld_d_cprecsetc(p,what,val,info,pos) type(mld_dprec_type), intent(inout) :: p character(len=*), intent(in) :: what character(len=*), intent(in) :: val integer(psb_ipk_), intent(out) :: info + character(len=*), optional, intent(in) :: pos call p%set(what,val,info) end subroutine mld_d_cprecsetc diff --git a/mlprec/mld_d_prec_type.f90 b/mlprec/mld_d_prec_type.f90 index 765f6999..3f5cd081 100644 --- a/mlprec/mld_d_prec_type.f90 +++ b/mlprec/mld_d_prec_type.f90 @@ -175,23 +175,25 @@ module mld_d_prec_type end interface interface - subroutine mld_dprecsetsm(prec,val,info,ilev) + subroutine mld_dprecsetsm(prec,val,info,ilev,pos) import :: psb_dspmat_type, psb_desc_type, psb_dpk_, & & mld_dprec_type, mld_d_base_smoother_type, psb_ipk_ - class(mld_dprec_type), intent(inout) :: prec + class(mld_dprec_type), target, intent(inout) :: prec class(mld_d_base_smoother_type), intent(in) :: val integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), optional, intent(in) :: ilev + character(len=*), optional, intent(in) :: pos end subroutine mld_dprecsetsm - subroutine mld_dprecsetsv(prec,val,info,ilev) + subroutine mld_dprecsetsv(prec,val,info,ilev,pos) import :: psb_dspmat_type, psb_desc_type, psb_dpk_, & & mld_dprec_type, mld_d_base_solver_type, psb_ipk_ class(mld_dprec_type), intent(inout) :: prec class(mld_d_base_solver_type), intent(in) :: val integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), optional, intent(in) :: ilev + character(len=*), optional, intent(in) :: pos end subroutine mld_dprecsetsv - subroutine mld_dprecseti(prec,what,val,info,ilev) + subroutine mld_dprecseti(prec,what,val,info,ilev,pos) import :: psb_dspmat_type, psb_desc_type, psb_dpk_, & & mld_dprec_type, psb_ipk_ class(mld_dprec_type), intent(inout) :: prec @@ -199,8 +201,9 @@ module mld_d_prec_type integer(psb_ipk_), intent(in) :: val integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), optional, intent(in) :: ilev + character(len=*), optional, intent(in) :: pos end subroutine mld_dprecseti - subroutine mld_dprecsetr(prec,what,val,info,ilev) + subroutine mld_dprecsetr(prec,what,val,info,ilev,pos) import :: psb_dspmat_type, psb_desc_type, psb_dpk_, & & mld_dprec_type, psb_ipk_ class(mld_dprec_type), intent(inout) :: prec @@ -208,8 +211,9 @@ module mld_d_prec_type real(psb_dpk_), intent(in) :: val integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), optional, intent(in) :: ilev + character(len=*), optional, intent(in) :: pos end subroutine mld_dprecsetr - subroutine mld_dprecsetc(prec,what,string,info,ilev) + subroutine mld_dprecsetc(prec,what,string,info,ilev,pos) import :: psb_dspmat_type, psb_desc_type, psb_dpk_, & & mld_dprec_type, psb_ipk_ class(mld_dprec_type), intent(inout) :: prec @@ -217,8 +221,9 @@ module mld_d_prec_type character(len=*), intent(in) :: string integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), optional, intent(in) :: ilev + character(len=*), optional, intent(in) :: pos end subroutine mld_dprecsetc - subroutine mld_dcprecseti(prec,what,val,info,ilev) + subroutine mld_dcprecseti(prec,what,val,info,ilev,pos) import :: psb_dspmat_type, psb_desc_type, psb_dpk_, & & mld_dprec_type, psb_ipk_ class(mld_dprec_type), intent(inout) :: prec @@ -226,8 +231,9 @@ module mld_d_prec_type integer(psb_ipk_), intent(in) :: val integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), optional, intent(in) :: ilev + character(len=*), optional, intent(in) :: pos end subroutine mld_dcprecseti - subroutine mld_dcprecsetr(prec,what,val,info,ilev) + subroutine mld_dcprecsetr(prec,what,val,info,ilev,pos) import :: psb_dspmat_type, psb_desc_type, psb_dpk_, & & mld_dprec_type, psb_ipk_ class(mld_dprec_type), intent(inout) :: prec @@ -235,8 +241,9 @@ module mld_d_prec_type real(psb_dpk_), intent(in) :: val integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), optional, intent(in) :: ilev + character(len=*), optional, intent(in) :: pos end subroutine mld_dcprecsetr - subroutine mld_dcprecsetc(prec,what,string,info,ilev) + subroutine mld_dcprecsetc(prec,what,string,info,ilev,pos) import :: psb_dspmat_type, psb_desc_type, psb_dpk_, & & mld_dprec_type, psb_ipk_ class(mld_dprec_type), intent(inout) :: prec @@ -244,6 +251,7 @@ module mld_d_prec_type character(len=*), intent(in) :: string integer(psb_ipk_), intent(out) :: info integer(psb_ipk_), optional, intent(in) :: ilev + character(len=*), optional, intent(in) :: pos end subroutine mld_dcprecsetc end interface @@ -478,6 +486,18 @@ contains write(iout_,*) return end if + if (allocated(p%precv(1)%sm2a)) then + write(iout_,*) 'Post smoother details' + call p%precv(1)%sm2a%descr(info,iout=iout_) + if (nlev == 1) then + if (p%precv(1)%parms%sweeps > 1) then + write(iout_,*) ' Number of smoother sweeps : ',& + & p%precv(1)%parms%sweeps + end if + write(iout_,*) + return + end if + end if end if ! diff --git a/tests/pdegen/Makefile b/tests/pdegen/Makefile index 7bbe3439..66a0b906 100644 --- a/tests/pdegen/Makefile +++ b/tests/pdegen/Makefile @@ -11,12 +11,16 @@ FINCLUDES=$(FMFLAG). $(FMFLAG)$(MLDINCDIR) $(FMFLAG)$(PSBINCDIR) $(FIFLAG). EXEDIR=./runs -all: spde3d ppde3d spde2d ppde2d +all: spde3d ppde3d spde2d ppde2d ppde3d-gs ppde3d: ppde3d.o data_input.o $(F90LINK) ppde3d.o data_input.o -o ppde3d $(MLD_LIB) $(PSBLAS_LIB) $(LDLIBS) /bin/mv ppde3d $(EXEDIR) +ppde3d-gs: ppde3d-gs.o data_input.o + $(F90LINK) ppde3d-gs.o data_input.o -o ppde3d-gs $(MLD_LIB) $(PSBLAS_LIB) $(LDLIBS) + /bin/mv ppde3d-gs $(EXEDIR) + spde3d: spde3d.o data_input.o $(F90LINK) spde3d.o data_input.o -o spde3d $(MLD_LIB) $(PSBLAS_LIB) $(LDLIBS) /bin/mv spde3d $(EXEDIR) @@ -30,15 +34,15 @@ spde2d: spde2d.o data_input.o $(F90LINK) spde2d.o data_input.o -o spde2d $(MLD_LIB) $(PSBLAS_LIB) $(LDLIBS) /bin/mv spde2d $(EXEDIR) -ppde3d.o spde3d.o ppde2d.o spde2d.o: data_input.o +ppde3d-gs.o ppde3d.o spde3d.o ppde2d.o spde2d.o: data_input.o check: all cd runs && ./ppde2d