Skip to content

Commit

Permalink
Merge pull request #9824 from iyamazaki/amesos2-superlu-dist
Browse files Browse the repository at this point in the history 
Amesos2::SuperLU_dist: add option to apply Equil
  • Loading branch information
@ndellingwood
ndellingwood authored Oct 26, 2021
2 parents 3bae707 + 98bc0c8 commit 0c84fa1
Show file tree
Hide file tree
Showing 2 changed files with 64 additions and 41 deletions.
4 changes: 2 additions & 2 deletions packages/amesos2/src/Amesos2_Superludist_FunctionMap.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -532,14 +532,14 @@ namespace Amesos2 {
}

static void gsequ_loc(SLUD::SuperMatrix* A, double* r, double* c,
double* rowcnd, double* colcnd, double* amax, int* info,
double* rowcnd, double* colcnd, double* amax, SLUD::int_t* info,
SLUD::gridinfo_t* grid)
{
SLUD::Z::pzgsequ(A, r, c, rowcnd, colcnd, amax, info, grid);
}

static void gsequ(SLUD::SuperMatrix* A, double* r, double* c,
double* rowcnd, double* colcnd, double* amax, int* info)
double* rowcnd, double* colcnd, double* amax, SLUD::int_t* info)
{
SLUD::Z::zgsequ_dist(A, r, c, rowcnd, colcnd, amax, info);
}
Expand Down
101 changes: 62 additions & 39 deletions packages/amesos2/src/Amesos2_Superludist_def.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -355,17 +355,17 @@ namespace Amesos2 {
free( data_.fstVtxSep );
#endif
}
float info = 0.0;
{
#ifdef HAVE_AMESOS2_TIMERS
Teuchos::TimeMonitor preOrderTime( this->timers_.preOrderTime_ );
Teuchos::TimeMonitor preOrderTime( this->timers_.preOrderTime_ );
#endif

float info = 0.0;
info = SLUD::get_perm_c_parmetis( &(data_.A),
data_.perm_r.getRawPtr(), data_.perm_c.getRawPtr(),
data_.grid.nprow * data_.grid.npcol, data_.domains,
&(data_.sizes), &(data_.fstVtxSep),
&(data_.grid), &(data_.symb_comm) );

info = SLUD::get_perm_c_parmetis( &(data_.A),
data_.perm_r.getRawPtr(), data_.perm_c.getRawPtr(),
data_.grid.nprow * data_.grid.npcol, data_.domains,
&(data_.sizes), &(data_.fstVtxSep),
&(data_.grid), &(data_.symb_comm) );
}
TEUCHOS_TEST_FOR_EXCEPTION( info > 0.0,
std::runtime_error,
"SuperLU_DIST pre-ordering ran out of memory after allocating "
Expand All @@ -389,18 +389,18 @@ namespace Amesos2 {

if( in_grid_ ){

float info = 0.0;
{
#ifdef HAVE_AMESOS2_TIMERS
Teuchos::TimeMonitor symFactTime( this->timers_.symFactTime_ );
Teuchos::TimeMonitor symFactTime( this->timers_.symFactTime_ );
#endif

float info = 0.0;
info = SLUD::symbfact_dist((data_.grid.nprow) * (data_.grid.npcol),
data_.domains, &(data_.A), data_.perm_c.getRawPtr(),
data_.perm_r.getRawPtr(), data_.sizes,
data_.fstVtxSep, &(data_.pslu_freeable),
&(data_.grid.comm), &(data_.symb_comm),
&(data_.mem_usage));

info = SLUD::symbfact_dist((data_.grid.nprow) * (data_.grid.npcol),
data_.domains, &(data_.A), data_.perm_c.getRawPtr(),
data_.perm_r.getRawPtr(), data_.sizes,
data_.fstVtxSep, &(data_.pslu_freeable),
&(data_.grid.comm), &(data_.symb_comm),
&(data_.mem_usage));
}
TEUCHOS_TEST_FOR_EXCEPTION( info > 0.0,
std::runtime_error,
"SuperLU_DIST symbolic factorization ran out of memory after"
Expand All @@ -420,17 +420,25 @@ namespace Amesos2 {

// loadA_impl(); // Refresh the matrix values

// if( data_.options.Equil == SLUD::YES ){
// // Apply the scalings computed in preOrdering
// function_map::laqgs(&(data_.A), data_.R.getRawPtr(),
// data_.C.getRawPtr(), data_.rowcnd, data_.colcnd,
// data_.amax, &(data_.equed));
if( in_grid_ ) {
if( data_.options.Equil == SLUD::YES ) {
SLUD::int_t info = 0;

// data_.rowequ = (data_.equed == SLUD::ROW) || (data_.equed == SLUD::BOTH);
// data_.colequ = (data_.equed == SLUD::COL) || (data_.equed == SLUD::BOTH);
// }
// Compute scaling
data_.R.resize(this->globalNumRows_);
data_.C.resize(this->globalNumCols_);
function_map::gsequ_loc(&(data_.A), data_.R.getRawPtr(), data_.C.getRawPtr(),
&(data_.rowcnd), &(data_.colcnd), &(data_.amax), &info, &(data_.grid));

// Apply the scalings
function_map::laqgs_loc(&(data_.A), data_.R.getRawPtr(), data_.C.getRawPtr(),
data_.rowcnd, data_.colcnd, data_.amax,
&(data_.equed));

data_.rowequ = (data_.equed == SLUD::ROW) || (data_.equed == SLUD::BOTH);
data_.colequ = (data_.equed == SLUD::COL) || (data_.equed == SLUD::BOTH);
}

if( in_grid_ ){
// Apply the column ordering, so that AC is the column-permuted A, and compute etree
size_t nnz_loc = ((SLUD::NRformat_loc*)data_.A.Store)->nnz_loc;
for( size_t i = 0; i < nnz_loc; ++i ) colind_[i] = data_.perm_c[colind_[i]];
Expand Down Expand Up @@ -462,7 +470,6 @@ namespace Amesos2 {
#ifdef HAVE_AMESOS2_TIMERS
Teuchos::TimeMonitor numFactTimer(this->timers_.numFactTime_);
#endif

function_map::gstrf(&(data_.options), this->globalNumRows_,
this->globalNumCols_, anorm, &(data_.lu),
&(data_.grid), &(data_.stat), &info);
Expand Down Expand Up @@ -511,7 +518,6 @@ namespace Amesos2 {
#ifdef HAVE_AMESOS2_TIMERS
Teuchos::TimeMonitor convTimer(this->timers_.vecConvTime_);
#endif

{
// The input dense matrix for B should be distributed in the
// same manner as the superlu_dist matrix. That is, if a
Expand All @@ -522,7 +528,6 @@ namespace Amesos2 {
#ifdef HAVE_AMESOS2_TIMERS
Teuchos::TimeMonitor redistTimer(this->timers_.vecRedistTime_);
#endif

// get grid-distributed mv data. The multivector data will be
// distributed across the processes in the SuperLU_DIST grid.
typedef Util::get_1d_copy_helper<MultiVecAdapter<Vector>,slu_type> copy_helper;
Expand Down Expand Up @@ -563,12 +568,22 @@ namespace Amesos2 {
same_solve_struct_ = true;
}

// Apply row-scaling if requested
if (data_.options.Equil == SLUD::YES && data_.rowequ) {
SLUD::int_t ld = as<SLUD::int_t>(local_len_rhs);
for(global_size_type j = 0; j < nrhs; ++j) {
for(size_t i = 0; i < local_len_rhs; ++i) {
bvals_[i + j*ld] *= data_.R[first_global_row_b + i];
}
}
}

// Solve
int ierr = 0; // returned error code
{
#ifdef HAVE_AMESOS2_TIMERS
Teuchos::TimeMonitor solveTimer(this->timers_.solveTime_);
#endif

function_map::gstrs(as<SLUD::int_t>(this->globalNumRows_), &(data_.lu),
&(data_.scale_perm), &(data_.grid), bvals_.getRawPtr(),
as<SLUD::int_t>(local_len_rhs), as<SLUD::int_t>(first_global_row_b),
Expand Down Expand Up @@ -606,14 +621,23 @@ namespace Amesos2 {
as<int>(nrhs),
&(data_.grid));
}

// Apply col-scaling if requested
if (data_.options.Equil == SLUD::YES && data_.colequ) {
SLUD::int_t ld = as<SLUD::int_t>(local_len_rhs);
for(global_size_type j = 0; j < nrhs; ++j) {
for(size_t i = 0; i < local_len_rhs; ++i) {
xvals_[i + j*ld] *= data_.C[first_global_row_b + i];
}
}
}
}

/* Update X's global values */
{
#ifdef HAVE_AMESOS2_TIMERS
Teuchos::TimeMonitor redistTimer(this->timers_.vecRedistTime_);
#endif

typedef Util::put_1d_data_helper<MultiVecAdapter<Vector>,slu_type> put_helper;
put_helper::do_put(X,
xvals_(),
Expand Down Expand Up @@ -672,10 +696,9 @@ namespace Amesos2 {

data_.options.Trans = SLUD::NOTRANS; // should always be set this way;

// TODO: Uncomment when supported
// bool equil = parameterList->get<bool>("Equil", true);
// data_.options.Equil = equil ? SLUD::YES : SLUD::NO;
data_.options.Equil = SLUD::NO;
// Equilbration option
bool equil = parameterList->get<bool>("Equil", false);
data_.options.Equil = equil ? SLUD::YES : SLUD::NO;

if( parameterList->isParameter("ColPerm") ){
RCP<const ParameterEntryValidator> colperm_validator = valid_params->getEntry("ColPerm").validator();
Expand Down Expand Up @@ -743,8 +766,8 @@ namespace Amesos2 {
tuple<SLUD::trans_t>(SLUD::NOTRANS),
pl.getRawPtr());

// TODO: uncomment when supported
// pl->set("Equil", false, "Whether to equilibrate the system before solve");
// Equillbration
pl->set("Equil", false, "Whether to equilibrate the system before solve");

// TODO: uncomment when supported
// setStringToIntegralParameter<SLUD::IterRefine_t>("IterRefine", "NOREFINE",
Expand Down

0 comments on commit 0c84fa1

Please sign in to comment.