More work on interfacing with new OpenMP matching

TestFerdous
Salvatore Filippone 3 years ago
parent 6f419a2210
commit 9779eeec9f

@ -17,7 +17,7 @@ DMODOBJS=amg_d_prec_type.o \
amg_d_ainv_solver.o amg_d_base_ainv_mod.o \
amg_d_invk_solver.o amg_d_invt_solver.o amg_d_krm_solver.o \
amg_d_matchboxp_mod.o amg_d_parmatch_aggregator_mod.o \
amg_d_newmatch_aggregator_mod.o
amg_d_newmatch_aggregator_mod.o amg_d_decmatch_mod.o
SMODOBJS=amg_s_prec_type.o amg_s_ilu_fact_mod.o \
amg_s_inner_mod.o amg_s_ilu_solver.o amg_s_diag_solver.o amg_s_jac_smoother.o amg_s_as_smoother.o \
@ -127,7 +127,8 @@ amg_d_base_aggregator_mod.o: amg_base_prec_type.o
amg_d_parmatch_aggregator_mod.o amg_d_dec_aggregator_mod.o: amg_d_base_aggregator_mod.o
amg_d_hybrid_aggregator_mod.o amg_d_symdec_aggregator_mod.o: amg_d_dec_aggregator_mod.o
amg_d_parmatch_aggregator_mod.o: amg_d_matchboxp_mod.o
amg_d_newmatch_aggregator_mod.o: amg_d_dec_aggregator_mod.o
amg_d_newmatch_aggregator_mod.o: amg_d_base_aggregator_mod.o
amg_d_newmatch_aggregator_mod.o: amg_d_decmatch_mod.o
amg_c_base_aggregator_mod.o: amg_base_prec_type.o
amg_c_parmatch_aggregator_mod.o amg_c_dec_aggregator_mod.o: amg_c_base_aggregator_mod.o

