Adjusted CBIND test program.

This reverts commit a17f503486.

Make.inc.in

@@ -75,7 +75,7 @@ CDEFINES=$(AMGCDEFINES)
 AMGFDEFINES=@AMGFDEFINES@ $(PSBFDEFINES)
 FDEFINES=$(AMGFDEFINES)
 
-CXXDEFINES=@AMGCXXDEFINES@
+CXXDEFINES=@AMGCXXDEFINES@ $(PSBCXXDEFINES)
 
 @COMPILERULES@
 

Make_n.inc.in

@@ -49,7 +49,7 @@ PSBLAS_INCLUDES=@PSBLAS_INCLUDES@
 PSBLAS_LIBS=@PSBLAS_LIBS@	
 PSBBASEMODNAME=psb_base_mod
 PSBPRECMODNAME=psb_prec_mod
-PSBMETHDMODNAME=psb_krylov_mod
+PSBMETHDMODNAME=psb_linsolve_mod
 PSBUTILMODNAME=psb_util_mod
 
 

Makefile

@@ -3,9 +3,9 @@ include Make.inc
 
 all:  objs lib
 
-objs: amgp cbnd
+objs: libdir amgp cbnd
 
-lib: libdir objs
+lib:  objs
 	cd amgprec && $(MAKE) lib
 	cd cbind && $(MAKE)  lib
 
@@ -46,6 +46,7 @@ cleanlib:
 
 veryclean: cleanlib
 	(cd amgprec && $(MAKE) veryclean)
+	(cd cbind && $(MAKE) veryclean)
 	(cd samples/simple/fileread && $(MAKE) clean)
 	(cd samples/simple/pdegen && $(MAKE) clean)
 	(cd samples/advanced/fileread && $(MAKE) clean)
@@ -56,3 +57,4 @@ check: all
 
 clean:
 	(cd amgprec && $(MAKE) clean)
+	(cd cbind && $(MAKE) clean)

amgprec/amg_base_prec_type.F90

@@ -323,9 +323,10 @@ module amg_base_prec_type
   !
   ! Legal values for entry: amg_poly_variant_
   !
-  integer(psb_ipk_), parameter :: amg_poly_lottes_      = 0
-  integer(psb_ipk_), parameter :: amg_poly_lottes_beta_ = 1
-  integer(psb_ipk_), parameter :: amg_poly_new_         = 2
+  integer(psb_ipk_), parameter :: amg_cheb_4_     = 0
+  integer(psb_ipk_), parameter :: amg_cheb_4_opt_ = 1
+  integer(psb_ipk_), parameter :: amg_cheb_1_opt_ = 2
+  integer(psb_ipk_), parameter :: amg_poly_dbg_   = 8
   
   integer(psb_ipk_), parameter :: amg_poly_rho_est_power_ = 0
 
@@ -569,12 +570,14 @@ contains
       val = amg_as_
     case('POLY')
       val = amg_poly_
-    case('POLY_LOTTES')
-      val = amg_poly_lottes_
-    case('POLY_LOTTES_BETA')
-      val = amg_poly_lottes_beta_
-    case('POLY_NEW')
-      val = amg_poly_new_
+    case('CHEB_4')
+      val = amg_cheb_4_
+    case('CHEB_4_OPT')
+      val = amg_cheb_4_opt_
+    case('CHEB_1_OPT')
+      val = amg_cheb_1_opt_
+    case('POLY_DBG')
+      val = amg_poly_dbg_
     case('POLY_RHO_EST_POWER')
       val =  amg_poly_rho_est_power_
     case('A_NORMI')

amgprec/amg_c_ainv_solver.F90

@@ -106,7 +106,7 @@ module amg_c_ainv_solver
       Implicit None
 
       ! Arguments
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_c_ainv_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_c_as_smoother.f90

@@ -230,7 +230,7 @@ module amg_c_as_smoother
            & psb_desc_type, psb_c_base_sparse_mat, psb_ipk_,&
            & psb_i_base_vect_type
       implicit none 
-      type(psb_cspmat_type), intent(inout), target     :: a
+      type(psb_cspmat_type), intent(in), target        :: a
       Type(psb_desc_type), Intent(inout)                 :: desc_a 
       class(amg_c_as_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                     :: info

amgprec/amg_c_base_smoother_mod.f90

@@ -237,7 +237,7 @@ module amg_c_base_smoother_mod
            & psb_c_vect_type, psb_c_base_vect_type, psb_spk_, &
            & amg_c_base_smoother_type, psb_ipk_, psb_i_base_vect_type
       ! Arguments
-      type(psb_cspmat_type), intent(inout), target  :: a
+      type(psb_cspmat_type), intent(in), target     :: a
       Type(psb_desc_type), Intent(inout)              :: desc_a 
       class(amg_c_base_smoother_type), intent(inout) :: sm 
       integer(psb_ipk_), intent(out)                  :: info

amgprec/amg_c_base_solver_mod.f90

@@ -170,7 +170,7 @@ module amg_c_base_solver_mod
       Implicit None
       
       ! Arguments
-      type(psb_cspmat_type), intent(inout), target          :: a
+      type(psb_cspmat_type), intent(in), target             :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_c_base_solver_type), intent(inout)          :: sv
       integer(psb_ipk_), intent(out)                        :: info

amgprec/amg_c_diag_solver.f90

@@ -119,7 +119,7 @@ module amg_c_diag_solver
       import :: psb_desc_type, psb_cspmat_type,  psb_c_base_sparse_mat, &
            & psb_c_vect_type, psb_c_base_vect_type, psb_spk_, &
            & amg_c_diag_solver_type, psb_ipk_, psb_i_base_vect_type      
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_c_diag_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                        :: info
@@ -331,7 +331,7 @@ module amg_c_l1_diag_solver
       import :: psb_desc_type, psb_cspmat_type,  psb_c_base_sparse_mat, &
            & psb_c_vect_type, psb_c_base_vect_type, psb_spk_, &
            & amg_c_l1_diag_solver_type, psb_ipk_, psb_i_base_vect_type      
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_c_l1_diag_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                        :: info

amgprec/amg_c_gs_solver.f90

@@ -181,7 +181,7 @@ module amg_c_gs_solver
            & psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_c_gs_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info
@@ -195,7 +195,7 @@ module amg_c_gs_solver
            & psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_c_bwgs_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_c_id_solver.f90

@@ -123,7 +123,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_cspmat_type), intent(inout), target        :: a
+    type(psb_cspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_c_id_solver_type), intent(inout)          :: sv
     integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_c_ilu_solver.f90

@@ -144,7 +144,7 @@ module amg_c_ilu_solver
            & psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_cspmat_type), intent(inout), target       :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_c_ilu_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_c_invk_solver.f90

@@ -94,7 +94,7 @@ module amg_c_invk_solver
       Implicit None
 
       ! Arguments
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_c_invk_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_c_invt_solver.f90

@@ -94,7 +94,7 @@ module amg_c_invt_solver
       Implicit None
 
       ! Arguments
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_c_invt_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_c_jac_smoother.f90

@@ -151,7 +151,7 @@ module amg_c_jac_smoother
       import :: psb_desc_type, amg_c_jac_smoother_type, psb_c_vect_type, psb_spk_, &
            & psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_c_jac_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                      :: info
@@ -274,7 +274,7 @@ module amg_c_jac_smoother
       import :: psb_desc_type, amg_c_l1_jac_smoother_type, psb_c_vect_type, &
            & psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_c_l1_jac_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_c_jac_solver.f90

@@ -143,7 +143,7 @@ module amg_c_jac_solver
            & psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_c_jac_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info
@@ -160,7 +160,7 @@ module amg_c_jac_solver
            & psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_c_l1_jac_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_c_krm_solver.f90

@@ -174,7 +174,7 @@ module amg_c_krm_solver
            & psb_cspmat_type, psb_c_base_sparse_mat, psb_c_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_c_krm_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_c_mumps_solver.F90

@@ -163,7 +163,7 @@ module amg_c_mumps_solver
       Implicit None
 
       ! Arguments
-      type(psb_cspmat_type), intent(inout), target        :: a
+      type(psb_cspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_c_mumps_solver_type), intent(inout)       :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_c_slu_solver.F90

@@ -256,7 +256,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_cspmat_type), intent(inout), target        :: a
+    type(psb_cspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_c_slu_solver_type), intent(inout)         :: sv
     integer, intent(out)                                :: info

amgprec/amg_d_ainv_solver.F90

@@ -106,7 +106,7 @@ module amg_d_ainv_solver
       Implicit None
 
       ! Arguments
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_d_ainv_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_d_as_smoother.f90

@@ -230,7 +230,7 @@ module amg_d_as_smoother
            & psb_desc_type, psb_d_base_sparse_mat, psb_ipk_,&
            & psb_i_base_vect_type
       implicit none 
-      type(psb_dspmat_type), intent(inout), target     :: a
+      type(psb_dspmat_type), intent(in), target        :: a
       Type(psb_desc_type), Intent(inout)                 :: desc_a 
       class(amg_d_as_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                     :: info

amgprec/amg_d_base_smoother_mod.f90

@@ -237,7 +237,7 @@ module amg_d_base_smoother_mod
            & psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, &
            & amg_d_base_smoother_type, psb_ipk_, psb_i_base_vect_type
       ! Arguments
-      type(psb_dspmat_type), intent(inout), target  :: a
+      type(psb_dspmat_type), intent(in), target     :: a
       Type(psb_desc_type), Intent(inout)              :: desc_a 
       class(amg_d_base_smoother_type), intent(inout) :: sm 
       integer(psb_ipk_), intent(out)                  :: info

amgprec/amg_d_base_solver_mod.f90

@@ -170,7 +170,7 @@ module amg_d_base_solver_mod
       Implicit None
       
       ! Arguments
-      type(psb_dspmat_type), intent(inout), target          :: a
+      type(psb_dspmat_type), intent(in), target             :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_d_base_solver_type), intent(inout)          :: sv
       integer(psb_ipk_), intent(out)                        :: info

amgprec/amg_d_diag_solver.f90

@@ -119,7 +119,7 @@ module amg_d_diag_solver
       import :: psb_desc_type, psb_dspmat_type,  psb_d_base_sparse_mat, &
            & psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, &
            & amg_d_diag_solver_type, psb_ipk_, psb_i_base_vect_type      
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_d_diag_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                        :: info
@@ -331,7 +331,7 @@ module amg_d_l1_diag_solver
       import :: psb_desc_type, psb_dspmat_type,  psb_d_base_sparse_mat, &
            & psb_d_vect_type, psb_d_base_vect_type, psb_dpk_, &
            & amg_d_l1_diag_solver_type, psb_ipk_, psb_i_base_vect_type      
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_d_l1_diag_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                        :: info

amgprec/amg_d_gs_solver.f90

@@ -181,7 +181,7 @@ module amg_d_gs_solver
            & psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_d_gs_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info
@@ -195,7 +195,7 @@ module amg_d_gs_solver
            & psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_d_bwgs_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_d_id_solver.f90

@@ -123,7 +123,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_dspmat_type), intent(inout), target        :: a
+    type(psb_dspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_d_id_solver_type), intent(inout)          :: sv
     integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_d_ilu_solver.f90

@@ -144,7 +144,7 @@ module amg_d_ilu_solver
            & psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_dspmat_type), intent(inout), target       :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_d_ilu_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_d_invk_solver.f90

@@ -94,7 +94,7 @@ module amg_d_invk_solver
       Implicit None
 
       ! Arguments
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_d_invk_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_d_invt_solver.f90

@@ -94,7 +94,7 @@ module amg_d_invt_solver
       Implicit None
 
       ! Arguments
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_d_invt_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_d_jac_smoother.f90

@@ -151,7 +151,7 @@ module amg_d_jac_smoother
       import :: psb_desc_type, amg_d_jac_smoother_type, psb_d_vect_type, psb_dpk_, &
            & psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_d_jac_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                      :: info
@@ -274,7 +274,7 @@ module amg_d_jac_smoother
       import :: psb_desc_type, amg_d_l1_jac_smoother_type, psb_d_vect_type, &
            & psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_d_l1_jac_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_d_jac_solver.f90

@@ -143,7 +143,7 @@ module amg_d_jac_solver
            & psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_d_jac_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info
@@ -160,7 +160,7 @@ module amg_d_jac_solver
            & psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_d_l1_jac_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_d_krm_solver.f90

@@ -174,7 +174,7 @@ module amg_d_krm_solver
            & psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_d_krm_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_d_mumps_solver.F90

@@ -163,7 +163,7 @@ module amg_d_mumps_solver
       Implicit None
 
       ! Arguments
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_d_mumps_solver_type), intent(inout)       :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_d_poly_smoother.f90

@@ -140,7 +140,7 @@ module amg_d_poly_smoother
       import :: psb_desc_type, amg_d_poly_smoother_type, psb_d_vect_type, psb_dpk_, &
            & psb_dspmat_type, psb_d_base_sparse_mat, psb_d_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
