Skip to content

Commit

Permalink
fix local and non-local sparsity for local interfaces
Browse files Browse the repository at this point in the history
  • Loading branch information
greole committed Nov 7, 2024
1 parent 504636f commit fab7b70
Show file tree
Hide file tree
Showing 6 changed files with 172 additions and 73 deletions.
40 changes: 26 additions & 14 deletions include/OGL/MatrixWrapper/HostMatrix.H
Original file line number Diff line number Diff line change
Expand Up @@ -178,6 +178,27 @@ private:
[[nodiscard]] std::vector<interface_locality> collect_cells_on_interfaces(
const lduInterfaceFieldPtrsList &interfaces) const;

/** Based on OpenFOAMs ldu matrix and interfaces this function computes
** the local sparsity pattern. Here, the sparsity pattern
** contains coefficients that resides on the owning rank. The coefficient
** can include local interfaces.
**/
[[nodiscard]] std::tuple<std::vector<label>, std::vector<label>,
std::vector<label>, std::vector<gko::span>>
compute_local_sparsity(std::shared_ptr<const gko::Executor> exec) const;

/** Based on OpenFOAMs interfaces this function computes
** the interface sparsity pattern.
** Here, the sparsity pattern
** only contains all coefficient that reside on an interface (local and
*non-local).
**/
[[nodiscard]] std::tuple<std::vector<label>, std::vector<label>,
std::vector<label>, std::vector<label>,
std::vector<gko::span>>
compute_interface_sparsity(std::shared_ptr<const gko::Executor> exec) const;


public:
// segregated matrix wrapper constructor
HostMatrixWrapper(const ExecutorHandler &exec, const objectRegistry &db,
Expand All @@ -201,13 +222,12 @@ public:
const dictionary &solverControls, const word &fieldName,
label verbose);

/** Based on OpenFOAMs ldu matrix format this function computes two
** consecutive index arrays (local_sparsisty_.row_idxs and col_idxs) in row
** major ordering and the permutation index (local_sparsity_.ldu_mapping),
** which are required to create to a ginkgo matrix
/** Based on OpenFOAMs ldu matrix and interfaces this function computes
** the local and non-local sparsity pattern.
**/
[[nodiscard]] std::shared_ptr<SparsityPattern> compute_local_sparsity(
std::shared_ptr<const gko::Executor> exec) const;
[[nodiscard]] std::pair<std::shared_ptr<SparsityPattern>,
std::shared_ptr<SparsityPattern>>
compute_sparsity_patterns(std::shared_ptr<const gko::Executor> exec) const;

/** Copies the LDU matrix coefficients to local_coeffs without changing or
** reinstantiating the sparsity pattern.
Expand All @@ -224,14 +244,6 @@ public:
std::shared_ptr<const SparsityPattern> sparsity,
gko::array<scalar> &target_coeffs) const;

/** Based on OpenFOAMs interfaces this function computes two
** consecutive index arrays (non_local_sparsisty_.row_idxs and col_idxs) in
** row *major ordering and the permutation index
** (non_local_sparsity_.ldu_mapping), which are required to create a ginkgo
** matrix
**/
[[nodiscard]] std::shared_ptr<SparsityPattern> compute_non_local_sparsity(
std::shared_ptr<const gko::Executor> exec) const;

