fdurastante/amg4psblas
1
0
Fork 0
mirror of https://github.com/sfilippone/amg4psblas.git synced 2026-10-06 22:55:12 +00:00
Code Issues Packages Projects Releases Wiki Activity

First version working with BootCMatch. To be tested in detail.

Browse Source
This commit is contained in:
Salvatore Filippone
2018-09-17 09:48:05 +01:00
parent c7ab19ec9b
commit 660d00d49b
10 changed files with 142 additions and 46 deletions
Download Patch File Download Diff File Expand all files Collapse all files
mlprec/impl/level/mld_d_base_onelev_cseti.F90
+2 -1
View File
@@ -260,7 +260,8 @@ subroutine mld_d_base_onelev_cseti(lv,what,val,info,pos)
call lv%sm2a%set(what,val,info)
end if
end if
if (allocated(lv%aggr)) call lv%aggr%set(what,val,info)
end select
if (info /= psb_success_) goto 9999
call psb_erractionrestore(err_act)
mlprec/impl/level/mld_d_base_onelev_setag.f90
+1
View File
@@ -76,6 +76,7 @@ subroutine mld_d_base_onelev_setag(lv,val,info,pos)
lv%parms%par_aggr_alg = mld_ext_aggr_
lv%parms%aggr_type = mld_noalg_
end if
call lv%aggr%default()
end subroutine mld_d_base_onelev_setag
mlprec/mld_base_prec_type.F90
+1 -1
View File
@@ -571,7 +571,7 @@ contains
write(iout,*) ' Parallel aggregation algorithm: ',&
& par_aggr_alg_names(pm%par_aggr_alg)
write(iout,*) ' Aggregation type: ',&
if (pm%aggr_type>0) write(iout,*) ' Aggregation type: ',&
& aggr_type_names(pm%aggr_type)
if (pm%par_aggr_alg /= mld_ext_aggr_) then
if ( pm%aggr_ord /= mld_aggr_ord_nat_) &
mlprec/mld_d_base_aggregator_mod.f90
+17 -1
View File
@@ -104,12 +104,27 @@ module mld_d_base_aggregator_mod
procedure, pass(ag) :: default => mld_d_base_aggregator_default
procedure, pass(ag) :: descr => mld_d_base_aggregator_descr
procedure, pass(ag) :: set_aggr_type => mld_d_base_aggregator_set_aggr_type
procedure, pass(ag) :: cseti => mld_d_base_aggregator_cseti
generic, public :: set => cseti
procedure, nopass :: fmt => mld_d_base_aggregator_fmt
end type mld_d_base_aggregator_type
contains
subroutine mld_d_base_aggregator_cseti(ag,what,val,info)
Implicit None
! Arguments
class(mld_d_base_aggregator_type), intent(inout) :: ag
character(len=*), intent(in) :: what
integer(psb_ipk_), intent(in) :: val
integer(psb_ipk_), intent(out) :: info
! Do nothing
info = 0
end subroutine mld_d_base_aggregator_cseti
subroutine mld_d_base_aggregator_update_next(ag,agnext,info)
implicit none
class(mld_d_base_aggregator_type), target, intent(inout) :: ag, agnext
@@ -159,7 +174,7 @@ contains
implicit none
character(len=32) :: val
val = "Null "
val = "Default aggregator "
end function mld_d_base_aggregator_fmt
subroutine mld_d_base_aggregator_descr(ag,parms,iout,info)
@@ -169,6 +184,7 @@ contains
integer(psb_ipk_), intent(in) :: iout
integer(psb_ipk_), intent(out) :: info
write(iout,*) 'Aggregator object type: ',ag%fmt()
call parms%mldescr(iout,info)
return
tests/Bcmatch/Makefile
+5 -5
View File
@@ -6,20 +6,20 @@ MLDLIBDIR=$(MLDDIR)/lib
MLD_LIBS=-L$(MLDLIBDIR) -lpsb_krylov -lmld_prec -lpsb_prec
FINCLUDES=$(FMFLAG). $(FMFLAG)$(MLDMODDIR) $(FMFLAG)$(MLDINCDIR) $(PSBLAS_INCLUDES) $(FIFLAG).
HSL_DIR=/opt/hsl/2.3.1/gnu/6.4.0
HSL_DIR=/opt/hsl/2.3.1/sys
HSL_INCDIR=$(HSL_DIR)/include
HSL_LIBDIR=$(HSL_DIR)/lib
HSL_LIBS=-lhsl_mc64 -L$(HSL_LIBDIR)
HSL_FLAGS= -DHAVE_HSL -I$(HSL_INCDIR)
# SPRAL package for auction algorithm
SPRAL_DIR=/opt/spral/2015.04.20/gnu/6.4.0
SPRAL_DIR=/opt/spral/2015.04.20/sys
SPRAL_INCDIR=$(SPRAL_DIR)/include
SPRAL_LIBDIR=$(SPRAL_DIR)/lib
SPRAL_LIBS=-lspral -L$(SPRAL_LIBDIR)
SPRAL_FLAGS=-DHAVE_SPRAL -I$(SPRAL_INCDIR)
BCM_DIR=/opt/bcm/0.9/gnu/6.4.0
BCM_DIR=/opt/bcm/0.9/sys
BCM_INCDIR=$(BCM_DIR)/include
BCM_LIBDIR=$(BCM_DIR)/lib
BCM_LDLIBS=-lBCM -L$(BCM_LIBDIR) $(HSL_LIBS) $(SPRAL_LIBS)
@@ -47,8 +47,8 @@ check: all
clean:
/bin/rm -f data_input.o mld_d_pde3d.o mld_s_pde3d.o mld_d_pde2d.o mld_s_pde2d.o *$(.mod)\
$(EXEDIR)/mld_d_pde3d $(EXEDIR)/mld_s_pde3d $(EXEDIR)/mld_d_pde2d $(EXEDIR)/mld_s_pde2d
/bin/rm -f data_input.o mld_d_pde3d.o *$(.mod) $(BCMOBJS)\
$(EXEDIR)/mld_d_pde3d
verycleanlib:
(cd ../..; make veryclean)
tests/Bcmatch/bootCMatch_interface.c
+7
View File
@@ -1,3 +1,4 @@
#include <string.h>
#include <stdio.h>
@@ -38,6 +39,12 @@ int mld_bootCMatch_if(bcm_CSRMatrix *C, int match_algorithm, int n_sweeps,
int ftcoarse=1;
int cr_it=0, cr_relax_type=0;
double cr_relax_weight=0.0;
// Sanity checks
nr = bcm_CSRMatrixNumRows(C);
nc = bcm_VectorSize(w);
// fprintf(stderr,"Sanity check: %d %d \n",nr,nc);
// Here I am building Ac but I won't use it.
Ac=bcm_CSRMatchingAgg(C, &w, &P, match_algorithm, n_sweeps, max_nlevels,max_csize , &ftcoarse,
cr_it, cr_relax_type, cr_relax_weight);
tests/Bcmatch/mld_d_bcmatch_aggregator_mod.F90
+86 -35
View File
@@ -113,17 +113,20 @@ module mld_d_bcmatch_aggregator_mod
type, extends(mld_d_base_aggregator_type) :: mld_d_bcmatch_aggregator_type
integer(psb_ipk_) :: matching_alg
integer(psb_ipk_) :: n_sweeps
real(psb_dpk_), allocatable :: w_tmp(:)
type(bcm_Vector) :: w_par
real(psb_dpk_), allocatable :: w_tmp(:), w_nxt(:)
type(bcm_Vector) :: w_par
integer(psb_ipk_) :: max_csize
integer(psb_ipk_) :: max_nlevels
!type(psb_d_vect_type) :: w
contains
procedure, pass(ag) :: bld_tprol => mld_d_bcmatch_aggregator_build_tprol
procedure, pass(ag) :: set => d_bcmatch_aggr_cseti
procedure, pass(ag) :: cseti => d_bcmatch_aggr_cseti
procedure, pass(ag) :: default => d_bcmatch_aggr_set_default
procedure, pass(ag) :: mat_asb => mld_d_bcmatch_aggregator_mat_asb
procedure, pass(ag) :: update_level => d_bcmatch_aggregator_update_level
procedure, pass(ag) :: update_next => d_bcmatch_aggregator_update_next
procedure, pass(ag) :: bld_wnxt => d_bcmatch_bld_wnxt
procedure, pass(ag) :: bld_default_w => d_bld_default_w
procedure, pass(ag) :: set_c_default_w => d_set_default_bcm_w
!!$ procedure, pass(ag) :: clone => mld_d_base_aggregator_clone
!!$ procedure, pass(ag) :: free => mld_d_bcmatch_aggregator_free
!!$ procedure, pass(ag) :: default => mld_d_base_aggregator_default
@@ -195,7 +198,60 @@ module mld_d_bcmatch_aggregator_mod
contains
subroutine d_bld_default_w(ag,nr)
use psb_realloc_mod
implicit none
class(mld_d_bcmatch_aggregator_type), target, intent(inout) :: ag
integer(psb_ipk_), intent(in) :: nr
integer(psb_ipk_) :: info
call psb_realloc(nr,ag%w_tmp,info)
if (info /= psb_success_) return
ag%w_tmp = done
call ag%set_c_default_w()
end subroutine d_bld_default_w
subroutine d_set_default_bcm_w(ag)
use psb_realloc_mod
use iso_c_binding
implicit none
class(mld_d_bcmatch_aggregator_type), target, intent(inout) :: ag
ag%w_par%size = psb_size(ag%w_tmp)
ag%w_par%owns_data = 0
if (ag%w_par%size > 0) call set_cloc(ag%w_tmp, ag%w_par)
end subroutine d_set_default_bcm_w
subroutine set_cloc(vect,w_par)
use iso_c_binding
real(psb_dpk_), target :: vect(:)
type(bcm_Vector) :: w_par
w_par%data = c_loc(vect)
end subroutine set_cloc
subroutine d_bcmatch_bld_wnxt(ag,ilaggr,valaggr,nx)
use psb_realloc_mod
implicit none
class(mld_d_bcmatch_aggregator_type), target, intent(inout) :: ag
integer(psb_ipk_), intent(in) :: ilaggr(:)
real(psb_dpk_), intent(in) :: valaggr(:)
integer(psb_ipk_), intent(in) :: nx
integer(psb_ipk_) :: info,i,j
call psb_realloc(nx,ag%w_nxt,info)
associate(w_nxt => ag%w_nxt, w_tmp=>ag%w_tmp)
w_nxt = dzero
do j=1, size(ilaggr)
i = ilaggr(j)
w_nxt(i) = w_nxt(i) + valaggr(j)*w_tmp(j)
end do
end associate
end subroutine d_bcmatch_bld_wnxt
function mld_d_bcmatch_aggregator_fmt() result(val)
implicit none
character(len=32) :: val
@@ -203,7 +259,8 @@ contains
val = "BootCMatch aggregation"
end function mld_d_bcmatch_aggregator_fmt
subroutine d_bcmatch_aggregator_update_level(ag,agnext,info)
subroutine d_bcmatch_aggregator_update_next(ag,agnext,info)
use psb_realloc_mod
implicit none
class(mld_d_bcmatch_aggregator_type), target, intent(inout) :: ag
class(mld_d_base_aggregator_type), target, intent(inout) :: agnext
@@ -212,17 +269,20 @@ contains
!
!
select type(agnext)
type is (mld_d_bcmatch_aggregator_type)
class is (mld_d_bcmatch_aggregator_type)
agnext%matching_alg = ag%matching_alg
agnext%n_sweeps = ag%n_sweeps
agnext%max_csize = ag%max_csize
agnext%max_nlevels = ag%max_nlevels
! Is this going to generate shallow copies/memory leaks/double frees?
! To be investigated further.
agnext%w_par = ag%w_par
call psb_safe_ab_cpy(ag%w_nxt,agnext%w_tmp,info)
call agnext%set_c_default_w()
class default
! What should we do here?
end select
info = 0
end subroutine d_bcmatch_aggregator_update_level
end subroutine d_bcmatch_aggregator_update_next
subroutine d_bcmatch_aggr_cseti(ag,what,val,info)
@@ -238,31 +298,22 @@ contains
info = psb_success_
select case(what)
case('BCM_MATCH_ALG')
ag%matching_alg=val
case('BCM_SWEEPS')
ag%n_sweeps=val
case('BCM_MAX_CSIZE')
ag%max_csize=val
case('BCM_MAX_NLEVELS')
ag%max_nlevels=val
case('BCM_W_SIZE')
ag%w_par%size=val
ag%w_par%owns_data=0
allocate(ag%w_tmp(val))
ag%w_tmp = 1.0_psb_dpk_
call set_cloc(ag%w_tmp, ag%w_par)
case default
case('BCM_MATCH_ALG')
ag%matching_alg=val
case('BCM_SWEEPS')
ag%n_sweeps=val
case('BCM_MAX_CSIZE')
ag%max_csize=val
case('BCM_MAX_NLEVELS')
ag%max_nlevels=val
case('BCM_W_SIZE')
!write(0,*) 'Setting W_SIZE'
call ag%bld_default_w(val)
case default
end select
return
contains
subroutine set_cloc(vect,w_par)
real(psb_dpk_), target :: vect(:)
type(bcm_Vector) :: w_par
w_par%data = c_loc(vect)
end subroutine set_cloc
end subroutine d_bcmatch_aggr_cseti
subroutine d_bcmatch_aggr_set_default(ag)
@@ -272,10 +323,10 @@ contains
! Arguments
class(mld_d_bcmatch_aggregator_type), intent(inout) :: ag
character(len=20) :: name='d_bcmatch_aggr_set_default'
ag%matching_alg=0
ag%n_sweeps=1
ag%max_nlevels=36
ag%max_csize=10
ag%matching_alg = 0
ag%n_sweeps = 1
ag%max_nlevels = 36
ag%max_csize = 10
return
tests/Bcmatch/mld_d_bcmatch_aggregator_tprol.f90
+8 -1
View File
@@ -292,7 +292,9 @@ subroutine mld_d_bcmatch_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr
call a%csclip(b=a_tmp, info=info, jmax=a%get_nrows(), imax=a%get_nrows())
call a_tmp%mv_to(acsr)
nr = a%get_nrows()
if (psb_size(ag%w_tmp) < nr) call ag%bld_default_w(nr)
!write(*,*) 'Build_tprol:',acsr%get_nrows(),acsr%get_ncols()
C%num_rows=acsr%get_nrows()
C%num_cols=acsr%get_ncols()
@@ -326,7 +328,12 @@ subroutine mld_d_bcmatch_aggregator_build_tprol(ag,parms,a,desc_a,ilaggr,nlaggr
call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_bootCMatch_if')
goto 9999
end if
!!$ write(0,*) 'On output from BootCMatch',nr,num_pcols,size(ilaggr),maxval(ilaggr),&
!!$ & minval(ilaggr),minval(ilaggr(1:nr)),a%get_nrows(),a%get_ncols()
! Prepare vector W for next level, just in case
call ag%bld_wnxt(ilaggr(1:nr),valaggr(1:nr),num_pcols)
call psb_realloc(np,nlaggr,info)
if (info /= psb_success_) then
info=psb_err_alloc_request_
tests/Bcmatch/mld_d_pde3d.f90
+14 -2
View File
@@ -671,6 +671,7 @@ program mld_d_pde3d
integer(psb_ipk_) :: cfill ! fill-in for incomplete LU factorization
real(psb_dpk_) :: cthres ! threshold for ILUT factorization
integer(psb_ipk_) :: cjswp ! sweeps for GS or JAC coarsest-lev subsolver
logical :: use_bcm ! Use BootCMatch aggregation
end type precdata
type(precdata) :: p_choice
@@ -816,13 +817,23 @@ program mld_d_pde3d
call prec%set('coarse_fillin', p_choice%cfill, info)
call prec%set('coarse_iluthrs', p_choice%cthres, info)
call prec%set('coarse_sweeps', p_choice%cjswp, info)
!call prec%set(bcmag,info)
if (p_choice%use_bcm) then
call prec%set(bcmag,info)
call prec%set('BCM_MATCH_ALG',2, info)
call prec%set('BCM_SWEEPS',3, info)
!!$ if (p_choice%csize>0) call prec%set('BCM_MAX_CSIZE',p_choice%csize, info)
call prec%set('BCM_MAX_NLEVELS',p_choice%maxlevs, info)
!call prec%set('BCM_W_SIZE',desc_a%get_local_rows(), info,ilev=2)
end if
end select
! build the preconditioner
call psb_barrier(ictxt)
t1 = psb_wtime()
!call psb_set_debug_level(9999)
call prec%hierarchy_build(a,desc_a,info)
!call psb_set_debug_level(0)
thier = psb_wtime()-t1
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='mld_hierarchy_bld')
@@ -1024,6 +1035,7 @@ contains
call read_data(prec%cfill,inp_unit) ! fill-in for incompl LU
call read_data(prec%cthres,inp_unit) ! Threshold for ILUT
call read_data(prec%cjswp,inp_unit) ! sweeps for GS/JAC subsolver
call read_data(prec%use_bcm,inp_unit) ! BootCMatch?
if (inp_unit /= psb_inp_unit) then
close(inp_unit)
end if
@@ -1085,7 +1097,7 @@ contains
call psb_bcast(icontxt,prec%cfill)
call psb_bcast(icontxt,prec%cthres)
call psb_bcast(icontxt,prec%cjswp)
call psb_bcast(icontxt,prec%use_bcm)
end subroutine get_parms
tests/Bcmatch/runs/mld-bcm.inp
+1
View File
@@ -47,3 +47,4 @@ DIST ! Coarsest-level matrix distribution: DIST REPL, DE
1 ! Coarsest-level fillin P for ILU(P) and ILU(T,P)
1.d-4 ! Coarsest-level threshold T for ILU(T,P)
1 ! Number of sweeps for JACOBI/GS/BJAC coarsest-level solver
T ! Use BootCMatch aggregator