Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Amesos2 : matrix-gather when colMap does not match rowMap #13790

Merged
merged 4 commits into from
Feb 7, 2025
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion packages/amesos2/src/Amesos2_KLU2_def.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -453,7 +453,7 @@ KLU2<Matrix,Vector>::loadA_impl(EPhase current_phase)
Kokkos::resize(host_col_ptr_view_, this->globalNumRows_ + 1);
}

local_ordinal_type nnz_ret = 0;
local_ordinal_type nnz_ret = -1;
bool gather_supported = (this->matrixA_->getComm()->getSize() > 1 && (std::is_same<scalar_type, float>::value || std::is_same<scalar_type, double>::value));
{
#ifdef HAVE_AMESOS2_TIMERS
Expand Down
2 changes: 1 addition & 1 deletion packages/amesos2/src/Amesos2_ShyLUBasker_def.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -587,7 +587,7 @@ ShyLUBasker<Matrix,Vector>::loadA_impl(EPhase current_phase)
Kokkos::resize(colptr_view_, this->globalNumCols_ + 1); //this will be wrong for case of gapped col ids, e.g. 0,2,4,9; num_cols = 10 ([0,10)) but num GIDs = 4...
}

local_ordinal_type nnz_ret = 0;
local_ordinal_type nnz_ret = -1;
bool gather_supported = (this->matrixA_->getComm()->getSize() > 1 && (std::is_same<scalar_type, float>::value || std::is_same<scalar_type, double>::value));
{
#ifdef HAVE_AMESOS2_TIMERS
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -98,16 +98,18 @@ namespace Amesos2 {
const EPhase current_phase) const
{
typedef Tpetra::Map< local_ordinal_t, global_ordinal_t, node_t> contiguous_map_type;
auto rowMap = this->mat_->getRowMap();
auto colMap = this->mat_->getColMap();
auto rowComm = rowMap->getComm();
auto colComm = colMap->getComm();

using Teuchos::RCP;
using Teuchos::rcp;
#ifdef HAVE_AMESOS2_TIMERS
auto reindexTimer = Teuchos::TimeMonitor::getNewTimer("Time to re-index matrix gids");
Teuchos::TimeMonitor ReindexTimer(*reindexTimer);
#endif

auto rowMap = this->mat_->getRowMap();
auto colMap = this->mat_->getColMap();
auto rowComm = rowMap->getComm();
auto colComm = colMap->getComm();

GlobalOrdinal indexBase = rowMap->getIndexBase();
GlobalOrdinal numDoFs = this->mat_->getGlobalNumRows();
LocalOrdinal nRows = this->mat_->getLocalNumRows();
Expand All @@ -120,12 +122,17 @@ namespace Amesos2 {
global_ordinal_t frow = tmpMap->getMinGlobalIndex();

// Create new GID list for RowMap
Kokkos::View<global_ordinal_t*, HostExecSpaceType> rowIndexList ("indexList", nRows);
Kokkos::View<global_ordinal_t*, HostExecSpaceType> rowIndexList ("rowIndexList", nRows);
for (local_ordinal_t k = 0; k < nRows; k++) {
rowIndexList(k) = frow+k; // based on index-base of rowMap
}
// Create new GID list for ColMap
Kokkos::View<global_ordinal_t*, HostExecSpaceType> colIndexList ("indexList", nCols);
Kokkos::View<global_ordinal_t*, HostExecSpaceType> colIndexList ("colIndexList", nCols);
// initialize to catch col GIDs that are not in row GIDs
// they will be all assigned to (n+1)th columns
for (local_ordinal_t k = 0; k < nCols; k++) {
colIndexList(k) = numDoFs+indexBase;
}
typedef Tpetra::MultiVector<global_ordinal_t,
local_ordinal_t,
global_ordinal_t,
Expand Down Expand Up @@ -307,9 +314,11 @@ namespace Amesos2 {
recvDispls(p+1) = recvDispls(p) + recvCounts(p);
}
}
// -- convert to global colids & convert to 0-base
// -- convert to global colids & ** convert to base-zero **
KV_GO lclColind_ ("localColind_", lclColind.extent(0));
for (int i = 0; i < int(lclColind.extent(0)); i++) lclColind_(i) = (colMap->getGlobalElement((lclColind(i))) - colIndexBase);
for (int i = 0; i < int(lclColind.extent(0)); i++) {
lclColind_(i) = (colMap->getGlobalElement((lclColind(i))) - colIndexBase);
}
if (column_major || need_to_perm) {
Kokkos::resize(indices_t, indices.extent(0));
Teuchos::gatherv<int, LocalOrdinal> (lclColind_.data(), lclColind_.extent(0), indices_t.data(),
Expand All @@ -326,6 +335,7 @@ namespace Amesos2 {
Teuchos::RCP< Teuchos::Time > gatherTime = Teuchos::TimeMonitor::getNewCounter ("Amesos2::gather(transpose index)");
Teuchos::TimeMonitor GatherTimer(*gatherTime);
#endif
// (note: column idexes are now in base-0)
if (column_major) {
// Map to transpose
Kokkos::resize(transpose_map, ret);
Expand All @@ -347,9 +357,14 @@ namespace Amesos2 {
int i = perm_l2g(row);
for (int k=pointers_t(i); k<pointers_t(i+1); k++) {
int col = indices_t(k);
transpose_map(k) = pointers(1+col);
indices(pointers(1+col)) = row;
pointers(1+col) ++;
if (col < nRows) {
transpose_map(k) = pointers(1+col);
indices(pointers(1+col)) = row;
pointers(1+col) ++;
} else {
// extra columns
transpose_map(k) = -1;
}
}
}
} else if (need_to_perm) {
Expand All @@ -366,9 +381,13 @@ namespace Amesos2 {
int row = perm_g2l(i);
for (int k=pointers_t(i); k<pointers_t(i+1); k++) {
int col = indices_t(k);
transpose_map(k) = pointers(1+row);
indices(pointers(1+row)) = col;
pointers(1+row) ++;
if (col < nRows) {
transpose_map(k) = pointers(1+row);
indices(pointers(1+row)) = col;
pointers(1+row) ++;
} else {
transpose_map(k) = -1;
}
}
}
} else {
Expand Down Expand Up @@ -409,7 +428,9 @@ namespace Amesos2 {
#endif
if (transpose_map.extent(0) > 0) {
for (int k=0; k<ret; k++) {
nzvals(transpose_map(k)) = nzvals_t(k);
if (transpose_map(k) >= 0) {
nzvals(transpose_map(k)) = nzvals_t(k);
}
}
}
}
Expand Down
Loading