Skip to content

Commit

Permalink
Merge Pull Request #9807 from trilinos/Trilinos/master_merge_20211013…
Browse files Browse the repository at this point in the history
…_000549

Automatically Merged using Trilinos Master Merge AutoTester
PR Title: Trilinos Master Merge PR Generator: Auto PR created to promote from master_merge_20211013_000549 branch to master
PR Author: trilinos-autotester
  • Loading branch information
trilinos-autotester authored Oct 13, 2021
2 parents 45f7e16 + eebd63c commit 5a67b70
Show file tree
Hide file tree
Showing 83 changed files with 4,451 additions and 4,507 deletions.
10 changes: 8 additions & 2 deletions .github/workflows/muelu_label_to_project.yml
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
name: Auto assign MueLu labelled issues to MueLu Project
name: Auto assign MueLu or Tpetra labelled issues to MueLu or Tpetra Projects

on:
issues:
Expand All @@ -9,7 +9,7 @@ env:

jobs:
assign_one_project:
name: Assign to MueLu Project
name: Assign to MueLu or Tpetra Project
runs-on: ubuntu-latest
steps:
- name: Add to MueLu Project
Expand All @@ -18,3 +18,9 @@ jobs:
with:
project: 'https://github.com/trilinos/Trilinos/projects/5'
column_name: 'Feature request'
- name: Add to Tpetra Project
uses: srggrs/[email protected]
if: contains(github.event.label.name, 'Tpetra')
with:
project: 'https://github.com/trilinos/Trilinos/projects/3'
column_name: 'Backlog'
64 changes: 64 additions & 0 deletions packages/amesos2/src/Amesos2_Superlu_FunctionMap.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -94,6 +94,11 @@ at_plus_a(

namespace S { // single-precision real definitions

extern float slangs (char *, SLU::SuperMatrix *);

extern void sgscon (char *, SuperMatrix *, SuperMatrix *,
float, float *, SuperLUStat_t*, int *);

#ifdef HAVE_AMESOS2_SUPERLU5_API
typedef struct {
int *xsup; /* supernode and column mapping */
Expand Down Expand Up @@ -176,6 +181,11 @@ at_plus_a(

namespace D { // double-precision real definitions

extern double dlangs (char *, SLU::SuperMatrix *);

extern void dgscon (char *, SuperMatrix *, SuperMatrix *,
double, double *, SuperLUStat_t*, int *);

#ifdef HAVE_AMESOS2_SUPERLU5_API
typedef struct {
int *xsup; /* supernode and column mapping */
Expand Down Expand Up @@ -259,6 +269,11 @@ at_plus_a(
#ifdef HAVE_TEUCHOS_COMPLEX
namespace C { // single-precision complex definitions

extern float clangs (char *, SLU::SuperMatrix *);

extern void cgscon (char *, SuperMatrix *, SuperMatrix *,
float, float *, SuperLUStat_t*, int *);

#ifdef HAVE_AMESOS2_SUPERLU5_API
typedef struct {
int *xsup; /* supernode and column mapping */
Expand Down Expand Up @@ -341,6 +356,11 @@ at_plus_a(

namespace Z { // double-precision complex definitions

extern double zlangs (char *, SLU::SuperMatrix *);

extern void zgscon (char *, SuperMatrix *, SuperMatrix *,
double, double *, SuperLUStat_t*, int *);

#ifdef HAVE_AMESOS2_SUPERLU5_API
typedef struct {
int *xsup; /* supernode and column mapping */
Expand Down Expand Up @@ -466,6 +486,17 @@ namespace Amesos2 {
typedef typename SLU::S::GlobalLU_t GlobalLU_type;
#endif

static float langs(char *norm, SLU::SuperMatrix *A)
{
return SLU::S::slangs(norm, A);
}

static void gscon (char *norm, SLU::SuperMatrix *L, SLU::SuperMatrix *U,
float anorm, float *rcond, SLU::SuperLUStat_t *stat, int *info)
{
SLU::S::sgscon (norm, L, U, anorm, rcond, stat, info);
}

/**
* \brief Binds to the appropriate Superlu solver driver based on data type
*/
Expand Down Expand Up @@ -651,6 +682,17 @@ namespace Amesos2 {
typedef typename SLU::D::GlobalLU_t GlobalLU_type;
#endif

static double langs(char *norm, SLU::SuperMatrix *A)
{
return SLU::D::dlangs(norm, A);
}

static void gscon (char *norm, SLU::SuperMatrix *L, SLU::SuperMatrix *U,
double anorm, double *rcond, SLU::SuperLUStat_t *stat, int *info)
{
SLU::D::dgscon (norm, L, U, anorm, rcond, stat, info);
}

static void gssvx(SLU::superlu_options_t* options, SLU::SuperMatrix* A,
int* perm_c, int* perm_r, int* etree, char* equed, double* R, double* C,
SLU::SuperMatrix* L, SLU::SuperMatrix* U, void* work, int lwork,
Expand Down Expand Up @@ -781,6 +823,17 @@ namespace Amesos2 {
typedef typename SLU::C::GlobalLU_t GlobalLU_type;
#endif

static float langs(char *norm, SLU::SuperMatrix *A)
{
return SLU::C::clangs(norm, A);
}

static void gscon (char *norm, SLU::SuperMatrix *L, SLU::SuperMatrix *U,
float anorm, float *rcond, SLU::SuperLUStat_t *stat, int *info)
{
SLU::C::cgscon (norm, L, U, anorm, rcond, stat, info);
}

static void gssvx(SLU::superlu_options_t* options, SLU::SuperMatrix* A,
int* perm_c, int* perm_r, int* etree, char* equed, float* R, float* C,
SLU::SuperMatrix* L, SLU::SuperMatrix* U, void* work, int lwork,
Expand Down Expand Up @@ -922,6 +975,17 @@ namespace Amesos2 {
typedef typename SLU::Z::GlobalLU_t GlobalLU_type;
#endif

static double langs(char *norm, SLU::SuperMatrix *A)
{
return SLU::Z::zlangs(norm, A);
}

static void gscon (char *norm, SLU::SuperMatrix *L, SLU::SuperMatrix *U,
double anorm, double *rcond, SLU::SuperLUStat_t *stat, int *info)
{
SLU::Z::zgscon (norm, L, U, anorm, rcond, stat, info);
}

static void gssvx(SLU::superlu_options_t* options, SLU::SuperMatrix* A,
int* perm_c, int* perm_r, int* etree, char* equed, double* R, double* C,
SLU::SuperMatrix* L, SLU::SuperMatrix* U, void* work, int lwork,
Expand Down
1 change: 1 addition & 0 deletions packages/amesos2/src/Amesos2_Superlu_decl.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -269,6 +269,7 @@ class Superlu : public SolverCore<Amesos2::Superlu, Matrix, Vector>
char equed;
bool rowequ, colequ; // flags what type of equilibration
// has been performed
magnitude_type anorm, rcond; // condition number estimate

int relax;
int panel_size;
Expand Down
93 changes: 80 additions & 13 deletions packages/amesos2/src/Amesos2_Superlu_def.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -79,17 +79,17 @@ Superlu<Matrix,Vector>::Superlu(
Teuchos::RCP<Vector> X,
Teuchos::RCP<const Vector> B )
: SolverCore<Amesos2::Superlu,Matrix,Vector>(A, X, B)
, is_contiguous_(true) // default is set by params
, use_triangular_solves_(false) // default is set by params
, use_metis_(false)
, symmetrize_metis_(true)
#if defined(KOKKOSKERNELS_ENABLE_SUPERNODAL_SPTRSV) && defined(KOKKOSKERNELS_ENABLE_TPL_SUPERLU)
, sptrsv_invert_diag_(true)
, sptrsv_invert_offdiag_ (false)
, sptrsv_u_in_csr_ (true)
, sptrsv_merge_supernodes_ (false)
, sptrsv_use_spmv_ (false)
#endif
, is_contiguous_(true) // default is set by params
, use_triangular_solves_(false) // default is set by params
, use_metis_(false)
, symmetrize_metis_(true)
{
// ilu_set_default_options is called later in set parameter list if required.
// This is not the ideal way, but the other option to always call
Expand Down Expand Up @@ -422,6 +422,19 @@ Superlu<Matrix,Vector>::numericFactorization_impl()
&(data_.stat), &info);
}

if (data_.options.ConditionNumber == SLU::YES) {
char norm[1];
if (data_.options.Trans == SLU::NOTRANS) {
*(unsigned char *)norm = '1';
} else {
*(unsigned char *)norm = 'I';
}

data_.anorm = function_map::langs(norm, &(data_.A));
function_map::gscon(norm, &(data_.L), &(data_.U),
data_.anorm, &(data_.rcond),
&(data_.stat), &info);
}
}
// Cleanup. AC data will be alloc'd again for next factorization (if at all)
SLU::Destroy_CompCol_Permuted( &(data_.AC) );
Expand All @@ -446,6 +459,24 @@ Superlu<Matrix,Vector>::numericFactorization_impl()
same_symbolic_ = true;

if(use_triangular_solves_) {
#ifdef HAVE_AMESOS2_VERBOSE_DEBUG
#if defined(KOKKOSKERNELS_ENABLE_SUPERNODAL_SPTRSV) && defined(KOKKOSKERNELS_ENABLE_TPL_SUPERLU)
if (this->getComm()->getRank()) {
std::cout << " > Metis : " << (use_metis_ ? "YES" : "NO") << std::endl;
std::cout << " > Equil : " << (data_.options.Equil == SLU::YES ? "YES" : "NO") << std::endl;
std::cout << " > Cond Number : " << (data_.options.ConditionNumber == SLU::YES ? "YES" : "NO") << std::endl;
std::cout << " > Invert diag : " << sptrsv_invert_diag_ << std::endl;
std::cout << " > Invert off-diag : " << sptrsv_invert_offdiag_ << std::endl;
std::cout << " > U in CSR : " << sptrsv_u_in_csr_ << std::endl;
std::cout << " > Merge : " << sptrsv_merge_supernodes_ << std::endl;
std::cout << " > Use SpMV : " << sptrsv_use_spmv_ << std::endl;
}
//std::cout << myRank << " : siize(A) " << (data_.A.nrow) << " x " << (data_.A.ncol) << std::endl;
//std::cout << myRank << " : nnz(A) " << ((SLU::NCformat*)data_.A.Store)->nnz << std::endl;
//std::cout << myRank << " : nnz(L) " << ((SLU::SCformat*)data_.L.Store)->nnz << std::endl;
//std::cout << myRank << " : nnz(U) " << ((SLU::SCformat*)data_.U.Store)->nnz << std::endl;
#endif
#endif
#ifdef HAVE_AMESOS2_TIMERS
Teuchos::RCP< Teuchos::Time > SpTrsvTimer_ = Teuchos::TimeMonitor::getNewCounter ("Time for SpTrsv setup");
Teuchos::TimeMonitor numFactTimer(*SpTrsvTimer_);
Expand Down Expand Up @@ -715,6 +746,9 @@ Superlu<Matrix,Vector>::setParameters_impl(const Teuchos::RCP<Teuchos::Parameter
bool equil = parameterList->get<bool>("Equil", true);
data_.options.Equil = equil ? SLU::YES : SLU::NO;

bool condNum = parameterList->get<bool>("ConditionNumber", false);
data_.options.ConditionNumber = condNum ? SLU::YES : SLU::NO;

bool symmetric_mode = parameterList->get<bool>("SymmetricMode", false);
data_.options.SymmetricMode = symmetric_mode ? SLU::YES : SLU::NO;

Expand Down Expand Up @@ -751,10 +785,10 @@ Superlu<Matrix,Vector>::setParameters_impl(const Teuchos::RCP<Teuchos::Parameter

is_contiguous_ = parameterList->get<bool>("IsContiguous", true);

use_triangular_solves_ = parameterList->get<bool>("Enable_KokkosKernels_TriangularSolves", false);

use_metis_ = parameterList->get<bool>("UseMetis", false);
symmetrize_metis_ = parameterList->get<bool>("SymmetrizeMetis", true);

use_triangular_solves_ = parameterList->get<bool>("Enable_KokkosKernels_TriangularSolves", false);
if(use_triangular_solves_) {
#if defined(KOKKOSKERNELS_ENABLE_SUPERNODAL_SPTRSV) && defined(KOKKOSKERNELS_ENABLE_TPL_SUPERLU)
// specify whether to invert diagonal blocks
Expand Down Expand Up @@ -836,6 +870,7 @@ Superlu<Matrix,Vector>::getValidParameters_impl() const
diag_pivot_thresh_validator); // partial pivoting

pl->set("Equil", true, "Whether to equilibrate the system before solve");
pl->set("ConditionNumber", false, "Whether to approximate condition number");

pl->set("SymmetricMode", false,
"Specifies whether to use the symmetric mode. "
Expand Down Expand Up @@ -896,14 +931,13 @@ Superlu<Matrix,Vector>::getValidParameters_impl() const

pl->set("ILU_Flag", false, "ILU flag: if true, run ILU routines");

pl->set("Enable_KokkosKernels_TriangularSolves", false, "Whether to use triangular solves.");

pl->set("IsContiguous", true, "Whether GIDs contiguous");

// call METIS before SuperLU
pl->set("UseMetis", false, "Whether to call METIS before SuperLU");
pl->set("SymmetrizeMetis", true, "Whether to symmetrize matrix before METIS");

pl->set("Enable_KokkosKernels_TriangularSolves", false, "Whether to use triangular solves.");
#if defined(KOKKOSKERNELS_ENABLE_SUPERNODAL_SPTRSV) && defined(KOKKOSKERNELS_ENABLE_TPL_SUPERLU)
pl->set("SpTRSV_Invert_Diag", true, "specify whether to invert diagonal blocks for supernodal sparse-trianguular solve");
pl->set("SpTRSV_Invert_OffDiag", false, "specify whether to apply diagonal-inversion to off-diagonal blocks for supernodal sparse-trianguular solve");
Expand Down Expand Up @@ -1022,8 +1056,37 @@ Superlu<Matrix,Vector>::triangular_solve_factor()
deep_copy_or_assign_view(device_trsv_C_, data_.C); // will use device to scale
}

bool condition_flag = false;
if (data_.options.ConditionNumber == SLU::YES) {
using STM = Teuchos::ScalarTraits<magnitude_type>;
const magnitude_type eps = STM::eps ();
int n = data_.perm_r.extent(0);

SCformat *Lstore = (SCformat*)(data_.L.Store);
int nsuper = 1 + Lstore->nsuper;
int *nb = Lstore->sup_to_col;
int max_cols = 0;
for (int i = 0; i < nsuper; i++) {
if (nb[i+1] - nb[i] > max_cols) {
max_cols = nb[i+1] - nb[i];
}
}

// when rcond is small, it is ill-conditioned and flag is true
const magnitude_type multiply_fact (10.0); // larger the value, more likely flag is true, and no invert
condition_flag = (((double)max_cols * nsuper) * eps * multiply_fact >= data_.rcond);

#ifdef HAVE_AMESOS2_VERBOSE_DEBUG
std::cout << this->getComm()->getRank()
<< " : anorm = " << data_.anorm << ", rcond = " << data_.rcond << ", n = " << n
<< ", num super cols = " << nsuper << ", max super cols = " << max_cols
<< " -> " << ((double)max_cols * nsuper) * eps / data_.rcond
<< (((double)max_cols * nsuper) * eps * multiply_fact < data_.rcond ? " (okay)" : " (warn)") << std::endl;
#endif
}

// Create handles for U and U^T solves
if (sptrsv_use_spmv_) {
if (sptrsv_use_spmv_ && !condition_flag) {
device_khL_.create_sptrsv_handle(
KokkosSparse::Experimental::SPTRSVAlgorithm::SUPERNODAL_SPMV_DAG, ld_rhs, true);
device_khU_.create_sptrsv_handle(
Expand All @@ -1042,13 +1105,17 @@ Superlu<Matrix,Vector>::triangular_solve_factor()
device_khL_.set_sptrsv_merge_supernodes (sptrsv_merge_supernodes_);
device_khU_.set_sptrsv_merge_supernodes (sptrsv_merge_supernodes_);

// invert only if flag is not on
bool sptrsv_invert_diag = (!condition_flag && sptrsv_invert_diag_);
bool sptrsv_invert_offdiag = (!condition_flag && sptrsv_invert_offdiag_);

// specify whether to invert diagonal blocks
device_khL_.set_sptrsv_invert_diagonal (sptrsv_invert_diag_);
device_khU_.set_sptrsv_invert_diagonal (sptrsv_invert_diag_);
device_khL_.set_sptrsv_invert_diagonal (sptrsv_invert_diag);
device_khU_.set_sptrsv_invert_diagonal (sptrsv_invert_diag);

// specify wheather to apply diagonal-inversion to off-diagonal blocks (optional, default is false)
device_khL_.set_sptrsv_invert_offdiagonal (sptrsv_invert_offdiag_);
device_khU_.set_sptrsv_invert_offdiagonal (sptrsv_invert_offdiag_);
device_khL_.set_sptrsv_invert_offdiagonal (sptrsv_invert_offdiag);
device_khU_.set_sptrsv_invert_offdiagonal (sptrsv_invert_offdiag);

// set etree
device_khL_.set_sptrsv_etree(data_.parents.data());
Expand Down
4 changes: 4 additions & 0 deletions packages/ifpack/src/Ifpack_Hypre.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -798,7 +798,9 @@ int Ifpack_Hypre::Compute(){
hypre_ParAMGData *amg_data = (hypre_ParAMGData*) Preconditioner_;
hypre_ParCSRMatrix **A_array = hypre_ParAMGDataAArray(amg_data);
hypre_ParCSRMatrix **P_array = hypre_ParAMGDataPArray(amg_data);
#if 0
HYPRE_Int **CF_marker_array = hypre_ParAMGDataCFMarkerArray(amg_data);
#endif
HYPRE_Int num_levels = hypre_ParAMGDataNumLevels(amg_data);

char ofs[80];
Expand All @@ -811,6 +813,7 @@ int Ifpack_Hypre::Compute(){
sprintf(ofs,"P_matrix.bmg.%d.dat",k);
HYPRE_ParCSRMatrixPrint(P_array[k], ofs);

#if 0
// CF
// Note: Hypre outputs "-1" for F Points and "1" for C Points
HYPRE_Int local_size = hypre_CSRMatrixNumRows(hypre_ParCSRMatrixDiag(A_array[k]));
Expand All @@ -821,6 +824,7 @@ int Ifpack_Hypre::Compute(){
for(int i=0; i<local_size; i++)
fprintf(f,"%d\n",(int)CF_marker_array[k][i]);
fclose(f);
#endif
}

}
Expand Down
Loading

0 comments on commit 5a67b70

Please sign in to comment.