-      type(psb_dspmat_type), intent(inout), target        :: a
+      type(psb_dspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_d_poly_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                      :: info
@@ -322,7 +322,7 @@ contains
     !
     sm%pdegree      = 1
     sm%rho_ba       = -done
-    sm%variant      = amg_poly_lottes_
+    sm%variant      = amg_cheb_4_
     sm%rho_estimate = amg_poly_rho_est_power_
     sm%rho_estimate_iterations = 20
     if (allocated(sm%sv)) then

amgprec/amg_d_slu_solver.F90

@@ -256,7 +256,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_dspmat_type), intent(inout), target        :: a
+    type(psb_dspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_d_slu_solver_type), intent(inout)         :: sv
     integer, intent(out)                                :: info

amgprec/amg_d_sludist_solver.F90

@@ -259,7 +259,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_dspmat_type), intent(inout), target        :: a
+    type(psb_dspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_d_sludist_solver_type), intent(inout)     :: sv
     integer, intent(out)                                :: info

amgprec/amg_d_umf_solver.F90

@@ -259,7 +259,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_dspmat_type), intent(inout), target        :: a
+    type(psb_dspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_d_umf_solver_type), intent(inout)         :: sv
     integer, intent(out)                                :: info

amgprec/amg_s_ainv_solver.F90

@@ -106,7 +106,7 @@ module amg_s_ainv_solver
       Implicit None
 
       ! Arguments
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_s_ainv_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_s_as_smoother.f90

@@ -230,7 +230,7 @@ module amg_s_as_smoother
            & psb_desc_type, psb_s_base_sparse_mat, psb_ipk_,&
            & psb_i_base_vect_type
       implicit none 
-      type(psb_sspmat_type), intent(inout), target     :: a
+      type(psb_sspmat_type), intent(in), target        :: a
       Type(psb_desc_type), Intent(inout)                 :: desc_a 
       class(amg_s_as_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                     :: info

amgprec/amg_s_base_smoother_mod.f90

@@ -237,7 +237,7 @@ module amg_s_base_smoother_mod
            & psb_s_vect_type, psb_s_base_vect_type, psb_spk_, &
            & amg_s_base_smoother_type, psb_ipk_, psb_i_base_vect_type
       ! Arguments
-      type(psb_sspmat_type), intent(inout), target  :: a
+      type(psb_sspmat_type), intent(in), target     :: a
       Type(psb_desc_type), Intent(inout)              :: desc_a 
       class(amg_s_base_smoother_type), intent(inout) :: sm 
       integer(psb_ipk_), intent(out)                  :: info

amgprec/amg_s_base_solver_mod.f90

@@ -170,7 +170,7 @@ module amg_s_base_solver_mod
       Implicit None
       
       ! Arguments
-      type(psb_sspmat_type), intent(inout), target          :: a
+      type(psb_sspmat_type), intent(in), target             :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_s_base_solver_type), intent(inout)          :: sv
       integer(psb_ipk_), intent(out)                        :: info

amgprec/amg_s_diag_solver.f90

@@ -119,7 +119,7 @@ module amg_s_diag_solver
       import :: psb_desc_type, psb_sspmat_type,  psb_s_base_sparse_mat, &
            & psb_s_vect_type, psb_s_base_vect_type, psb_spk_, &
            & amg_s_diag_solver_type, psb_ipk_, psb_i_base_vect_type      
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_s_diag_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                        :: info
@@ -331,7 +331,7 @@ module amg_s_l1_diag_solver
       import :: psb_desc_type, psb_sspmat_type,  psb_s_base_sparse_mat, &
            & psb_s_vect_type, psb_s_base_vect_type, psb_spk_, &
            & amg_s_l1_diag_solver_type, psb_ipk_, psb_i_base_vect_type      
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_s_l1_diag_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                        :: info

amgprec/amg_s_gs_solver.f90

@@ -181,7 +181,7 @@ module amg_s_gs_solver
            & psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_s_gs_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info
@@ -195,7 +195,7 @@ module amg_s_gs_solver
            & psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_s_bwgs_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_s_id_solver.f90

@@ -123,7 +123,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_sspmat_type), intent(inout), target        :: a
+    type(psb_sspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_s_id_solver_type), intent(inout)          :: sv
     integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_s_ilu_solver.f90

@@ -144,7 +144,7 @@ module amg_s_ilu_solver
            & psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_sspmat_type), intent(inout), target       :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_s_ilu_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_s_invk_solver.f90

@@ -94,7 +94,7 @@ module amg_s_invk_solver
       Implicit None
 
       ! Arguments
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_s_invk_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_s_invt_solver.f90

@@ -94,7 +94,7 @@ module amg_s_invt_solver
       Implicit None
 
       ! Arguments
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_s_invt_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_s_jac_smoother.f90

@@ -151,7 +151,7 @@ module amg_s_jac_smoother
       import :: psb_desc_type, amg_s_jac_smoother_type, psb_s_vect_type, psb_spk_, &
            & psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_s_jac_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                      :: info
@@ -274,7 +274,7 @@ module amg_s_jac_smoother
       import :: psb_desc_type, amg_s_l1_jac_smoother_type, psb_s_vect_type, &
            & psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_s_l1_jac_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_s_jac_solver.f90

@@ -143,7 +143,7 @@ module amg_s_jac_solver
            & psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_s_jac_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info
@@ -160,7 +160,7 @@ module amg_s_jac_solver
            & psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_s_l1_jac_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_s_krm_solver.f90

@@ -174,7 +174,7 @@ module amg_s_krm_solver
            & psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_s_krm_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_s_mumps_solver.F90

@@ -163,7 +163,7 @@ module amg_s_mumps_solver
       Implicit None
 
       ! Arguments
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_s_mumps_solver_type), intent(inout)       :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_s_poly_smoother.f90

@@ -140,7 +140,7 @@ module amg_s_poly_smoother
       import :: psb_desc_type, amg_s_poly_smoother_type, psb_s_vect_type, psb_spk_, &
            & psb_sspmat_type, psb_s_base_sparse_mat, psb_s_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
-      type(psb_sspmat_type), intent(inout), target        :: a
+      type(psb_sspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_s_poly_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                      :: info
@@ -322,7 +322,7 @@ contains
     !
     sm%pdegree      = 1
     sm%rho_ba       = -sone
-    sm%variant      = amg_poly_lottes_
+    sm%variant      = amg_cheb_4_
     sm%rho_estimate = amg_poly_rho_est_power_
     sm%rho_estimate_iterations = 20
     if (allocated(sm%sv)) then

amgprec/amg_s_slu_solver.F90

@@ -256,7 +256,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_sspmat_type), intent(inout), target        :: a
+    type(psb_sspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_s_slu_solver_type), intent(inout)         :: sv
     integer, intent(out)                                :: info

amgprec/amg_z_ainv_solver.F90

@@ -106,7 +106,7 @@ module amg_z_ainv_solver
       Implicit None
 
       ! Arguments
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_z_ainv_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_z_as_smoother.f90

@@ -230,7 +230,7 @@ module amg_z_as_smoother
            & psb_desc_type, psb_z_base_sparse_mat, psb_ipk_,&
            & psb_i_base_vect_type
       implicit none 
-      type(psb_zspmat_type), intent(inout), target     :: a
+      type(psb_zspmat_type), intent(in), target        :: a
       Type(psb_desc_type), Intent(inout)                 :: desc_a 
       class(amg_z_as_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                     :: info

amgprec/amg_z_base_smoother_mod.f90

@@ -237,7 +237,7 @@ module amg_z_base_smoother_mod
            & psb_z_vect_type, psb_z_base_vect_type, psb_dpk_, &
            & amg_z_base_smoother_type, psb_ipk_, psb_i_base_vect_type
       ! Arguments
-      type(psb_zspmat_type), intent(inout), target  :: a
+      type(psb_zspmat_type), intent(in), target     :: a
       Type(psb_desc_type), Intent(inout)              :: desc_a 
       class(amg_z_base_smoother_type), intent(inout) :: sm 
       integer(psb_ipk_), intent(out)                  :: info

amgprec/amg_z_base_solver_mod.f90

@@ -170,7 +170,7 @@ module amg_z_base_solver_mod
       Implicit None
       
       ! Arguments
-      type(psb_zspmat_type), intent(inout), target          :: a
+      type(psb_zspmat_type), intent(in), target             :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_z_base_solver_type), intent(inout)          :: sv
       integer(psb_ipk_), intent(out)                        :: info

amgprec/amg_z_diag_solver.f90

@@ -119,7 +119,7 @@ module amg_z_diag_solver
       import :: psb_desc_type, psb_zspmat_type,  psb_z_base_sparse_mat, &
            & psb_z_vect_type, psb_z_base_vect_type, psb_dpk_, &
            & amg_z_diag_solver_type, psb_ipk_, psb_i_base_vect_type      
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_z_diag_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                        :: info
@@ -331,7 +331,7 @@ module amg_z_l1_diag_solver
       import :: psb_desc_type, psb_zspmat_type,  psb_z_base_sparse_mat, &
            & psb_z_vect_type, psb_z_base_vect_type, psb_dpk_, &
            & amg_z_l1_diag_solver_type, psb_ipk_, psb_i_base_vect_type      
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                    :: desc_a 
       class(amg_z_l1_diag_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                        :: info

amgprec/amg_z_gs_solver.f90

@@ -181,7 +181,7 @@ module amg_z_gs_solver
            & psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_z_gs_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info
@@ -195,7 +195,7 @@ module amg_z_gs_solver
            & psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_z_bwgs_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_z_id_solver.f90

@@ -123,7 +123,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_zspmat_type), intent(inout), target        :: a
+    type(psb_zspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_z_id_solver_type), intent(inout)          :: sv
     integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_z_ilu_solver.f90

@@ -144,7 +144,7 @@ module amg_z_ilu_solver
            & psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_zspmat_type), intent(inout), target       :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_z_ilu_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_z_invk_solver.f90

@@ -94,7 +94,7 @@ module amg_z_invk_solver
       Implicit None
 
       ! Arguments
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_z_invk_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_z_invt_solver.f90

@@ -94,7 +94,7 @@ module amg_z_invt_solver
       Implicit None
 
       ! Arguments
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_z_invt_solver_type), intent(inout)        :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_z_jac_smoother.f90

@@ -151,7 +151,7 @@ module amg_z_jac_smoother
       import :: psb_desc_type, amg_z_jac_smoother_type, psb_z_vect_type, psb_dpk_, &
            & psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_z_jac_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                      :: info
@@ -274,7 +274,7 @@ module amg_z_jac_smoother
       import :: psb_desc_type, amg_z_l1_jac_smoother_type, psb_z_vect_type, &
            & psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_z_l1_jac_smoother_type), intent(inout)       :: sm
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_z_jac_solver.f90

@@ -143,7 +143,7 @@ module amg_z_jac_solver
            & psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_z_jac_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info
@@ -160,7 +160,7 @@ module amg_z_jac_solver
            & psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none 
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_z_l1_jac_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_z_krm_solver.f90

@@ -174,7 +174,7 @@ module amg_z_krm_solver
            & psb_zspmat_type, psb_z_base_sparse_mat, psb_z_base_vect_type,&
            & psb_ipk_, psb_i_base_vect_type
       implicit none
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a
       class(amg_z_krm_solver_type), intent(inout)         :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_z_mumps_solver.F90

@@ -163,7 +163,7 @@ module amg_z_mumps_solver
       Implicit None
 
       ! Arguments
-      type(psb_zspmat_type), intent(inout), target        :: a
+      type(psb_zspmat_type), intent(in), target           :: a
       Type(psb_desc_type), Intent(inout)                  :: desc_a 
       class(amg_z_mumps_solver_type), intent(inout)       :: sv
       integer(psb_ipk_), intent(out)                      :: info

amgprec/amg_z_slu_solver.F90

@@ -256,7 +256,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_zspmat_type), intent(inout), target        :: a
+    type(psb_zspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_z_slu_solver_type), intent(inout)         :: sv
     integer, intent(out)                                :: info

amgprec/amg_z_sludist_solver.F90

@@ -259,7 +259,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_zspmat_type), intent(inout), target        :: a
+    type(psb_zspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_z_sludist_solver_type), intent(inout)     :: sv
     integer, intent(out)                                :: info

amgprec/amg_z_umf_solver.F90

@@ -259,7 +259,7 @@ contains
     Implicit None
 
     ! Arguments
-    type(psb_zspmat_type), intent(inout), target        :: a
+    type(psb_zspmat_type), intent(in), target           :: a
     Type(psb_desc_type), Intent(inout)                  :: desc_a 
     class(amg_z_umf_solver_type), intent(inout)         :: sv
     integer, intent(out)                                :: info

amgprec/impl/aggregator/MatchBoxPC.cpp