@ -73,16 +73,25 @@ module amg_d_newmatch_aggregator_mod
! the W argument with its update. Hence, copy it in w_nxt
! before passing it to the matching
!
integer(psb_ipk_) :: orig_aggr_size
integer(psb_ipk_) :: jacobi_sweeps
real(psb_dpk_), allocatable :: w(:), w_nxt(:)
type(psb_dspmat_type), allocatable :: prol, restr
type(psb_dspmat_type), allocatable :: ac, base_a, rwa
type(psb_desc_type), allocatable :: desc_ac, desc_ax, base_desc, rwdesc
type(nwm_Vector) :: w_c_nxt
integer(psb_ipk_) :: max_csize
integer(psb_ipk_) :: max_nlevels
logical :: reproducible_matching = .false.
logical :: need_symmetrize = .false.
logical :: unsmoothed_hierarchy = .true.
contains
procedure, pass(ag) :: bld_tprol => amg_d_newmatch_aggregator_build_tprol
procedure, pass(ag) :: cseti => d_newmatch_aggr_cseti
procedure, pass(ag) :: default => d_newmatch_aggr_set_default
procedure, pass(ag) :: mat_asb => amg_d_newmatch_aggregator_mat_asb
procedure, pass(ag) :: mat_bld => amg_d_newmatch_aggregator_mat_bld
procedure, pass(ag) :: inner_mat_asb => amg_d_newmatch_aggregator_inner_mat_asb
procedure, pass(ag) :: update_next => d_newmatch_aggregator_update_next
procedure, pass(ag) :: bld_wnxt => d_newmatch_bld_wnxt
procedure, pass(ag) :: bld_default_w => d_bld_default_w
@ -179,6 +188,24 @@ module amg_d_newmatch_aggregator_mod
integer(psb_ipk_), intent(out) :: info
end subroutine amg_daggrmat_unsmth_spmm_asb
end interface
interface
subroutine amg_d_newmatch_aggregator_inner_mat_asb(ag,parms,a,desc_a,&
& ac,desc_ac, op_prol,op_restr,info)
import :: amg_d_newmatch_aggregator_type, psb_desc_type, psb_dspmat_type,&
& psb_ldspmat_type, psb_dpk_, psb_ipk_, psb_lpk_, amg_dml_parms, amg_daggr_data
implicit none
class(amg_d_newmatch_aggregator_type), target, intent(inout) :: ag
type(amg_dml_parms), intent(inout) :: parms
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dspmat_type), intent(inout) :: op_prol,op_restr
type(psb_dspmat_type), intent(inout) :: ac
type(psb_desc_type), intent(inout) :: desc_ac
integer(psb_ipk_), intent(out) :: info
end subroutine amg_d_newmatch_aggregator_inner_mat_asb
end interface
!!$ interface
!!$ subroutine amg_d_newmatch_unsmth_spmm_bld(a,desc_a,ilaggr,nlaggr,parms,&

@ -12,6 +12,7 @@ subroutine amg_d_newmatch_aggregator_build_tprol(ag,parms,ag_data,&
use psb_base_mod
use amg_base_prec_type
use amg_d_inner_mod
use amg_d_decmatch_mod
#if defined(SERIAL_MPI)
use amg_d_newmatch_aggregator_mod
#else
@ -55,43 +56,6 @@ subroutine amg_d_newmatch_aggregator_build_tprol(ag,parms,ag_data,&
integer(psb_ipk_), save :: idx_mboxp=-1, idx_spmmbld=-1, idx_sweeps_mult=-1
logical, parameter :: dump=.false., do_timings=.true., debug=.false., &
& dump_prol_restr=.false.
interface
function bootCMatch(C,match_alg,n_sweeps,max_nlevels,max_csize,w)&
& bind(c,name='bootCMatch') result(P)
use iso_c_binding
import
implicit none
type(nwm_CSRMatrix) :: C, P
type(nwm_Vector) :: w
integer(c_int) :: match_alg
integer(c_int) :: n_sweeps
integer(c_int) :: max_nlevels
integer(c_int) :: max_csize
end function bootCMatch
end interface
interface
function amg_bootCMatch_if(C,match_alg,n_sweeps,max_nlevels,max_csize,&
& w,isz,ilaggr,valaggr, num_cols) &
& bind(c,name='amg_bootCMatch_if') result(iret)
use iso_c_binding
import
implicit none
type(nwm_CSRMatrix) :: C, P
type(nwm_Vector) :: w
integer(c_int), value :: match_alg
integer(c_int), value :: n_sweeps
integer(c_int), value :: max_nlevels
integer(c_int), value :: max_csize
integer(c_int), value :: isz
integer(c_int) :: num_cols
integer(c_int) :: ilaggr(*)
real(c_double) :: valaggr(*)
integer(c_int) :: iret
end function amg_bootCMatch_if
end interface
name='d_newmatch_tprol'
ictxt = desc_a%get_context()
call psb_info(ictxt,me,np)
@ -121,16 +85,6 @@ subroutine amg_d_newmatch_aggregator_build_tprol(ag,parms,ag_data,&
& amg_aggr_ord_nat_,is_legal_ml_aggr_ord)
call amg_check_def(parms%aggr_thresh,'Aggr_Thresh',dzero,is_legal_d_aggr_thrs)
!write(*,*) 'Build_tprol:',acsr%get_nrows(),acsr%get_ncols()
C%num_rows = acsr%get_nrows()
C%num_cols = acsr%get_ncols()
C%num_nonzeros = acsr%get_nzeros()
C%owns_data = 0
acsr%irp = acsr%irp - 1
acsr%ja = acsr%ja - 1
C%i = c_loc(acsr%irp)
C%j = c_loc(acsr%ja)
C%data = c_loc(acsr%val)
#if !defined(SERIAL_MPI)
match_algorithm = ag%matching_alg
@ -166,8 +120,9 @@ subroutine amg_d_newmatch_aggregator_build_tprol(ag,parms,ag_data,&
if (n_sweeps /= ag%n_sweeps) then
write(0,*) me,' Inconsistent N_SWEEPS ',n_sweeps,ag%n_sweeps
end if
!!$ if (me==0) write(0,*) 'Matching sweeps: ',n_sweeps
n_sweeps = max(1,n_sweeps)
if (debug) write(0,*) me,' Copies, with n_sweeps: ',n_sweeps,max_csize
if (ag%unsmoothed_hierarchy.and.allocated(ag%base_a)) then
call ag%base_a%cp_to(acsr)
@ -276,8 +231,8 @@ subroutine amg_d_newmatch_aggregator_build_tprol(ag,parms,ag_data,&
!
if (debug) write(0,*) me,' Into matchbox_build_prol ',info
if (do_timings) call psb_tic(idx_mboxp)
!!$ call amg_dmatchboxp_build_prol(tmpw,acv(i-1),desc_acv(i-1),ixaggr,nxaggr,tmp_prol,info,&
!!$ & symmetrize=ag%need_symmetrize,reproducible=ag%reproducible_matching)
call amg_ddecmatch_build_prol(tmpw,acv(i-1),desc_acv(i-1),ixaggr,nxaggr,tmp_prol,info,&
& symmetrize=ag%need_symmetrize,reproducible=ag%reproducible_matching)
if (do_timings) call psb_toc(idx_mboxp)
if (debug) write(0,*) me,' Out from matchbox_build_prol ',info
if (psb_errstatus_fatal()) write(0,*)me,trim(name),'Error fatal on exit bld_tprol',info

@ -1,72 +1,76 @@
#include <string.h>
#include <stdio.h>
#include "psb_base_cbind.h"
#include "MatchingAlgorithms.h"
#include "bcm.h"
bcm_CSRMatrix bootCMatch(bcm_CSRMatrix *C, int *match_algorithm, int *n_sweeps, int *max_nlevels, int *max_csize, bcm_Vector *w);
bcm_CSRMatrix bootCMatch(bcm_CSRMatrix *C, int *match_algorithm, int *n_sweeps, int *max_nlevels, int *max_csize, bcm_Vector *w){
bcm_Vector *w_temp;
int info;
//double *w_inp;
//w_inp=bcm_VectorData(w);
bcm_CSRMatrix *P;
bcm_CSRMatrix *Ac;
int ftcoarse=1;
int cr_it=0, cr_relax_type=0;
double cr_relax_weight=0.0;
// Here I am building Ac but I won't use it.
Ac=bcm_CSRMatchingAgg(C, &w, &P, *match_algorithm, *n_sweeps,*max_csize, *max_nlevels, &ftcoarse,
cr_it, cr_relax_type, cr_relax_weight);
//w_inp=bcm_VectorData(w);
bcm_CSRMatrixDestroy(Ac);
return *P;
}
#ifdef __cplusplus
extern "C" {
#endif
int mld_bootCMatch_if(bcm_CSRMatrix *C, int match_algorithm, int n_sweeps,
int max_nlevels, int max_csize, bcm_Vector *w,
int isz, int ilaggr[], double valaggr[], int *num_cols){
bcm_Vector *w_temp;
int info;
//double *w_inp;
//w_inp=bcm_VectorData(w);
psb_i_t dnew_Match_If(psb_i_t nr, psb_i_t irp[], psb_i_t ja[],
psb_d_t val[], psb_d_t diag[],
psb_d_t w[], psb_i_t mate[])
bcm_CSRMatrix *P;
bcm_CSRMatrix *Ac;
int *irp, *ja, nr, nz, nc,i,j;
double *val;
int ftcoarse=1;
int cr_it=0, cr_relax_type=0;
double cr_relax_weight=0.0;
#ifdef __cplusplus
}
#endif
psb_i_t dnew_Match_If(psb_i_t nr, psb_i_t irp[], psb_i_t ja[],
psb_d_t val[], psb_d_t diag[], psb_d_t w[],
psb_i_t mate[])
{
psb_i_t info;
psb_i_t i,j;
psb_i_t ftcoarse=1;
psb_i_t cr_it=0, cr_relax_type=0;
psb_d_t cr_relax_weight=0.0;
// Sanity checks
nr = bcm_CSRMatrixNumRows(C);
nc = bcm_VectorSize(w);
// fprintf(stderr,"Sanity check: %d %d \n",nr,nc);
vector<NODE_T> s;
vector<NODE_T> t;
vector<VAL_T> weights;
vector<NODE_T> mateNode;
NODE_T u,v;
VAL_T weight;
psb_i_t preprocess = atoi(argv[2]);
psb_i_t romaInput = atoi(argv[3]);
VAL_T lambda = atof(argv[4]);
psb_d_t aii, ajj, aij, wii, wjj, tmp1, tmp2, minabs, edgnrm;
psb_i_t nt;
psb_d_t timeDiff;
MatchStat pstat;
// Here I am building Ac but I won't use it.
Ac=bcm_CSRMatchingAgg(C, &w, &P, match_algorithm, n_sweeps, max_csize, max_nlevels, &ftcoarse,
cr_it, cr_relax_type, cr_relax_weight);
irp = bcm_CSRMatrixI(P);
ja = bcm_CSRMatrixJ(P);
val = bcm_CSRMatrixData(P);
nr = bcm_CSRMatrixNumRows(P);
nc = bcm_CSRMatrixNumCols(P);
nz = bcm_CSRMatrixNumNonzeros(P);
// fprintf(stderr,"Sanity check: %d %d \n",nr,nc);
for (i=1; i<nr; i++) {
for (j=irp[i-1]; j<irp[i]; j++) {
v = i-1; // I
u = ja[j-1]; // J
if (v>u) {
// Define Ahat entry
aij = val[j-1];
aii = diag[v];
ajj = diag[u];
wii = w[v];
wjj = w[u];
edgnrm = aii*(wii*wii) + ajj*(wjj*wjj);
if (edgnrm > eps) {
weight = 1.0 - (2*1.0*aij*wii*wjj)/(aii*(wii*wii) + ajj*(wjj*wjj));
} else {
weight = 1e-16;
}
//
s.push_back(u);
t.push_back(v);
weights.push_back(weight);
}
}
}
runRomaWrapper(s,t,weights, nr, mateNode,preprocess,romaInput,lambda ,nt, pstat, timeDiff);
if (isz < nr) return(-1);
if (nz != nr) return(-2);
/* loop here only makes sense when nr==nz */
for (i=0; i< nr; i++) {
for (j=irp[i]; j<irp[i+1]; j++) {
ilaggr[i] = ja[j] + 1;
valaggr[i] = val[j];
}
mate[i] = mateNode[i]+1;
}
*num_cols = nc;
//w_inp=bcm_VectorData(w);
bcm_CSRMatrixDestroy(Ac);
bcm_CSRMatrixDestroy(P);
return(0);
}

Loading…
Cancel
Save