/** Iterates all interfaces and counts the number of unique neighbour
** processors and number of interfaces in total for this processor
Expand Down
3 changes: 1 addition & 2 deletions src/MatrixWrapper/Distributed.C
Original file line number Diff line number Diff line change
Expand Up @@ -614,8 +614,7 @@ std::shared_ptr<RepartDistMatrix> create_impl(
auto exec = exec_handler.get_ref_exec();
auto comm = *exec_handler.get_communicator().get();

auto local_sparsity = host_A->compute_local_sparsity(exec);
auto non_local_sparsity = host_A->compute_non_local_sparsity(exec);
auto [local_sparsity,non_local_sparsity] = host_A->compute_sparsity_patterns(exec);

auto src_comm_pattern = host_A->create_communication_pattern();
if (non_local_sparsity->spans.size() !=
Expand Down
109 changes: 96 additions & 13 deletions src/MatrixWrapper/HostMatrix.C
Original file line number Diff line number Diff line change
Expand Up @@ -362,6 +362,7 @@ std::vector<interface_locality> HostMatrixWrapper::collect_cells_on_interfaces(
interface_ctr++;
}
}

if (isA<cyclicFvPatch>(iface->interface())) {
const cyclicFvPatch &patch =
refCast<const cyclicFvPatch>(iface->interface());
Expand All @@ -373,13 +374,11 @@ std::vector<interface_locality> HostMatrixWrapper::collect_cells_on_interfaces(
const label neighbPatchId = patch.nbrPatchID();
#endif
const labelUList &cols = addr_.patchAddr(neighbPatchId);

for (label cellI = 0; cellI < interface_size; cellI++) {
interface_idxs.emplace_back(interface_id, cols[cellI],
face_cells[cellI], rank, local_ctr,
interface_ctr);
}

local_ctr++;
}
interface_id++;
Expand All @@ -390,14 +389,17 @@ std::vector<interface_locality> HostMatrixWrapper::collect_cells_on_interfaces(
return interface_idxs;
}

std::shared_ptr<SparsityPattern> HostMatrixWrapper::compute_non_local_sparsity(
std::tuple<std::vector<label>, std::vector<label>, std::vector<label>,
std::vector<label>, std::vector<gko::span>>
HostMatrixWrapper::compute_interface_sparsity(
std::shared_ptr<const gko::Executor> exec) const
{
auto interface_loc_vec = collect_cells_on_interfaces(interfaces_);
label total_interface_nnz = interface_loc_vec.back().non_local_nnz_ctr + 1;
std::vector<label> rows_vec(total_interface_nnz);
std::vector<label> cols_vec(total_interface_nnz);
std::vector<label> mapping_vec(total_interface_nnz);
std::vector<label> ranks{};
std::vector<gko::span> spans{};

auto rows = rows_vec.data();
Expand All @@ -406,6 +408,7 @@ std::shared_ptr<SparsityPattern> HostMatrixWrapper::compute_non_local_sparsity(
label prev_interface_ctr{0};
label end{0};
label start{0};
label rank_;

size_t element_ctr = 0;

Expand All @@ -415,34 +418,112 @@ std::shared_ptr<SparsityPattern> HostMatrixWrapper::compute_non_local_sparsity(
cols[element_ctr] = col;
permute[element_ctr] = element_ctr;

// a new interface started been reached
// a new interface start has been reached
if (interface_idx > prev_interface_ctr) {
ranks.push_back(rank_);
end = element_ctr;
spans.emplace_back(start, element_ctr);
start = end;
prev_interface_ctr = interface_idx;
}
rank_ = rank;
element_ctr++;
}
spans.emplace_back(start, element_ctr);
ranks.push_back(rank_);
return {rows_vec, cols_vec, mapping_vec, ranks, spans};
}


std::pair<std::shared_ptr<SparsityPattern>, std::shared_ptr<SparsityPattern>>
HostMatrixWrapper::compute_sparsity_patterns(
std::shared_ptr<const gko::Executor> exec) const
{
auto rank = get_exec_handler().get_rank();
auto [local_rows, local_cols, local_mapping, local_spans] =
compute_local_sparsity(exec);
auto [non_local_rows, non_local_cols, non_local_mapping, non_local_ranks,
non_local_spans] = compute_interface_sparsity(exec);

// move all local interfaces to local rows and cols
std::vector<size_t> erase_non_local{};
std::vector<size_t> keep_non_local{};
for (int interface_ctr = 0; interface_ctr < non_local_spans.size();
interface_ctr++) {
if (non_local_ranks[interface_ctr] == rank) {
auto [begin, end] = non_local_spans[interface_ctr];
size_t start = local_rows.size();
local_rows.insert(local_rows.end(), non_local_rows.data() + begin,
non_local_rows.data() + end);
local_cols.insert(local_cols.end(), non_local_cols.data() + begin,
non_local_cols.data() + end);
local_mapping.insert(local_mapping.end(),
non_local_mapping.data() + begin,
non_local_mapping.data() + end);
local_spans.emplace_back(start, local_rows.size());
erase_non_local.push_back(interface_ctr);
} else {
keep_non_local.push_back(interface_ctr);
}
}

// for (size_t erase_ctr=0; erase_ctr<erase_non_local.size();erase_ctr++) {
// auto erase_id {erase_non_local[erase_non_local.size() - erase_ctr
// - 1]}; auto [begin, end] = non_local_spans[erase_id];
// non_local_rows.erase(non_local_rows.begin() + begin,
// non_local_rows.begin() + end);
// non_local_cols.erase(non_local_cols.begin() + begin,
// non_local_cols.begin() + end);
// non_local_mapping.erase(non_local_mapping.begin() + begin,
// non_local_mapping.begin() + end);
// }

std::vector<label> non_local_rows_copy, non_local_cols_copy,
non_local_mapping_copy;
std::vector<gko::span> non_local_spans_copy{};
size_t begin = 0;
for (auto keep : keep_non_local) {
auto span = non_local_spans[keep];
size_t length{span.end - span.begin};
non_local_spans_copy.emplace_back(begin, begin + length);
non_local_rows_copy.insert(non_local_rows_copy.end(),
non_local_rows.data() + span.begin,
non_local_rows.data() + span.end);

non_local_cols_copy.insert(non_local_cols_copy.end(),
non_local_cols.data() + span.begin,
non_local_cols.data() + span.end);

non_local_mapping_copy.insert(non_local_mapping_copy.end(),
non_local_mapping.data() + span.begin,
non_local_mapping.data() + span.end);

begin += length;
}

gko::dim<2> dim{static_cast<gko::size_type>(nrows_),
static_cast<gko::size_type>(non_local_matrix_nnz_)};
return std::make_shared<SparsityPattern>(exec->get_master(), dim, rows_vec,
cols_vec, mapping_vec, spans);
auto local_sparsity = std::make_shared<SparsityPattern>(
exec->get_master(), get_size(), local_rows, local_cols, local_mapping,
local_spans);
auto non_local_sparsity = std::make_shared<SparsityPattern>(
exec->get_master(), gko::dim<2>{nrows_, non_local_rows_copy.size()},
non_local_rows_copy, non_local_cols_copy, non_local_mapping_copy,
non_local_spans_copy);
return {local_sparsity, non_local_sparsity};
}


std::shared_ptr<SparsityPattern> HostMatrixWrapper::compute_local_sparsity(
std::tuple<std::vector<label>, std::vector<label>, std::vector<label>,
std::vector<gko::span>>
HostMatrixWrapper::compute_local_sparsity(
std::shared_ptr<const gko::Executor> exec) const
{
LOG_1(verbose_, "start init host sparsity pattern")

std::vector<label> rows_vec(local_matrix_w_interfaces_nnz_);
std::vector<label> cols_vec(local_matrix_w_interfaces_nnz_);
std::vector<label> mapping_vec(local_matrix_w_interfaces_nnz_);
std::vector<gko::span> spans{gko::span{
0, static_cast<gko::size_type>(local_matrix_w_interfaces_nnz_)}};
std::vector<gko::span> spans{
gko::span{0, static_cast<gko::size_type>(local_matrix_nnz_)}};

auto rows = rows_vec.data();
auto cols = cols_vec.data();
Expand Down Expand Up @@ -477,8 +558,10 @@ std::shared_ptr<SparsityPattern> HostMatrixWrapper::compute_local_sparsity(
// i_j,k j=interface index and k cell index on the interface

LOG_1(verbose_, "done init host sparsity pattern")
return std::make_shared<SparsityPattern>(
exec->get_master(), get_size(), rows_vec, cols_vec, mapping_vec, spans);
// return std::make_shared<SparsityPattern>(
// exec->get_master(), get_size(), rows_vec, cols_vec, mapping_vec,
// spans);
return {rows_vec, cols_vec, mapping_vec, spans};
}

} // namespace Foam
2 changes: 2 additions & 0 deletions src/Repartitioner.C
Original file line number Diff line number Diff line change
Expand Up @@ -362,6 +362,8 @@ std::vector<InterfaceLocality> Repartitioner::build_non_local_interfaces(
auto [begin, end] = non_local_spans[i];
bool local = reparts_to_local(exec_handler, comm_target_ids[i]);

// TODO this seems to depend on same copying functionality like HostMatrix.C to
// copy local interfaces over.
if (local) {
local_ctr++; // local interface starts counting at 1;
gko::size_type rows_start = local_rows.size();
Expand Down
8 changes: 4 additions & 4 deletions unitTests/MatrixWrapper/HostMatrix_l2d.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -183,8 +183,8 @@ TEST(HostMatrix, canGenerateLocalSparsityPattern)
auto hostMatrix = ((HostMatrixEnvironment *)global_env)->hostMatrix;
auto exec = ((HostMatrixEnvironment *)global_env)->exec;

auto localSparsity =
hostMatrix->compute_local_sparsity(exec->get_device_exec());
auto [localSparsity, nonLocalSparsity] =
hostMatrix->compute_sparsity_patterns(exec->get_device_exec());
std::vector<label> rows_expected({0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 3,
3, 3, 3, 4, 4, 4, 4, 4, 5, 5, 5,
5, 6, 6, 6, 7, 7, 7, 7, 8, 8, 8});
Expand Down Expand Up @@ -245,8 +245,8 @@ TEST(HostMatrix, canGenerateNonLocalSparsityPattern)
auto exec = ((HostMatrixEnvironment *)global_env)->exec;
auto comm = exec->get_gko_mpi_device_comm();

auto nonLocalSparsity =
hostMatrix->compute_non_local_sparsity(exec->get_device_exec());
auto [localSparsity, nonLocalSparsity] =
hostMatrix->compute_sparsity_patterns(exec->get_device_exec());

// corresponds to cell ids
std::vector<std::vector<label>> rows_expected({{2, 5, 8, 6, 7, 8},
Expand Down
83 changes: 43 additions & 40 deletions unitTests/MatrixWrapper/HostMatrix_p2d.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -250,46 +250,49 @@ TEST(HostMatrix, canCreateCommunicationPattern)
// EXPECT_EQ(mapping_expected, mapping_res);
// }

// TEST(HostMatrix, canGenerateNonLocalSparsityPattern)
// {
// auto hostMatrix = ((HostMatrixEnvironment *)global_env)->hostMatrix;
// auto exec = ((HostMatrixEnvironment *)global_env)->exec;
// auto comm = exec->get_gko_mpi_device_comm();

// auto nonLocalSparsity =
// hostMatrix->compute_non_local_sparsity(exec->get_device_exec());

// // corresponds to cell ids
// std::vector<std::vector<label>> rows_expected({{2, 5, 8, 6, 7, 8},
// {0, 3, 6, 6, 7, 8},
// {0, 1, 2, 2, 5, 8},
// {0, 1, 2, 0, 3, 6}});
// // cols expected
// std::vector<std::vector<label>> cols_expected({{0, 3, 6, 0, 1, 2},
// {2, 5, 8, 0, 1, 2},
// {6, 7, 8, 0, 3, 6},
// {6, 7, 8, 2, 5, 8}});

// std::vector<label> mapping_expected({0, 1, 2, 3, 4, 5});

// // we dont test the cols expected for now,
// // as they are in compressed format
// EXPECT_EQ(nonLocalSparsity->num_nnz, 6);
// EXPECT_EQ(nonLocalSparsity->spans.size(), 2);

// EXPECT_EQ(nonLocalSparsity->spans[0].begin, 0);
// EXPECT_EQ(nonLocalSparsity->spans[0].end, 3);
// EXPECT_EQ(nonLocalSparsity->spans[1].begin, 3);
// EXPECT_EQ(nonLocalSparsity->spans[1].end, 6);

// auto res_size{nonLocalSparsity->row_idxs.get_size()};
// auto rows_res = convert_to_vector(nonLocalSparsity->row_idxs);
// auto cols_res = convert_to_vector(nonLocalSparsity->col_idxs);
// auto mapping_res = convert_to_vector(nonLocalSparsity->ldu_mapping);
// EXPECT_EQ(rows_expected[comm->rank()], rows_res);
// EXPECT_EQ(cols_expected[comm->rank()], cols_res);
// EXPECT_EQ(mapping_expected, mapping_res);
// }
TEST(HostMatrix, canGenerateNonLocalSparsityPattern)
{
auto hostMatrix = ((HostMatrixEnvironment *)global_env)->hostMatrix;
auto exec = ((HostMatrixEnvironment *)global_env)->exec;
auto comm = exec->get_gko_mpi_device_comm();
auto rank = exec->get_rank();

auto [localSparsity, nonLocalSparsity] =
hostMatrix->compute_sparsity_patterns(exec->get_device_exec());

// corresponds to cell ids
std::vector<std::vector<label>> rows_expected({{2, 5, 8, 6, 7, 8},
{0, 3, 6, 6, 7, 8},
{0, 1, 2, 2, 5, 8},
{0, 1, 2, 0, 3, 6}});
// cols expected
std::vector<std::vector<label>> cols_expected({{0, 3, 6, 0, 1, 2},
{2, 5, 8, 0, 1, 2},
{6, 7, 8, 0, 3, 6},
{6, 7, 8, 2, 5, 8}});

std::vector<label> mapping_expected({0, 1, 2, 3, 4, 5});

// we dont test the cols expected for now,
// as they are in compressed format
std::vector<label> exp_send_idx_size{8, 16, 16, 8};
EXPECT_EQ(nonLocalSparsity->num_nnz, exp_send_idx_size[rank]);
std::vector<label> exp_spans_size{1, 2, 2, 1};
EXPECT_EQ(nonLocalSparsity->spans.size(), exp_spans_size[rank]);

// EXPECT_EQ(nonLocalSparsity->spans[0].begin, 0);
// EXPECT_EQ(nonLocalSparsity->spans[0].end, 3);
// EXPECT_EQ(nonLocalSparsity->spans[1].begin, 3);
// EXPECT_EQ(nonLocalSparsity->spans[1].end, 6);

// auto res_size{nonLocalSparsity->row_idxs.get_size()};
// auto rows_res = convert_to_vector(nonLocalSparsity->row_idxs);
// auto cols_res = convert_to_vector(nonLocalSparsity->col_idxs);
// auto mapping_res = convert_to_vector(nonLocalSparsity->ldu_mapping);
// EXPECT_EQ(rows_expected[comm->rank()], rows_res);
// EXPECT_EQ(cols_expected[comm->rank()], cols_res);
// EXPECT_EQ(mapping_expected, mapping_res);
}


int main(int argc, char *argv[])
Expand Down

0 comments on commit fab7b70

Please sign in to comment.