@@ -67,14 +67,13 @@ void dMatchBoxPC(MilanLongInt NLVer, MilanLongInt NLEdge,
 #endif
 
 
-#define TIME_TRACKER
-    #ifdef TIME_TRACKER
-        double tmr = MPI_Wtime();
-    #endif
+#undef TIME_TRACKER
+#ifdef TIME_TRACKER
+  double tmr = MPI_Wtime();
+#endif
 
-// Rimosso per tornare al vecchio matching #define OMP
-#ifdef OMP
-				fprintf(stderr,"Warning: using buggy OpenMP matching!\n");
+#if defined(OPENMP)
+	//fprintf(stderr,"Warning: using buggy OpenMP matching!\n");
         dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(NLVer, NLEdge,
 							   verLocPtr, verLocInd, edgeLocWeight,
 							   verDistance,  Mate,
@@ -93,11 +92,11 @@ void dMatchBoxPC(MilanLongInt NLVer, MilanLongInt NLEdge,
 #endif
 
 
-  #ifdef TIME_TRACKER
-    tmr = MPI_Wtime() - tmr;
-    fprintf(stderr, "Elaboration time: %f for %ld nodes\n", tmr, NLVer);
-  #endif
-
+#ifdef TIME_TRACKER
+	tmr = MPI_Wtime() - tmr;
+	fprintf(stderr, "Elaboration time: %f for %ld nodes\n", tmr, NLVer);
+#endif
+	
 #endif
 }
 
@@ -115,13 +114,24 @@ void sMatchBoxPC(MilanLongInt NLVer, MilanLongInt NLEdge,
   fprintf(stderr,"MatchBoxPC: rank %d nlver %ld nledge %ld [ %ld %ld ]\n",
 	  myRank,NLVer, NLEdge,verDistance[0],verDistance[1]);
 #endif
+#if defined(OPENMP)
+	//fprintf(stderr,"Warning: using buggy OpenMP matching!\n");
+  salgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(NLVer, NLEdge,
+							      verLocPtr, verLocInd, edgeLocWeight,
+							      verDistance,  Mate,
+							      myRank, numProcs, C_comm,
+							      msgIndSent, msgActualSent, msgPercent,
+							      ph0_time, ph1_time, ph2_time,
+							      ph1_card, ph2_card );
+#else
   salgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateC(NLVer, NLEdge,
-							   verLocPtr, verLocInd, edgeLocWeight,
-							   verDistance,  Mate,
-							   myRank, numProcs, C_comm,
-							   msgIndSent, msgActualSent, msgPercent,
-							   ph0_time, ph1_time, ph2_time,
-							   ph1_card, ph2_card );
+							    verLocPtr, verLocInd, edgeLocWeight,
+							    verDistance,  Mate,
+							    myRank, numProcs, C_comm,
+							    msgIndSent, msgActualSent, msgPercent,
+							    ph0_time, ph1_time, ph2_time,
+							    ph1_card, ph2_card );
+#endif
 #endif
 }
 

amgprec/impl/aggregator/MatchBoxPC.h

@@ -59,7 +59,7 @@
 #include <assert.h>
 #include <map>
 #include <vector>
-#ifdef OMP
+#ifdef OPENMP
 // OpenMP is included and used if and only if the OpenMP version of the matching
 // is required
 #include "omp.h"
@@ -82,6 +82,8 @@ const int BundleTag = 9;  // Predefined tag
 
 static vector<MilanLongInt> DEFAULT_VECTOR;
 
+#if !defined(SERIAL_MPI)
+
 // MPI type map
 template <typename T>
 MPI_Datatype TypeMap();
@@ -93,6 +95,7 @@ template <>
 inline MPI_Datatype TypeMap<double>() { return MPI_DOUBLE; }
 template <>
 inline MPI_Datatype TypeMap<float>() { return MPI_FLOAT; }
+#endif
 
 #ifdef __cplusplus
 extern "C"
@@ -178,264 +181,416 @@ extern "C"
 #define MilanRealMin MINUS_INFINITY
 #endif
 
-#ifdef OMP
 /* These functions are only used in the experimental OMP implementation, if that
 is disabled there is no reason to actually compile or reference them. */
 
-    // Function of find the owner of a ghost vertex using binary search:
-    MilanInt findOwnerOfGhost(MilanLongInt vtxIndex, MilanLongInt *mVerDistance,
-                              MilanInt myRank, MilanInt numProcs);
+  // Function of find the owner of a ghost vertex using binary search:
+  MilanInt findOwnerOfGhost(MilanLongInt vtxIndex, MilanLongInt *mVerDistance,
+			    MilanInt myRank, MilanInt numProcs);
+  
+  MilanLongInt firstComputeCandidateMateD(MilanLongInt adj1,
+					 MilanLongInt adj2,
+					 MilanLongInt *verLocInd,
+					 MilanReal *edgeLocWeight);
 
-    MilanLongInt firstComputeCandidateMate(MilanLongInt adj1,
-                                           MilanLongInt adj2,
-                                           MilanLongInt *verLocInd,
-                                           MilanReal *edgeLocWeight);
+  
+  void queuesTransfer(vector<MilanLongInt> &U,
+		      vector<MilanLongInt> &privateU,
+		      vector<MilanLongInt> &QLocalVtx,
+		      vector<MilanLongInt> &QGhostVtx,
+		      vector<MilanLongInt> &QMsgType,
+		      vector<MilanInt> &QOwner,
+		      vector<MilanLongInt> &privateQLocalVtx,
+		      vector<MilanLongInt> &privateQGhostVtx,
+		      vector<MilanLongInt> &privateQMsgType,
+		      vector<MilanInt> &privateQOwner);
+  
+  bool isAlreadyMatched(MilanLongInt node,
+			MilanLongInt StartIndex,
+			MilanLongInt EndIndex,
+			vector<MilanLongInt> &GMate,
+			MilanLongInt *Mate,
+			map<MilanLongInt, MilanLongInt> &Ghost2LocalMap);
+  
+  MilanLongInt computeCandidateMateD(MilanLongInt adj1,
+				    MilanLongInt adj2,
+				    MilanReal *edgeLocWeight,
+				    MilanLongInt k,
+				    MilanLongInt *verLocInd,
+				    MilanLongInt StartIndex,
+				    MilanLongInt EndIndex,
+				    vector<MilanLongInt> &GMate,
+				    MilanLongInt *Mate,
+				    map<MilanLongInt, MilanLongInt> &Ghost2LocalMap);
+  
+  void initialize(MilanLongInt NLVer, MilanLongInt NLEdge,
+		  MilanLongInt StartIndex, MilanLongInt EndIndex,
+		  MilanLongInt *numGhostEdgesPtr,
+		  MilanLongInt *numGhostVerticesPtr,
+		  MilanLongInt *S,
+		  MilanLongInt *verLocInd,
+		  MilanLongInt *verLocPtr,
+		  map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
+		  vector<MilanLongInt> &Counter,
+		  vector<MilanLongInt> &verGhostPtr,
+		  vector<MilanLongInt> &verGhostInd,
+		  vector<MilanLongInt> &tempCounter,
+		  vector<MilanLongInt> &GMate,
+		  vector<MilanLongInt> &Message,
+		  vector<MilanLongInt> &QLocalVtx,
+		  vector<MilanLongInt> &QGhostVtx,
+		  vector<MilanLongInt> &QMsgType,
+		  vector<MilanInt> &QOwner,
+		  MilanLongInt *&candidateMate,
+		  vector<MilanLongInt> &U,
+		  vector<MilanLongInt> &privateU,
+		  vector<MilanLongInt> &privateQLocalVtx,
+		  vector<MilanLongInt> &privateQGhostVtx,
+		  vector<MilanLongInt> &privateQMsgType,
+		  vector<MilanInt> &privateQOwner);
+  
+  void clean(MilanLongInt NLVer,
+	     MilanInt myRank,
+	     MilanLongInt MessageIndex,
+	     vector<MPI_Request> &SRequest,
+	     vector<MPI_Status> &SStatus,
+	     MilanInt BufferSize,
+	     MilanLongInt *Buffer,
+	     MilanLongInt msgActual,
+	     MilanLongInt *msgActualSent,
+	     MilanLongInt msgInd,
+	     MilanLongInt *msgIndSent,
+	     MilanLongInt NumMessagesBundled,
+	     MilanReal *msgPercent);
+  
+  void PARALLEL_COMPUTE_CANDIDATE_MATE_BD(MilanLongInt NLVer,
+					 MilanLongInt *verLocPtr,
+					 MilanLongInt *verLocInd,
+					 MilanInt myRank,
+					 MilanReal *edgeLocWeight,
+					 MilanLongInt *candidateMate);
+  
+  void PARALLEL_PROCESS_EXPOSED_VERTEX_BD(MilanLongInt NLVer,
+					 MilanLongInt *candidateMate,
+					 MilanLongInt *verLocInd,
+					 MilanLongInt *verLocPtr,
+					 MilanLongInt StartIndex,
+					 MilanLongInt EndIndex,
+					 MilanLongInt *Mate,
+					 vector<MilanLongInt> &GMate,
+					 map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
+					 MilanReal *edgeLocWeight,
+					 MilanLongInt *myCardPtr,
+					 MilanLongInt *msgIndPtr,
+					 MilanLongInt *NumMessagesBundledPtr,
+					 MilanLongInt *SPtr,
+					 MilanLongInt *verDistance,
+					 MilanLongInt *PCounter,
+					 vector<MilanLongInt> &Counter,
+					 MilanInt myRank,
+					 MilanInt numProcs,
+					 vector<MilanLongInt> &U,
+					 vector<MilanLongInt> &privateU,
+					 vector<MilanLongInt> &QLocalVtx,
+					 vector<MilanLongInt> &QGhostVtx,
+					 vector<MilanLongInt> &QMsgType,
+					 vector<MilanInt> &QOwner,
+					 vector<MilanLongInt> &privateQLocalVtx,
+					 vector<MilanLongInt> &privateQGhostVtx,
+					 vector<MilanLongInt> &privateQMsgType,
+					 vector<MilanInt> &privateQOwner);
+  
+  void PROCESS_CROSS_EDGE(MilanLongInt *edge,
+			  MilanLongInt *SPtr);
+  
+  void processMatchedVerticesD(
+			      MilanLongInt NLVer,
+			      vector<MilanLongInt> &UChunkBeingProcessed,
+			      vector<MilanLongInt> &U,
+			      vector<MilanLongInt> &privateU,
+			      MilanLongInt StartIndex,
+			      MilanLongInt EndIndex,
+			      MilanLongInt *myCardPtr,
+			      MilanLongInt *msgIndPtr,
+			      MilanLongInt *NumMessagesBundledPtr,
+			      MilanLongInt *SPtr,
+			      MilanLongInt *verLocPtr,
+			      MilanLongInt *verLocInd,
+			      MilanLongInt *verDistance,
+			      MilanLongInt *PCounter,
+			      vector<MilanLongInt> &Counter,
+			      MilanInt myRank,
+			      MilanInt numProcs,
+			      MilanLongInt *candidateMate,
+			      vector<MilanLongInt> &GMate,
+			      MilanLongInt *Mate,
+			      map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
+			      MilanReal *edgeLocWeight,
+			      vector<MilanLongInt> &QLocalVtx,
+			      vector<MilanLongInt> &QGhostVtx,
+			      vector<MilanLongInt> &QMsgType,
+			      vector<MilanInt> &QOwner,
+			      vector<MilanLongInt> &privateQLocalVtx,
+			      vector<MilanLongInt> &privateQGhostVtx,
+			      vector<MilanLongInt> &privateQMsgType,
+			      vector<MilanInt> &privateQOwner);
+  
+  void processMatchedVerticesAndSendMessagesD(
+					     MilanLongInt NLVer,
+					     vector<MilanLongInt> &UChunkBeingProcessed,
+					     vector<MilanLongInt> &U,
+					     vector<MilanLongInt> &privateU,
+					     MilanLongInt StartIndex,
+					     MilanLongInt EndIndex,
+					     MilanLongInt *myCardPtr,
+					     MilanLongInt *msgIndPtr,
+					     MilanLongInt *NumMessagesBundledPtr,
+					     MilanLongInt *SPtr,
+					     MilanLongInt *verLocPtr,
+					     MilanLongInt *verLocInd,
+					     MilanLongInt *verDistance,
+					     MilanLongInt *PCounter,
+					     vector<MilanLongInt> &Counter,
+					     MilanInt myRank,
+					     MilanInt numProcs,
+					     MilanLongInt *candidateMate,
+					     vector<MilanLongInt> &GMate,
+					     MilanLongInt *Mate,
+					     map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
+					     MilanReal *edgeLocWeight,
+					     vector<MilanLongInt> &QLocalVtx,
+					     vector<MilanLongInt> &QGhostVtx,
+					     vector<MilanLongInt> &QMsgType,
+					     vector<MilanInt> &QOwner,
+					     vector<MilanLongInt> &privateQLocalVtx,
+					     vector<MilanLongInt> &privateQGhostVtx,
+					     vector<MilanLongInt> &privateQMsgType,
+					     vector<MilanInt> &privateQOwner,
+					     MPI_Comm comm,
+					     MilanLongInt *msgActual,
+					     vector<MilanLongInt> &Message);
+  
+  void sendBundledMessages(MilanLongInt *numGhostEdgesPtr,
+			   MilanInt *BufferSizePtr,
+			   MilanLongInt *Buffer,
+			   vector<MilanLongInt> &PCumulative,
+			   vector<MilanLongInt> &PMessageBundle,
+			   vector<MilanLongInt> &PSizeInfoMessages,
+			   MilanLongInt *PCounter,
+			   MilanLongInt NumMessagesBundled,
+			   MilanLongInt *msgActualPtr,
+			   MilanLongInt *MessageIndexPtr,
+			   MilanInt numProcs,
+			   MilanInt myRank,
+			   MPI_Comm comm,
+			   vector<MilanLongInt> &QLocalVtx,
+			   vector<MilanLongInt> &QGhostVtx,
+			   vector<MilanLongInt> &QMsgType,
+			   vector<MilanInt> &QOwner,
+			   vector<MPI_Request> &SRequest,
+			   vector<MPI_Status> &SStatus);
+  
+  void processMessagesD(
+		       MilanLongInt NLVer,
+		       MilanLongInt *Mate,
+		       MilanLongInt *candidateMate,
+		       map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
+		       vector<MilanLongInt> &GMate,
+		       vector<MilanLongInt> &Counter,
+		       MilanLongInt StartIndex,
+		       MilanLongInt EndIndex,
+		       MilanLongInt *myCardPtr,
+		       MilanLongInt *msgIndPtr,
+		       MilanLongInt *msgActualPtr,
+		       MilanReal *edgeLocWeight,
+		       MilanLongInt *verDistance,
+		       MilanLongInt *verLocPtr,
+		       MilanLongInt k,
+		       MilanLongInt *verLocInd,
+		       MilanInt numProcs,
+		       MilanInt myRank,
+		       MPI_Comm comm,
+		       vector<MilanLongInt> &Message,
+		       MilanLongInt numGhostEdges,
+		       MilanLongInt u,
+		       MilanLongInt v,
+		       MilanLongInt *SPtr,
+		       vector<MilanLongInt> &U);
+  
+  void extractUChunk(
+		     vector<MilanLongInt> &UChunkBeingProcessed,
+		     vector<MilanLongInt> &U,
+		     vector<MilanLongInt> &privateU);
 
-    void queuesTransfer(vector<MilanLongInt> &U,
-                        vector<MilanLongInt> &privateU,
-                        vector<MilanLongInt> &QLocalVtx,
-                        vector<MilanLongInt> &QGhostVtx,
-                        vector<MilanLongInt> &QMsgType,
-                        vector<MilanInt> &QOwner,
-                        vector<MilanLongInt> &privateQLocalVtx,
-                        vector<MilanLongInt> &privateQGhostVtx,
-                        vector<MilanLongInt> &privateQMsgType,
-                        vector<MilanInt> &privateQOwner);
+  MilanLongInt firstComputeCandidateMateS(MilanLongInt adj1,
+					 MilanLongInt adj2,
+					 MilanLongInt *verLocInd,
+					  MilanFloat *edgeLocWeight);
 
-    bool isAlreadyMatched(MilanLongInt node,
-                          MilanLongInt StartIndex,
-                          MilanLongInt EndIndex,
-                          vector<MilanLongInt> &GMate,
-                          MilanLongInt *Mate,
-                          map<MilanLongInt, MilanLongInt> &Ghost2LocalMap);
+  MilanLongInt computeCandidateMateS(MilanLongInt adj1,
+				    MilanLongInt adj2,
+				    MilanFloat *edgeLocWeight,
+				    MilanLongInt k,
+				    MilanLongInt *verLocInd,
+				    MilanLongInt StartIndex,
+				    MilanLongInt EndIndex,
+				    vector<MilanLongInt> &GMate,
+				    MilanLongInt *Mate,
+				    map<MilanLongInt, MilanLongInt> &Ghost2LocalMap);
+  
+  void PARALLEL_COMPUTE_CANDIDATE_MATE_BS(MilanLongInt NLVer,
+					 MilanLongInt *verLocPtr,
+					 MilanLongInt *verLocInd,
+					 MilanInt myRank,
+					 MilanFloat *edgeLocWeight,
+					 MilanLongInt *candidateMate);
+  
+  void PARALLEL_PROCESS_EXPOSED_VERTEX_BS(MilanLongInt NLVer,
+					 MilanLongInt *candidateMate,
+					 MilanLongInt *verLocInd,
+					 MilanLongInt *verLocPtr,
+					 MilanLongInt StartIndex,
+					 MilanLongInt EndIndex,
+					 MilanLongInt *Mate,
+					 vector<MilanLongInt> &GMate,
+					 map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
+					 MilanFloat *edgeLocWeight,
+					 MilanLongInt *myCardPtr,
+					 MilanLongInt *msgIndPtr,
+					 MilanLongInt *NumMessagesBundledPtr,
+					 MilanLongInt *SPtr,
+					 MilanLongInt *verDistance,
+					 MilanLongInt *PCounter,
+					 vector<MilanLongInt> &Counter,
+					 MilanInt myRank,
+					 MilanInt numProcs,
+					 vector<MilanLongInt> &U,
+					 vector<MilanLongInt> &privateU,
+					 vector<MilanLongInt> &QLocalVtx,
+					 vector<MilanLongInt> &QGhostVtx,
+					 vector<MilanLongInt> &QMsgType,
+					 vector<MilanInt> &QOwner,
+					 vector<MilanLongInt> &privateQLocalVtx,
+					 vector<MilanLongInt> &privateQGhostVtx,
+					 vector<MilanLongInt> &privateQMsgType,
+					 vector<MilanInt> &privateQOwner);
+  void processMatchedVerticesS(
+			      MilanLongInt NLVer,
+			      vector<MilanLongInt> &UChunkBeingProcessed,
+			      vector<MilanLongInt> &U,
+			      vector<MilanLongInt> &privateU,
+			      MilanLongInt StartIndex,
+			      MilanLongInt EndIndex,
+			      MilanLongInt *myCardPtr,
+			      MilanLongInt *msgIndPtr,
+			      MilanLongInt *NumMessagesBundledPtr,
+			      MilanLongInt *SPtr,
+			      MilanLongInt *verLocPtr,
+			      MilanLongInt *verLocInd,
+			      MilanLongInt *verDistance,
+			      MilanLongInt *PCounter,
+			      vector<MilanLongInt> &Counter,
+			      MilanInt myRank,
+			      MilanInt numProcs,
+			      MilanLongInt *candidateMate,
+			      vector<MilanLongInt> &GMate,
+			      MilanLongInt *Mate,
+			      map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
+			      MilanFloat *edgeLocWeight,
+			      vector<MilanLongInt> &QLocalVtx,
+			      vector<MilanLongInt> &QGhostVtx,
+			      vector<MilanLongInt> &QMsgType,
+			      vector<MilanInt> &QOwner,
+			      vector<MilanLongInt> &privateQLocalVtx,
+			      vector<MilanLongInt> &privateQGhostVtx,
+			      vector<MilanLongInt> &privateQMsgType,
+			      vector<MilanInt> &privateQOwner);
+  
+  void processMatchedVerticesAndSendMessagesS(
+					     MilanLongInt NLVer,
+					     vector<MilanLongInt> &UChunkBeingProcessed,
+					     vector<MilanLongInt> &U,
+					     vector<MilanLongInt> &privateU,
+					     MilanLongInt StartIndex,
+					     MilanLongInt EndIndex,
+					     MilanLongInt *myCardPtr,
+					     MilanLongInt *msgIndPtr,
+					     MilanLongInt *NumMessagesBundledPtr,
+					     MilanLongInt *SPtr,
+					     MilanLongInt *verLocPtr,
+					     MilanLongInt *verLocInd,
+					     MilanLongInt *verDistance,
+					     MilanLongInt *PCounter,
+					     vector<MilanLongInt> &Counter,
+					     MilanInt myRank,
+					     MilanInt numProcs,
+					     MilanLongInt *candidateMate,
+					     vector<MilanLongInt> &GMate,
+					     MilanLongInt *Mate,
+					     map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
+					     MilanFloat *edgeLocWeight,
+					     vector<MilanLongInt> &QLocalVtx,
+					     vector<MilanLongInt> &QGhostVtx,
+					     vector<MilanLongInt> &QMsgType,
+					     vector<MilanInt> &QOwner,
+					     vector<MilanLongInt> &privateQLocalVtx,
+					     vector<MilanLongInt> &privateQGhostVtx,
+					     vector<MilanLongInt> &privateQMsgType,
+					     vector<MilanInt> &privateQOwner,
+					     MPI_Comm comm,
+					     MilanLongInt *msgActual,
+					     vector<MilanLongInt> &Message);
+  
+  void processMessagesS(
+		       MilanLongInt NLVer,
+		       MilanLongInt *Mate,
+		       MilanLongInt *candidateMate,
+		       map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
+		       vector<MilanLongInt> &GMate,
+		       vector<MilanLongInt> &Counter,
+		       MilanLongInt StartIndex,
+		       MilanLongInt EndIndex,
+		       MilanLongInt *myCardPtr,
+		       MilanLongInt *msgIndPtr,
+		       MilanLongInt *msgActualPtr,
+		       MilanFloat *edgeLocWeight,
+		       MilanLongInt *verDistance,
+		       MilanLongInt *verLocPtr,
+		       MilanLongInt k,
+		       MilanLongInt *verLocInd,
+		       MilanInt numProcs,
+		       MilanInt myRank,
+		       MPI_Comm comm,
+		       vector<MilanLongInt> &Message,
+		       MilanLongInt numGhostEdges,
+		       MilanLongInt u,
+		       MilanLongInt v,
+		       MilanLongInt *SPtr,
+		       vector<MilanLongInt> &U);
+  
+  
+  void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
+			 MilanLongInt NLVer, MilanLongInt NLEdge,
+			 MilanLongInt *verLocPtr, MilanLongInt *verLocInd, MilanReal *edgeLocWeight,
+			 MilanLongInt *verDistance,
+			 MilanLongInt *Mate,
+			 MilanInt myRank, MilanInt numProcs, MPI_Comm comm,
+			 MilanLongInt *msgIndSent, MilanLongInt *msgActualSent, MilanReal *msgPercent,
+			 MilanReal *ph0_time, MilanReal *ph1_time, MilanReal *ph2_time,
+			 MilanLongInt *ph1_card, MilanLongInt *ph2_card);
 
-    MilanLongInt computeCandidateMate(MilanLongInt adj1,
-                                      MilanLongInt adj2,
-                                      MilanReal *edgeLocWeight,
-                                      MilanLongInt k,
-                                      MilanLongInt *verLocInd,
-                                      MilanLongInt StartIndex,
-                                      MilanLongInt EndIndex,
-                                      vector<MilanLongInt> &GMate,
-                                      MilanLongInt *Mate,
-                                      map<MilanLongInt, MilanLongInt> &Ghost2LocalMap);
-
-    void initialize(MilanLongInt NLVer, MilanLongInt NLEdge,
-                    MilanLongInt StartIndex, MilanLongInt EndIndex,
-                    MilanLongInt *numGhostEdgesPtr,
-                    MilanLongInt *numGhostVerticesPtr,
-                    MilanLongInt *S,
-                    MilanLongInt *verLocInd,
-                    MilanLongInt *verLocPtr,
-                    map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
-                    vector<MilanLongInt> &Counter,
-                    vector<MilanLongInt> &verGhostPtr,
-                    vector<MilanLongInt> &verGhostInd,
-                    vector<MilanLongInt> &tempCounter,
-                    vector<MilanLongInt> &GMate,
-                    vector<MilanLongInt> &Message,
-                    vector<MilanLongInt> &QLocalVtx,
-                    vector<MilanLongInt> &QGhostVtx,
-                    vector<MilanLongInt> &QMsgType,
-                    vector<MilanInt> &QOwner,
-                    MilanLongInt *&candidateMate,
-                    vector<MilanLongInt> &U,
-                    vector<MilanLongInt> &privateU,
-                    vector<MilanLongInt> &privateQLocalVtx,
-                    vector<MilanLongInt> &privateQGhostVtx,
-                    vector<MilanLongInt> &privateQMsgType,
-                    vector<MilanInt> &privateQOwner);
-
-    void clean(MilanLongInt NLVer,
-               MilanInt myRank,
-               MilanLongInt MessageIndex,
-               vector<MPI_Request> &SRequest,
-               vector<MPI_Status> &SStatus,
-               MilanInt BufferSize,
-               MilanLongInt *Buffer,
-               MilanLongInt msgActual,
-               MilanLongInt *msgActualSent,
-               MilanLongInt msgInd,
-               MilanLongInt *msgIndSent,
-               MilanLongInt NumMessagesBundled,
-               MilanReal *msgPercent);
-
-    void PARALLEL_COMPUTE_CANDIDATE_MATE_B(MilanLongInt NLVer,
-                                           MilanLongInt *verLocPtr,
-                                           MilanLongInt *verLocInd,
-                                           MilanInt myRank,
-                                           MilanReal *edgeLocWeight,
-                                           MilanLongInt *candidateMate);
-
-    void PARALLEL_PROCESS_EXPOSED_VERTEX_B(MilanLongInt NLVer,
-                                           MilanLongInt *candidateMate,
-                                           MilanLongInt *verLocInd,
-                                           MilanLongInt *verLocPtr,
-                                           MilanLongInt StartIndex,
-                                           MilanLongInt EndIndex,
-                                           MilanLongInt *Mate,
-                                           vector<MilanLongInt> &GMate,
-                                           map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
-                                           MilanReal *edgeLocWeight,
-                                           MilanLongInt *myCardPtr,
-                                           MilanLongInt *msgIndPtr,
-                                           MilanLongInt *NumMessagesBundledPtr,
-                                           MilanLongInt *SPtr,
-                                           MilanLongInt *verDistance,
-                                           MilanLongInt *PCounter,
-                                           vector<MilanLongInt> &Counter,
-                                           MilanInt myRank,
-                                           MilanInt numProcs,
-                                           vector<MilanLongInt> &U,
-                                           vector<MilanLongInt> &privateU,
-                                           vector<MilanLongInt> &QLocalVtx,
-                                           vector<MilanLongInt> &QGhostVtx,
-                                           vector<MilanLongInt> &QMsgType,
-                                           vector<MilanInt> &QOwner,
-                                           vector<MilanLongInt> &privateQLocalVtx,
-                                           vector<MilanLongInt> &privateQGhostVtx,
-                                           vector<MilanLongInt> &privateQMsgType,
-                                           vector<MilanInt> &privateQOwner);
-
-    void PROCESS_CROSS_EDGE(MilanLongInt *edge,
-                            MilanLongInt *SPtr);
-
-    void processMatchedVertices(
-        MilanLongInt NLVer,
-        vector<MilanLongInt> &UChunkBeingProcessed,
-        vector<MilanLongInt> &U,
-        vector<MilanLongInt> &privateU,
-        MilanLongInt StartIndex,
-        MilanLongInt EndIndex,
-        MilanLongInt *myCardPtr,
-        MilanLongInt *msgIndPtr,
-        MilanLongInt *NumMessagesBundledPtr,
-        MilanLongInt *SPtr,
-        MilanLongInt *verLocPtr,
-        MilanLongInt *verLocInd,
-        MilanLongInt *verDistance,
-        MilanLongInt *PCounter,
-        vector<MilanLongInt> &Counter,
-        MilanInt myRank,
-        MilanInt numProcs,
-        MilanLongInt *candidateMate,
-        vector<MilanLongInt> &GMate,
-        MilanLongInt *Mate,
-        map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
-        MilanReal *edgeLocWeight,
-        vector<MilanLongInt> &QLocalVtx,
-        vector<MilanLongInt> &QGhostVtx,
-        vector<MilanLongInt> &QMsgType,
-        vector<MilanInt> &QOwner,
-        vector<MilanLongInt> &privateQLocalVtx,
-        vector<MilanLongInt> &privateQGhostVtx,
-        vector<MilanLongInt> &privateQMsgType,
-        vector<MilanInt> &privateQOwner);
-
-    void processMatchedVerticesAndSendMessages(
-        MilanLongInt NLVer,
-        vector<MilanLongInt> &UChunkBeingProcessed,
-        vector<MilanLongInt> &U,
-        vector<MilanLongInt> &privateU,
-        MilanLongInt StartIndex,
-        MilanLongInt EndIndex,
-        MilanLongInt *myCardPtr,
-        MilanLongInt *msgIndPtr,
-        MilanLongInt *NumMessagesBundledPtr,
-        MilanLongInt *SPtr,
-        MilanLongInt *verLocPtr,
-        MilanLongInt *verLocInd,
-        MilanLongInt *verDistance,
-        MilanLongInt *PCounter,
-        vector<MilanLongInt> &Counter,
-        MilanInt myRank,
-        MilanInt numProcs,
-        MilanLongInt *candidateMate,
-        vector<MilanLongInt> &GMate,
-        MilanLongInt *Mate,
-        map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
-        MilanReal *edgeLocWeight,
-        vector<MilanLongInt> &QLocalVtx,
-        vector<MilanLongInt> &QGhostVtx,
-        vector<MilanLongInt> &QMsgType,
-        vector<MilanInt> &QOwner,
-        vector<MilanLongInt> &privateQLocalVtx,
-        vector<MilanLongInt> &privateQGhostVtx,
-        vector<MilanLongInt> &privateQMsgType,
-        vector<MilanInt> &privateQOwner,
-        MPI_Comm comm,
-        MilanLongInt *msgActual,
-        vector<MilanLongInt> &Message);
-
-    void sendBundledMessages(MilanLongInt *numGhostEdgesPtr,
-                             MilanInt *BufferSizePtr,
-                             MilanLongInt *Buffer,
-                             vector<MilanLongInt> &PCumulative,
-                             vector<MilanLongInt> &PMessageBundle,
-                             vector<MilanLongInt> &PSizeInfoMessages,
-                             MilanLongInt *PCounter,
-                             MilanLongInt NumMessagesBundled,
-                             MilanLongInt *msgActualPtr,
-                             MilanLongInt *MessageIndexPtr,
-                             MilanInt numProcs,
-                             MilanInt myRank,
-                             MPI_Comm comm,
-                             vector<MilanLongInt> &QLocalVtx,
-                             vector<MilanLongInt> &QGhostVtx,
-                             vector<MilanLongInt> &QMsgType,
-                             vector<MilanInt> &QOwner,
-                             vector<MPI_Request> &SRequest,
-                             vector<MPI_Status> &SStatus);
-
-    void processMessages(
-        MilanLongInt NLVer,
-        MilanLongInt *Mate,
-        MilanLongInt *candidateMate,
-        map<MilanLongInt, MilanLongInt> &Ghost2LocalMap,
-        vector<MilanLongInt> &GMate,
-        vector<MilanLongInt> &Counter,
-        MilanLongInt StartIndex,
-        MilanLongInt EndIndex,
-        MilanLongInt *myCardPtr,
-        MilanLongInt *msgIndPtr,
-        MilanLongInt *msgActualPtr,
-        MilanReal *edgeLocWeight,
-        MilanLongInt *verDistance,
-        MilanLongInt *verLocPtr,
-        MilanLongInt k,
-        MilanLongInt *verLocInd,
-        MilanInt numProcs,
-        MilanInt myRank,
-        MPI_Comm comm,
-        vector<MilanLongInt> &Message,
-        MilanLongInt numGhostEdges,
-        MilanLongInt u,
-        MilanLongInt v,
-        MilanLongInt *SPtr,
-        vector<MilanLongInt> &U);
-
-    void extractUChunk(
-        vector<MilanLongInt> &UChunkBeingProcessed,
-        vector<MilanLongInt> &U,
-        vector<MilanLongInt> &privateU);
-
-    void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
-        MilanLongInt NLVer, MilanLongInt NLEdge,
-        MilanLongInt *verLocPtr, MilanLongInt *verLocInd, MilanReal *edgeLocWeight,
-        MilanLongInt *verDistance,
-        MilanLongInt *Mate,
-        MilanInt myRank, MilanInt numProcs, MPI_Comm comm,
-        MilanLongInt *msgIndSent, MilanLongInt *msgActualSent, MilanReal *msgPercent,
-        MilanReal *ph0_time, MilanReal *ph1_time, MilanReal *ph2_time,
-        MilanLongInt *ph1_card, MilanLongInt *ph2_card);
-#endif
+  
+  void salgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
+			 MilanLongInt NLVer, MilanLongInt NLEdge,
+			 MilanLongInt *verLocPtr, MilanLongInt *verLocInd, MilanFloat *edgeLocWeight,
+			 MilanLongInt *verDistance,
+			 MilanLongInt *Mate,
+			 MilanInt myRank, MilanInt numProcs, MPI_Comm comm,
+			 MilanLongInt *msgIndSent, MilanLongInt *msgActualSent, MilanReal *msgPercent,
+			 MilanReal *ph0_time, MilanReal *ph1_time, MilanReal *ph2_time,
+			 MilanLongInt *ph1_card, MilanLongInt *ph2_card);
 
 
-#ifndef OMP
-   //Function of find the owner of a ghost vertex using binary search:
-   inline MilanInt findOwnerOfGhost(MilanLongInt vtxIndex, MilanLongInt *mVerDistance,
-                                        MilanInt myRank, MilanInt numProcs);
-#endif
 
     void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateC(
         MilanLongInt NLVer, MilanLongInt NLEdge,

amgprec/impl/aggregator/algoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateC.cpp

@@ -1303,16 +1303,16 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateC(
 // SINGLE PRECISION VERSION
 
 void salgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateC(
-                                                              MilanLongInt NLVer, MilanLongInt NLEdge,
-                                                              MilanLongInt* verLocPtr, MilanLongInt* verLocInd,
-							      MilanFloat* edgeLocWeight,
-                                                              MilanLongInt* verDistance,
-                                                              MilanLongInt* Mate,
-                                                              MilanInt myRank, MilanInt numProcs, MPI_Comm comm,
-                                                              MilanLongInt* msgIndSent, MilanLongInt* msgActualSent,
-							      MilanReal* msgPercent,
-                                                              MilanReal* ph0_time, MilanReal* ph1_time, MilanReal* ph2_time,
-                                                              MilanLongInt* ph1_card, MilanLongInt* ph2_card ) {
+      MilanLongInt NLVer, MilanLongInt NLEdge,
+      MilanLongInt* verLocPtr, MilanLongInt* verLocInd,
+      MilanFloat* edgeLocWeight,
+      MilanLongInt* verDistance,
+      MilanLongInt* Mate,
+      MilanInt myRank, MilanInt numProcs, MPI_Comm comm,
+      MilanLongInt* msgIndSent, MilanLongInt* msgActualSent,
+      MilanReal* msgPercent,
+      MilanReal* ph0_time, MilanReal* ph1_time, MilanReal* ph2_time,
+      MilanLongInt* ph1_card, MilanLongInt* ph2_card ) {
 #if !defined(SERIAL_MPI)
 #ifdef PRINT_DEBUG_INFO_
     cout<<"\n("<<myRank<<")Within algoEdgeApproxDominatingEdgesLinearSearchMessageBundling()"; fflush(stdout);

amgprec/impl/aggregator/algoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP.cpp

@@ -1,5 +1,4 @@
 #include "MatchBoxPC.h"
-#ifdef OMP
 // ***********************************************************************
 //
 //        MatchboxP: A C++ library for approximate weighted matching
@@ -126,8 +125,10 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
     fflush(stdout);
 #endif
 
-    MilanLongInt StartIndex = verDistance[myRank];       // The starting vertex owned by the current rank
-    MilanLongInt EndIndex   = verDistance[myRank + 1] - 1; // The ending vertex owned by the current rank
+    // The starting vertex owned by the current rank
+    MilanLongInt StartIndex = verDistance[myRank];
+    // The ending vertex owned by the current rank
+    MilanLongInt EndIndex   = verDistance[myRank + 1] - 1;
 
     MPI_Status computeStatus;
 
@@ -145,7 +146,8 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
     // only one message will be sent in the initialization phase -
     // one of: REQUEST/FAILURE/SUCCESS
     vector<MilanLongInt> QLocalVtx, QGhostVtx, QMsgType;
-    vector<MilanInt> QOwner; // Changed by Fabio to be an integer, addresses needs to be integers!
+    // Changed by Fabio to be an integer, addresses needs to be integers!
+    vector<MilanInt> QOwner; 
 
     MilanLongInt *PCounter = new MilanLongInt[numProcs];
     for (int i = 0; i < numProcs; i++)
@@ -153,7 +155,8 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
 
     MilanLongInt NumMessagesBundled = 0;
     // TODO when the last computational section will be refactored this could be eliminated
-    MilanInt ghostOwner = 0; // Changed by Fabio to be an integer, addresses needs to be integers!
+    // Changed by Fabio to be an integer, addresses needs to be integers!
+    MilanInt ghostOwner = 0; 
     MilanLongInt *candidateMate = nullptr;
 #ifdef PRINT_DEBUG_INFO_
     cout << "\n(" << myRank << ")NV: " << NLVer << "  Edges: " << NLEdge;
@@ -168,9 +171,12 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
     MilanLongInt myCard = 0;
 
     // Build the Ghost Vertex Set: Vg
-    map<MilanLongInt, MilanLongInt> Ghost2LocalMap;       // Map each ghost vertex to a local vertex
-    vector<MilanLongInt> Counter;                         // Store the edge count for each ghost vertex
-    MilanLongInt numGhostVertices = 0, numGhostEdges = 0; // Number of Ghost vertices
+    // Map each ghost vertex to a local vertex
+    map<MilanLongInt, MilanLongInt> Ghost2LocalMap;
+    // Store the edge count for each ghost vertex
+    vector<MilanLongInt> Counter;
+    // Number of Ghost vertices
+    MilanLongInt numGhostVertices = 0, numGhostEdges = 0; 
 
 #ifdef PRINT_DEBUG_INFO_
     cout << "\n(" << myRank << ")About to compute Ghost Vertices...";
@@ -237,7 +243,7 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
      * PARALLEL_COMPUTE_CANDIDATE_MATE_B is now totally parallel.
      */
 
-    PARALLEL_COMPUTE_CANDIDATE_MATE_B(NLVer,
+    PARALLEL_COMPUTE_CANDIDATE_MATE_BD(NLVer,
                                       verLocPtr,
                                       verLocInd,
                                       myRank,
@@ -262,7 +268,7 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
      * TODO: Test when it's actually more efficient to execute this code
      *       in parallel.
      */
-    PARALLEL_PROCESS_EXPOSED_VERTEX_B(NLVer,
+    PARALLEL_PROCESS_EXPOSED_VERTEX_BD(NLVer,
                                       candidateMate,
                                       verLocInd,
                                       verLocPtr,
@@ -314,7 +320,7 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
     vector<MilanLongInt> UChunkBeingProcessed;
     UChunkBeingProcessed.reserve(UCHUNK);
 
-    processMatchedVertices(NLVer,
+    processMatchedVerticesD(NLVer,
                            UChunkBeingProcessed,
                            U,
                            privateU,
@@ -423,7 +429,7 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
         /////////////////////////// PROCESS MATCHED VERTICES //////////////////////////////
         ///////////////////////////////////////////////////////////////////////////////////
 
-      processMatchedVerticesAndSendMessages(NLVer,
+      processMatchedVerticesAndSendMessagesD(NLVer,
 					    UChunkBeingProcessed,
 					    U,
 					    privateU,
@@ -483,8 +489,8 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
         ///////////////////////////////////////////////////////////////////////////////////
         /////////////////////////// PROCESS MESSAGES //////////////////////////////////////
         ///////////////////////////////////////////////////////////////////////////////////
-
-        processMessages(NLVer,
+	//startTime = MPI_Wtime();
+        processMessagesD(NLVer,
                         Mate,
                         candidateMate,
                         Ghost2LocalMap,
@@ -549,6 +555,488 @@ void dalgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
     *ph2_card = myCard;                 // Cardinality at the end of Phase-2
 }
 // End of algoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMate
+
+void salgoDistEdgeApproxDomEdgesLinearSearchMesgBndlSmallMateCMP(
+    MilanLongInt NLVer, MilanLongInt NLEdge,
+    MilanLongInt *verLocPtr, MilanLongInt *verLocInd,
+    MilanFloat *edgeLocWeight,
+    MilanLongInt *verDistance,
+    MilanLongInt *Mate,
+    MilanInt myRank, MilanInt numProcs, MPI_Comm comm,
+    MilanLongInt *msgIndSent, MilanLongInt *msgActualSent,
+    MilanReal *msgPercent,
+    MilanReal *ph0_time, MilanReal *ph1_time, MilanReal *ph2_time,
+    MilanLongInt *ph1_card, MilanLongInt *ph2_card)
+{
+
+    /*
+     * verDistance: it's a vector long as the number of processors.
+     *              verDistance[i] contains the first node index of the i-th processor
+     *              verDistance[i + 1] contains the last node index of the i-th processor
+     * NLVer: number of elements in the LocPtr
+     * NLEdge: number of edges assigned to the current processor
+     *
+     * Contains the portion of matrix assigned to the processor in
+     * Yale notation
+     * verLocInd: contains the positions on row of the matrix
+     * verLocPtr: i-th value is the position of the first element on the i-th row and
+     *            i+1-th value is the position of the first element on the i+1-th row
+     */
+
+#if !defined(SERIAL_MPI)
+#ifdef PRINT_DEBUG_INFO_
+    cout << "\n(" << myRank << ")Within algoEdgeApproxDominatingEdgesLinearSearchMessageBundling()";
+    fflush(stdout);
+#endif
+
+#ifdef PRINT_DEBUG_INFO_
+    cout << "\n(" << myRank << ") verDistance [" ;
+    for (int i = 0; i < numProcs; i++)
+      cout << verDistance[i] << "," << verDistance[i+1];
+    cout  << "]\n";
+    fflush(stdout);
+#endif
+#ifdef DEBUG_HANG_
+    if (myRank == 0) {
+      cout << "\n(" << myRank << ") verDistance [" ;
+      for (int i = 0; i < numProcs; i++)
+	cout << verDistance[i] << "," ;
+      cout  << verDistance[numProcs]<< "]\n";
+    }
+    fflush(stdout);
+#endif
+
+    // The starting vertex owned by the current rank
+    MilanLongInt StartIndex = verDistance[myRank];
+    // The ending vertex owned by the current rank
+    MilanLongInt EndIndex   = verDistance[myRank + 1] - 1;
+
+    MPI_Status computeStatus;
+
+    MilanLongInt msgActual = 0, msgInd = 0;
+    MilanFloat heaviestEdgeWt = 0.0f; // Assumes positive weight
+    MilanReal startTime, finishTime;
+
+    startTime = MPI_Wtime();
+
+    // Data structures for sending and receiving messages:
+    vector<MilanLongInt> Message; // [ u, v, message_type ]
+    Message.resize(3, -1);
+    // Data structures for Message Bundling:
+    // Although up to two messages can be sent along any cross edge,
+    // only one message will be sent in the initialization phase -
+    // one of: REQUEST/FAILURE/SUCCESS
+    vector<MilanLongInt> QLocalVtx, QGhostVtx, QMsgType;
+    // Changed by Fabio to be an integer, addresses needs to be integers!
+    vector<MilanInt> QOwner; 
+
+    MilanLongInt *PCounter = new MilanLongInt[numProcs];
+    for (int i = 0; i < numProcs; i++)
+        PCounter[i] = 0;
+
+    MilanLongInt NumMessagesBundled = 0;
+    // TODO when the last computational section will be refactored this could be eliminated
+    // Changed by Fabio to be an integer, addresses needs to be integers!
+    MilanInt ghostOwner = 0; 
+    MilanLongInt *candidateMate = nullptr;
+#ifdef PRINT_DEBUG_INFO_
+    cout << "\n(" << myRank << ")NV: " << NLVer << "  Edges: " << NLEdge;
+    fflush(stdout);
+    cout << "\n(" << myRank << ")StartIndex: " << StartIndex << "  EndIndex: " << EndIndex;
+    fflush(stdout);
+#endif
+    // Other Variables:
+    MilanLongInt u = -1, v = -1, w = -1, i = 0;
+    MilanLongInt k = -1, adj1 = -1, adj2 = -1;
+    MilanLongInt k1 = -1, adj11 = -1, adj12 = -1;
+    MilanLongInt myCard = 0;
+
+    // Build the Ghost Vertex Set: Vg
+    // Map each ghost vertex to a local vertex
+    map<MilanLongInt, MilanLongInt> Ghost2LocalMap;
+    // Store the edge count for each ghost vertex
+    vector<MilanLongInt> Counter;
+    // Number of Ghost vertices
+    MilanLongInt numGhostVertices = 0, numGhostEdges = 0; 
+
+#ifdef PRINT_DEBUG_INFO_
+    cout << "\n(" << myRank << ")About to compute Ghost Vertices...";
+    fflush(stdout);
+#endif
+#ifdef DEBUG_HANG_
+    if (myRank == 0)
+        cout << "\n(" << myRank << ")About to compute Ghost Vertices...";
+    fflush(stdout);
+#endif
+
+    // Define Adjacency Lists for Ghost Vertices:
+    // cout<<"Building Ghost data structures ... \n\n";
+    vector<MilanLongInt> verGhostPtr, verGhostInd, tempCounter;
+    // Mate array for ghost vertices:
+    vector<MilanLongInt> GMate; // Proportional to the number of ghost vertices
+    MilanLongInt S;
+    MilanLongInt privateMyCard = 0;
+    vector<MilanLongInt> PCumulative, PMessageBundle, PSizeInfoMessages;
+    vector<MPI_Request> SRequest;  // Requests that are used for each send message
+    vector<MPI_Status> SStatus;    // Status of sent messages, used in MPI_Wait
+    MilanLongInt MessageIndex = 0; // Pointer for current message
+    MilanInt BufferSize;
+    MilanLongInt *Buffer;
+
+    vector<MilanLongInt> privateQLocalVtx, privateQGhostVtx, privateQMsgType;
+    vector<MilanInt> privateQOwner;
+    vector<MilanLongInt> U, privateU;
+
+    
+    initialize(NLVer, NLEdge, StartIndex,
+               EndIndex, &numGhostEdges,
+               &numGhostVertices, &S,
+               verLocInd, verLocPtr,
+               Ghost2LocalMap, Counter,
+               verGhostPtr, verGhostInd,
+               tempCounter, GMate,
+               Message, QLocalVtx,
+               QGhostVtx, QMsgType, QOwner,
+               candidateMate, U,
+               privateU,
+               privateQLocalVtx,
+               privateQGhostVtx,
+               privateQMsgType,
+               privateQOwner);
+
+    finishTime = MPI_Wtime();
+    *ph0_time = finishTime - startTime; // Time taken for Phase-0: Initialization
+#ifdef DEBUG_HANG_
+    cout << myRank << " Finished initialization" << endl;
+    fflush(stdout);
+#endif
+
+    startTime = MPI_Wtime();
+
+    /////////////////////////////////////////////////////////////////////////////////////////
+    //////////////////////////////////// INITIALIZATION /////////////////////////////////////
+    /////////////////////////////////////////////////////////////////////////////////////////
+    // Compute the Initial Matching Set:
+
+    /*
+     * OMP PARALLEL_COMPUTE_CANDIDATE_MATE_B has been splitted from
+     * PARALLEL_PROCESS_EXPOSED_VERTEX_B in order to better parallelize
+     * the two.
+     * PARALLEL_COMPUTE_CANDIDATE_MATE_B is now totally parallel.
+     */
+
+    PARALLEL_COMPUTE_CANDIDATE_MATE_BS(NLVer,
+                                      verLocPtr,
+                                      verLocInd,
+                                      myRank,
+                                      edgeLocWeight,
+                                      candidateMate);
+
+#ifdef DEBUG_HANG_
+    cout << myRank << " Finished Exposed Vertex" << endl;
+    fflush(stdout);
+#if 0
+    cout << myRank << " candidateMate after parallelCompute " <<endl;
+    for (int i=0; i<NLVer; i++) {
+      cout << candidateMate[i] << " " ;
+    }
+    cout << endl;
+#endif
+#endif
+    /*
+     * PARALLEL_PROCESS_EXPOSED_VERTEX_B
+     * TODO: write comment
+     *
+     * TODO: Test when it's actually more efficient to execute this code
+     *       in parallel.
+     */
+    PARALLEL_PROCESS_EXPOSED_VERTEX_BS(NLVer,
+                                      candidateMate,
+                                      verLocInd,
+                                      verLocPtr,
+                                      StartIndex,
+                                      EndIndex,
+                                      Mate,
+                                      GMate,
+                                      Ghost2LocalMap,
+                                      edgeLocWeight,
+                                      &myCard,
+                                      &msgInd,
+                                      &NumMessagesBundled,
+                                      &S,
+                                      verDistance,
+                                      PCounter,
+                                      Counter,
+                                      myRank,
+                                      numProcs,
+                                      U,
+                                      privateU,
+                                      QLocalVtx,
+                                      QGhostVtx,
+                                      QMsgType,
+                                      QOwner,
+                                      privateQLocalVtx,
+                                      privateQGhostVtx,
+                                      privateQMsgType,
+                                      privateQOwner);
+
+    tempCounter.clear(); // Do not need this any more
+
+#ifdef DEBUG_HANG_
+    cout << myRank << " Finished Exposed Vertex" << endl;
+    fflush(stdout);
+#if 0
+    cout << myRank << " Mate after Exposed Vertices " <<endl;
+    for (int i=0; i<NLVer; i++) {
+      cout << Mate[i] << " " ;
+    }
+    cout << endl;
+#endif
+#endif
+
+    ///////////////////////////////////////////////////////////////////////////////////
+    /////////////////////////// PROCESS MATCHED VERTICES //////////////////////////////
+    ///////////////////////////////////////////////////////////////////////////////////
+
+    // TODO what would be the optimal UCHUNK
+    vector<MilanLongInt> UChunkBeingProcessed;
+    UChunkBeingProcessed.reserve(UCHUNK);
+
+    processMatchedVerticesS(NLVer,
+                           UChunkBeingProcessed,
+                           U,
+                           privateU,
+                           StartIndex,
+                           EndIndex,
+                           &myCard,
+                           &msgInd,
+                           &NumMessagesBundled,
+                           &S,
+                           verLocPtr,
+                           verLocInd,
+                           verDistance,
+                           PCounter,
+                           Counter,
+                           myRank,
+                           numProcs,
+                           candidateMate,
+                           GMate,
+                           Mate,
+                           Ghost2LocalMap,
+                           edgeLocWeight,
+                           QLocalVtx,
+                           QGhostVtx,
+                           QMsgType,
+                           QOwner,
+                           privateQLocalVtx,
+                           privateQGhostVtx,
+                           privateQMsgType,
+                           privateQOwner);
+
+
+#ifdef DEBUG_HANG_
+    cout << myRank << " Finished Process Vertices" << endl;
+    fflush(stdout);
+#if 0
+    cout << myRank << " Mate after Matched Vertices " <<endl;
+    for (int i=0; i<NLVer; i++) {
+      cout << Mate[i] << " " ;
+    }
+    cout << endl;
+#endif
+#endif
+
+    /////////////////////////////////////////////////////////////////////////////////////////
+    ///////////////////////////// SEND BUNDLED MESSAGES /////////////////////////////////////
+    /////////////////////////////////////////////////////////////////////////////////////////
+
+    sendBundledMessages(&numGhostEdges,
+                        &BufferSize,
+                        Buffer,
+                        PCumulative,
+                        PMessageBundle,
+                        PSizeInfoMessages,
+                        PCounter,
+                        NumMessagesBundled,
+                        &msgActual,
+                        &MessageIndex,
+                        numProcs,
+                        myRank,
+                        comm,
+                        QLocalVtx,
+                        QGhostVtx,
+                        QMsgType,
+                        QOwner,
+                        SRequest,
+                        SStatus);
+
+    ///////////////////////// END OF SEND BUNDLED MESSAGES //////////////////////////////////
+
+    finishTime = MPI_Wtime();
+    *ph1_time = finishTime - startTime; // Time taken for Phase-1
+
+#ifdef DEBUG_HANG_
+    cout << myRank << " Finished sendBundles" << endl;
+    fflush(stdout);
+#endif
+
+    *ph1_card = myCard;                 // Cardinality at the end of Phase-1
+    startTime = MPI_Wtime();
+    /////////////////////////////////////////////////////////////////////////////////////////
+    //////////////////////////////////////// MAIN LOOP //////////////////////////////////////
+    /////////////////////////////////////////////////////////////////////////////////////////
+    // Main While Loop:
+#ifdef PRINT_DEBUG_INFO_
+    cout << "\n(" << myRank << "=========================************===============================" << endl;
+    fflush(stdout);
+    fflush(stdout);
+#endif
+#ifdef PRINT_DEBUG_INFO_
+    cout << "\n(" << myRank << ")Entering While(true) loop..";
+    fflush(stdout);
+#endif
+#ifdef PRINT_DEBUG_INFO_
+    cout << "\n(" << myRank << "=========================************===============================" << endl;
+    fflush(stdout);
+    fflush(stdout);
+#endif
+
+    while (true)  {
+#ifdef DEBUG_HANG_
+      //if (myRank == 0)
+      cout << "\n(" << myRank << ") Main loop" << endl;
+      fflush(stdout);
+#endif
+        ///////////////////////////////////////////////////////////////////////////////////
+        /////////////////////////// PROCESS MATCHED VERTICES //////////////////////////////
+        ///////////////////////////////////////////////////////////////////////////////////
+
+      processMatchedVerticesAndSendMessagesS(NLVer,
+					    UChunkBeingProcessed,
+					    U,
+					    privateU,
+					    StartIndex,
+					    EndIndex,
+					    &myCard,
+					    &msgInd,
+					    &NumMessagesBundled,
+					    &S,
+					    verLocPtr,
+					    verLocInd,
+					    verDistance,
+					    PCounter,
+					    Counter,
+					    myRank,
+					    numProcs,
+					    candidateMate,
+					    GMate,
+					    Mate,
+					    Ghost2LocalMap,
+					    edgeLocWeight,
+					    QLocalVtx,
+					    QGhostVtx,
+					    QMsgType,
+					    QOwner,
+					    privateQLocalVtx,
+					    privateQGhostVtx,
+					    privateQMsgType,
+					    privateQOwner,
+					    comm,
+					    &msgActual,
+					    Message);
+
+      ///////////////////////// END OF PROCESS MATCHED VERTICES /////////////////////////
+
+        //// BREAK IF NO MESSAGES EXPECTED /////////
+#ifdef DEBUG_HANG_
+#if 0
+      cout << myRank << " Mate after ProcessMatchedAndSend phase "<<S <<endl;
+      for (int i=0; i<NLVer; i++) {
+	cout << Mate[i] << " " ;
+      }
+      cout << endl;
+#endif
+#endif
+#ifdef PRINT_DEBUG_INFO_
+        cout << "\n(" << myRank << ")Deciding whether to break: S= " << S << endl;
+#endif
+
+        if (S == 0) {
+#ifdef DEBUG_HANG_
+	  cout << "\n(" << myRank << ") Breaking out" << endl;
+	  fflush(stdout);
+#endif
+	  break;
+        }
+        ///////////////////////////////////////////////////////////////////////////////////
+        /////////////////////////// PROCESS MESSAGES //////////////////////////////////////
+        ///////////////////////////////////////////////////////////////////////////////////
+
+        processMessagesS(NLVer,
+                        Mate,
+                        candidateMate,
+                        Ghost2LocalMap,
+                        GMate,
+                        Counter,
+                        StartIndex,
+                        EndIndex,
+                        &myCard,
+                        &msgInd,
+                        &msgActual,
+                        edgeLocWeight,
+                        verDistance,
+                        verLocPtr,
+                        k,
+                        verLocInd,
+                        numProcs,
+                        myRank,
+                        comm,
+                        Message,
+                        numGhostEdges,
+                        u,
+                        v,
+                        &S,
+                        U);
+
+        ///////////////////////// END OF PROCESS MESSAGES /////////////////////////////////
+#ifdef DEBUG_HANG_
+#if 0
+      cout << myRank << " Mate after ProcessMessages phase "<<S <<endl;
+      for (int i=0; i<NLVer; i++) {
+	cout << Mate[i] << " " ;
+      }
+      cout << endl;
+#endif
+#endif
+#ifdef PRINT_DEBUG_INFO_
+        cout << "\n(" << myRank << ")Finished Message processing phase: S= " << S;
+        fflush(stdout);
+        cout << "\n(" << myRank << ")** SENT     : ACTUAL= " << msgActual;
+        fflush(stdout);
+        cout << "\n(" << myRank << ")** SENT     : INDIVIDUAL= " << msgInd << endl;
+        fflush(stdout);
+#endif
+    } // End of while (true)
+
+    clean(NLVer,
+          myRank,
+          MessageIndex,
+          SRequest,
+          SStatus,
+          BufferSize,
+          Buffer,
+          msgActual,
+          msgActualSent,
+          msgInd,
+          msgIndSent,
+          NumMessagesBundled,
+          msgPercent);
+
+    finishTime = MPI_Wtime();
+    *ph2_time = finishTime - startTime; // Time taken for Phase-2
+    *ph2_card = myCard;                 // Cardinality at the end of Phase-2
+}
+
 #endif
 
 #endif

amgprec/impl/aggregator/amg_c_soc1_map_bld.F90

@@ -275,7 +275,8 @@ subroutine amg_c_soc1_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,in
             disjoint = all(ilaggr(icol(1:ip)) == -(nr+1)).or.(ip==0)
             if (disjoint) then
               locnaggr(kk)     = locnaggr(kk) + 1
-              itmp = (bnds(kk)-1+locnaggr(kk))*nths+kk
+              itmp = (bnds(kk)-1+locnaggr(kk)) !be careful about overflow
+              itmp = itmp*nths+kk
               if (itmp < (bnds(kk)-1+locnaggr(kk))) then
                 !$omp atomic update
                 info = max(12345678,info)

amgprec/impl/aggregator/amg_c_soc2_map_bld.F90

@@ -309,7 +309,8 @@ subroutine amg_c_soc2_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,in
             !
             if (disjoint) then
               locnaggr(kk)     = locnaggr(kk) + 1
-              itmp = (bnds(kk)-1+locnaggr(kk))*nths+kk
+              itmp = (bnds(kk)-1+locnaggr(kk)) !be careful about overflow
+              itmp = itmp*nths+kk
               if (itmp < (bnds(kk)-1+locnaggr(kk))) then
                 !$omp atomic update
                 info = max(12345678,info)

amgprec/impl/aggregator/amg_caggrmat_minnrg_bld.f90

@@ -300,7 +300,7 @@ subroutine amg_caggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr,parms,&
 !!$        endif
 !!$      enddo
 !!$      if (jd == -1) then 
-!!$        write(0,*) 'Wrong input: we need the diagonal!!!!', i
+!!$                write(0,*) name,': Warning: there is no diagonal element', i
 !!$      else
 !!$        acsrf%val(jd)=acsrf%val(jd)-tmp
 !!$      end if

amgprec/impl/aggregator/amg_caggrmat_smth_bld.f90

@@ -230,7 +230,7 @@ subroutine amg_caggrmat_smth_bld(a,desc_a,ilaggr,nlaggr,parms,&
 
       enddo
       if (jd == -1) then 
-        write(0,*) 'Wrong input: we need the diagonal!!!!', i
+        write(0,*) name,': Warning: there is no diagonal element', i
       else
         acsrf%val(jd)=acsrf%val(jd)-tmp
       end if

amgprec/impl/aggregator/amg_d_parmatch_smth_bld.F90

@@ -246,7 +246,7 @@ subroutine amg_d_parmatch_smth_bld(ag,a,desc_a,ilaggr,nlaggr,parms,&
 
       enddo
       if (jd == -1) then
-        write(0,*) 'Wrong input: we need the diagonal!!!!', i
+        write(0,*) name,': Warning: there is no diagonal element', i
       else
         acsrf%val(jd)=acsrf%val(jd)-tmp
       end if

amgprec/impl/aggregator/amg_d_soc1_map_bld.F90

@@ -275,7 +275,8 @@ subroutine amg_d_soc1_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,in
             disjoint = all(ilaggr(icol(1:ip)) == -(nr+1)).or.(ip==0)
             if (disjoint) then
               locnaggr(kk)     = locnaggr(kk) + 1
-              itmp = (bnds(kk)-1+locnaggr(kk))*nths+kk
+              itmp = (bnds(kk)-1+locnaggr(kk)) !be careful about overflow
+              itmp = itmp*nths+kk
               if (itmp < (bnds(kk)-1+locnaggr(kk))) then
                 !$omp atomic update
                 info = max(12345678,info)

amgprec/impl/aggregator/amg_d_soc2_map_bld.F90

@@ -309,7 +309,8 @@ subroutine amg_d_soc2_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,in
             !
             if (disjoint) then
               locnaggr(kk)     = locnaggr(kk) + 1
-              itmp = (bnds(kk)-1+locnaggr(kk))*nths+kk
+              itmp = (bnds(kk)-1+locnaggr(kk)) !be careful about overflow
+              itmp = itmp*nths+kk
               if (itmp < (bnds(kk)-1+locnaggr(kk))) then
                 !$omp atomic update
                 info = max(12345678,info)

amgprec/impl/aggregator/amg_daggrmat_minnrg_bld.f90

@@ -300,7 +300,7 @@ subroutine amg_daggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr,parms,&
 !!$        endif
 !!$      enddo
 !!$      if (jd == -1) then 
-!!$        write(0,*) 'Wrong input: we need the diagonal!!!!', i
+!!$                write(0,*) name,': Warning: there is no diagonal element', i
 !!$      else
 !!$        acsrf%val(jd)=acsrf%val(jd)-tmp
 !!$      end if

amgprec/impl/aggregator/amg_daggrmat_smth_bld.f90

@@ -230,7 +230,7 @@ subroutine amg_daggrmat_smth_bld(a,desc_a,ilaggr,nlaggr,parms,&
 
       enddo
       if (jd == -1) then 
-        write(0,*) 'Wrong input: we need the diagonal!!!!', i
+        write(0,*) name,': Warning: there is no diagonal element', i
       else
         acsrf%val(jd)=acsrf%val(jd)-tmp
       end if

amgprec/impl/aggregator/amg_s_parmatch_smth_bld.F90

@@ -246,7 +246,7 @@ subroutine amg_s_parmatch_smth_bld(ag,a,desc_a,ilaggr,nlaggr,parms,&
 
       enddo
       if (jd == -1) then
-        write(0,*) 'Wrong input: we need the diagonal!!!!', i
+        write(0,*) name,': Warning: there is no diagonal element', i
       else
         acsrf%val(jd)=acsrf%val(jd)-tmp
       end if

amgprec/impl/aggregator/amg_s_soc1_map_bld.F90

@@ -275,7 +275,8 @@ subroutine amg_s_soc1_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,in
             disjoint = all(ilaggr(icol(1:ip)) == -(nr+1)).or.(ip==0)
             if (disjoint) then
               locnaggr(kk)     = locnaggr(kk) + 1
-              itmp = (bnds(kk)-1+locnaggr(kk))*nths+kk
+              itmp = (bnds(kk)-1+locnaggr(kk)) !be careful about overflow
+              itmp = itmp*nths+kk
               if (itmp < (bnds(kk)-1+locnaggr(kk))) then
                 !$omp atomic update
                 info = max(12345678,info)

amgprec/impl/aggregator/amg_s_soc2_map_bld.F90

@@ -309,7 +309,8 @@ subroutine amg_s_soc2_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,in
             !
             if (disjoint) then
               locnaggr(kk)     = locnaggr(kk) + 1
-              itmp = (bnds(kk)-1+locnaggr(kk))*nths+kk
+              itmp = (bnds(kk)-1+locnaggr(kk)) !be careful about overflow
+              itmp = itmp*nths+kk
               if (itmp < (bnds(kk)-1+locnaggr(kk))) then
                 !$omp atomic update
                 info = max(12345678,info)

amgprec/impl/aggregator/amg_saggrmat_minnrg_bld.f90

@@ -300,7 +300,7 @@ subroutine amg_saggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr,parms,&
 !!$        endif
 !!$      enddo
 !!$      if (jd == -1) then 
-!!$        write(0,*) 'Wrong input: we need the diagonal!!!!', i
+!!$                write(0,*) name,': Warning: there is no diagonal element', i
 !!$      else
 !!$        acsrf%val(jd)=acsrf%val(jd)-tmp
 !!$      end if

amgprec/impl/aggregator/amg_saggrmat_smth_bld.f90

@@ -230,7 +230,7 @@ subroutine amg_saggrmat_smth_bld(a,desc_a,ilaggr,nlaggr,parms,&
 
       enddo
       if (jd == -1) then 
-        write(0,*) 'Wrong input: we need the diagonal!!!!', i
+        write(0,*) name,': Warning: there is no diagonal element', i
       else
         acsrf%val(jd)=acsrf%val(jd)-tmp
       end if

amgprec/impl/aggregator/amg_z_soc1_map_bld.F90

@@ -275,7 +275,8 @@ subroutine amg_z_soc1_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,in
             disjoint = all(ilaggr(icol(1:ip)) == -(nr+1)).or.(ip==0)
             if (disjoint) then
               locnaggr(kk)     = locnaggr(kk) + 1
-              itmp = (bnds(kk)-1+locnaggr(kk))*nths+kk
+              itmp = (bnds(kk)-1+locnaggr(kk)) !be careful about overflow
+              itmp = itmp*nths+kk
               if (itmp < (bnds(kk)-1+locnaggr(kk))) then
                 !$omp atomic update
                 info = max(12345678,info)

amgprec/impl/aggregator/amg_z_soc2_map_bld.F90

@@ -309,7 +309,8 @@ subroutine amg_z_soc2_map_bld(iorder,theta,clean_zeros,a,desc_a,nlaggr,ilaggr,in
             !
             if (disjoint) then
               locnaggr(kk)     = locnaggr(kk) + 1
-              itmp = (bnds(kk)-1+locnaggr(kk))*nths+kk
+              itmp = (bnds(kk)-1+locnaggr(kk)) !be careful about overflow
+              itmp = itmp*nths+kk
               if (itmp < (bnds(kk)-1+locnaggr(kk))) then
                 !$omp atomic update
                 info = max(12345678,info)

amgprec/impl/aggregator/amg_zaggrmat_minnrg_bld.f90

@@ -300,7 +300,7 @@ subroutine amg_zaggrmat_minnrg_bld(a,desc_a,ilaggr,nlaggr,parms,&
 !!$        endif
 !!$      enddo
 !!$      if (jd == -1) then 
-!!$        write(0,*) 'Wrong input: we need the diagonal!!!!', i
+!!$                write(0,*) name,': Warning: there is no diagonal element', i
 !!$      else
 !!$        acsrf%val(jd)=acsrf%val(jd)-tmp
 !!$      end if

amgprec/impl/aggregator/amg_zaggrmat_smth_bld.f90

@@ -230,7 +230,7 @@ subroutine amg_zaggrmat_smth_bld(a,desc_a,ilaggr,nlaggr,parms,&
 
       enddo
       if (jd == -1) then 
-        write(0,*) 'Wrong input: we need the diagonal!!!!', i
+        write(0,*) name,': Warning: there is no diagonal element', i
       else
         acsrf%val(jd)=acsrf%val(jd)-tmp
       end if

amgprec/impl/aggregator/clean.cpp

@@ -1,6 +1,6 @@
 #include "MatchBoxPC.h"
-#ifdef OMP
 // TODO comment
+#if !defined(SERIAL_MPI)
 
 void clean(MilanLongInt NLVer,
            MilanInt myRank,

amgprec/impl/aggregator/computeCandidateMate.cpp

@@ -1,5 +1,4 @@
 #include "MatchBoxPC.h"
-#ifdef OMP
 /**
  * Execute the research fr the Candidate Mate without controlling if the vertices are already matched.
  * Returns the vertices with the highest weight
@@ -9,7 +8,8 @@
  * @param edgeLocWeight
  * @return
  */
-MilanLongInt firstComputeCandidateMate(MilanLongInt adj1,
+
+MilanLongInt firstComputeCandidateMateD(MilanLongInt adj1,
                                        MilanLongInt adj2,
                                        MilanLongInt *verLocInd,
                                        MilanReal *edgeLocWeight)
@@ -42,7 +42,7 @@ MilanLongInt firstComputeCandidateMate(MilanLongInt adj1,
  * @param Ghost2LocalMap
  * @return
  */
-MilanLongInt computeCandidateMate(MilanLongInt adj1,
+MilanLongInt computeCandidateMateD(MilanLongInt adj1,
                                   MilanLongInt adj2,
                                   MilanReal *edgeLocWeight,
                                   MilanLongInt k,
@@ -71,4 +71,68 @@ MilanLongInt computeCandidateMate(MilanLongInt adj1,
 
     return w;
 }
-#endif
+
+
+MilanLongInt firstComputeCandidateMateS(MilanLongInt adj1,
+                                       MilanLongInt adj2,
+                                       MilanLongInt *verLocInd,
+                                       MilanFloat *edgeLocWeight)
+{
+    MilanInt w = -1;
+    MilanFloat heaviestEdgeWt = 0.0f; // Assign the smallest 
+    int finalK;
+    for (int k = adj1; k < adj2; k++) {
+      if ((edgeLocWeight[k] > heaviestEdgeWt) ||
+	  ((edgeLocWeight[k] == heaviestEdgeWt) && (w < verLocInd[k])))     {
+	heaviestEdgeWt = edgeLocWeight[k];
+	w = verLocInd[k];
+	finalK = k;
+      }
+    } // End of for loop
+    return finalK;
+}
+
+/**
+ * //TODO documentation
+ * @param adj1
+ * @param adj2
+ * @param edgeLocWeight
+ * @param k
+ * @param verLocInd
+ * @param StartIndex
+ * @param EndIndex
+ * @param GMate
+ * @param Mate
+ * @param Ghost2LocalMap
+ * @return
+ */
+MilanLongInt computeCandidateMateS(MilanLongInt adj1,
+                                  MilanLongInt adj2,
+				  MilanFloat *edgeLocWeight,
+                                  MilanLongInt k,
+                                  MilanLongInt *verLocInd,
+                                  MilanLongInt StartIndex,
+                                  MilanLongInt EndIndex,
+                                  vector<MilanLongInt> &GMate,
+                                  MilanLongInt *Mate,
+                                  map<MilanLongInt, MilanLongInt> &Ghost2LocalMap)
+{
+    // Start: PARALLEL_COMPUTE_CANDIDATE_MATE_B(v)
+
+    MilanInt w = -1;
+    MilanFloat heaviestEdgeWt = 0.0f; // Assign the smallest Value
+    for (k = adj1; k < adj2; k++)   {
+      if (isAlreadyMatched(verLocInd[k], StartIndex, EndIndex, GMate, Mate, Ghost2LocalMap))
+	continue;
+
+      if ((edgeLocWeight[k] > heaviestEdgeWt) ||
+	  ((edgeLocWeight[k] == heaviestEdgeWt) && (w < verLocInd[k])))     {
+	heaviestEdgeWt = edgeLocWeight[k];
+	w = verLocInd[k];
+      }
+    } // End of for loop
+      //  End: PARALLEL_COMPUTE_CANDIDATE_MATE_B(v)
+
+    return w;
+}
+

amgprec/impl/aggregator/extractUChunk.cpp

@@ -1,5 +1,4 @@
 #include "MatchBoxPC.h"
-#ifdef OMP
 void extractUChunk(
     vector<MilanLongInt> &UChunkBeingProcessed,
     vector<MilanLongInt> &U,
@@ -29,4 +28,3 @@ void extractUChunk(
 
     } // End of critical U // End of critical U
 }
-#endif

amgprec/impl/aggregator/findOwnerOfGhost.cpp

@@ -1,5 +1,4 @@
 #include "MatchBoxPC.h"
-#ifdef OMP
 /// Find the owner of a ghost node:
 MilanInt findOwnerOfGhost(MilanLongInt vtxIndex, MilanLongInt *mVerDistance,
                           MilanInt myRank, MilanInt numProcs)
@@ -27,4 +26,3 @@ MilanInt findOwnerOfGhost(MilanLongInt vtxIndex, MilanLongInt *mVerDistance,
 
   return Current;
 } // End of findOwnerOfGhost()
-#endif

amgprec/impl/aggregator/initialize.cpp

@@ -1,5 +1,4 @@
 #include "MatchBoxPC.h"
-#ifdef OMP
 void initialize(MilanLongInt NLVer, MilanLongInt NLEdge,
                 MilanLongInt StartIndex, MilanLongInt EndIndex,
                 MilanLongInt *numGhostEdges,
@@ -302,4 +301,3 @@ void initialize(MilanLongInt NLVer, MilanLongInt NLEdge,
         } // End of single region
     }     // End of parallel region
 }
-#endif

amgprec/impl/aggregator/isAlreadyMatched.cpp

@@ -1,5 +1,4 @@
 #include "MatchBoxPC.h"
-#ifdef OMP
 /**
  * //TODO documentation
  * @param k
@@ -44,4 +43,3 @@ bool isAlreadyMatched(MilanLongInt node,
 
     return val >= 0; // Already matched
 }
-#endif

amgprec/impl/aggregator/parallelComputeCandidateMateB.cpp

@@ -1,6 +1,6 @@
 #include "MatchBoxPC.h"
-#ifdef OMP
-void PARALLEL_COMPUTE_CANDIDATE_MATE_B(MilanLongInt NLVer,
+
+void PARALLEL_COMPUTE_CANDIDATE_MATE_BD(MilanLongInt NLVer,
                                               MilanLongInt *verLocPtr,
                                               MilanLongInt *verLocInd,
                                               MilanInt myRank,
@@ -20,9 +20,37 @@ void PARALLEL_COMPUTE_CANDIDATE_MATE_B(MilanLongInt NLVer,
             fflush(stdout);
 #endif
             // Start: PARALLEL_COMPUTE_CANDIDATE_MATE_B(v)
-            candidateMate[v] = firstComputeCandidateMate(verLocPtr[v], verLocPtr[v + 1], verLocInd, edgeLocWeight);
+            candidateMate[v] = firstComputeCandidateMateD(verLocPtr[v], verLocPtr[v + 1],
+							 verLocInd, edgeLocWeight);
             // End: PARALLEL_COMPUTE_CANDIDATE_MATE_B(v)
         }
     }
 }
+
+void PARALLEL_COMPUTE_CANDIDATE_MATE_BS(MilanLongInt NLVer,
+                                              MilanLongInt *verLocPtr,
+                                              MilanLongInt *verLocInd,
+                                              MilanInt myRank,
+                                              MilanFloat *edgeLocWeight,
+                                              MilanLongInt *candidateMate)
+{
+
+    MilanLongInt v = -1;
+
+#pragma omp parallel private(v) default(shared) num_threads(NUM_THREAD)
+    {
+
+#pragma omp for schedule(static)
+        for (v = 0; v < NLVer; v++) {
+#ifdef PRINT_DEBUG_INFO_
+            cout << "\n(" << myRank << ")Processing: " << v + StartIndex << endl;
+            fflush(stdout);
 #endif
+            // Start: PARALLEL_COMPUTE_CANDIDATE_MATE_B(v)
+            candidateMate[v] = firstComputeCandidateMateS(verLocPtr[v], verLocPtr[v + 1],
+							 verLocInd, edgeLocWeight);
+            // End: PARALLEL_COMPUTE_CANDIDATE_MATE_B(v)
+        }
+    }
+}
+

amgprec/impl/aggregator/processCrossEdge.cpp

@@ -1,5 +1,4 @@
 #include "MatchBoxPC.h"
-#ifdef OMP
 void PROCESS_CROSS_EDGE(MilanLongInt *edge,
                         MilanLongInt *S)
 {
@@ -22,4 +21,3 @@ void PROCESS_CROSS_EDGE(MilanLongInt *edge,
 
     // End: PARALLEL_PROCESS_CROSS_EDGE_B
 }
-#endif