Skip to content

Commit

Permalink
Use OpenMP in MatchingMeshPartitionerLonLatPolygon
Browse files Browse the repository at this point in the history
  • Loading branch information
wdeconinck committed Nov 9, 2022
1 parent 379de3d commit d30709c
Show file tree
Hide file tree
Showing 2 changed files with 70 additions and 33 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,7 @@
#include "atlas/mesh/Nodes.h"
#include "atlas/parallel/mpi/mpi.h"
#include "atlas/parallel/omp/fill.h"
#include "atlas/parallel/omp/omp.h"
#include "atlas/runtime/Exception.h"
#include "atlas/runtime/Log.h"
#include "atlas/util/CoordinateEnums.h"
Expand Down Expand Up @@ -51,27 +52,65 @@ void MatchingMeshPartitionerLonLatPolygon::partition(const Grid& grid, int parti

double west = poly.coordinatesMin().x();
double east = poly.coordinatesMax().x();
comm.allReduceInPlace(west, eckit::mpi::Operation::MIN);
comm.allReduceInPlace(east, eckit::mpi::Operation::MAX);
ATLAS_TRACE_MPI(ALLREDUCE) {
comm.allReduceInPlace(west, eckit::mpi::Operation::MIN);
comm.allReduceInPlace(east, eckit::mpi::Operation::MAX);
}

Projection projection = prePartitionedMesh_.projection();
omp::fill(partitioning, partitioning + grid.size(), -1);

auto compute = [&](double west) {
size_t i = 0;

for (PointLonLat P : grid.lonlat()) {
if (partitioning[i] < 0) {
projection.lonlat2xy(P);
P.normalise(west);
partitioning[i] = poly.contains(P) ? mpi_rank : -1;
#if !defined(__NVCOMPILER)
atlas_omp_parallel
#elif (__NVCOMPILER_MAJOR__ >= 21) && (__NVCOMPILER_MINOR__ > 9 )
// Internal compiler error with nvhpc 21.9:
//
// NVC++-S-0000-Internal compiler error. BAD sptr in var_refsym 0 (MatchingMeshPartitionerLonLatPolygon.cc: 64) (=following line)
// NVC++/x86-64 Linux 21.9-0: compilation completed with severe errors
atlas_omp_parallel
#endif
{
const idx_t num_threads = atlas_omp_get_num_threads();
const idx_t thread_num = atlas_omp_get_thread_num();
const idx_t begin = static_cast<idx_t>(thread_num * size_t(grid.size()) / num_threads);
const idx_t end =
static_cast<idx_t>((thread_num + 1) * size_t(grid.size()) / num_threads);
size_t i = begin;
auto it = grid.lonlat().begin() + i;
for(; i < end; ++i, ++it) {
PointLonLat P = *it;
if (partitioning[i] < 0) {
projection.lonlat2xy(P);
P.normalise(west);
partitioning[i] = poly.contains(P) ? mpi_rank : -1;
}
}
++i;
}
// Synchronize partitioning
comm.allReduceInPlace(partitioning, grid.size(), eckit::mpi::Operation::MAX);
ATLAS_TRACE_MPI(ALLREDUCE) {
comm.allReduceInPlace(partitioning, grid.size(), eckit::mpi::Operation::MAX);
}

return *std::min_element(partitioning, partitioning + grid.size());
std::vector<int> thread_min(atlas_omp_get_max_threads(),std::numeric_limits<int>::max());
#if !defined(__NVCOMPILER)
atlas_omp_parallel
#elif (__NVCOMPILER_MAJOR__ >= 21) && (__NVCOMPILER_MINOR__ > 9 )
// Internal compiler error with nvhpc 21.9:
//
// NVC++-S-0000-Internal compiler error. BAD sptr in var_refsym 0 (MatchingMeshPartitionerLonLatPolygon.cc: 64) (=following line)
// NVC++/x86-64 Linux 21.9-0: compilation completed with severe errors
atlas_omp_parallel
#endif
{
const idx_t num_threads = atlas_omp_get_num_threads();
const idx_t thread_num = atlas_omp_get_thread_num();
const idx_t begin = static_cast<idx_t>(thread_num * size_t(grid.size()) / num_threads);
const idx_t end =
static_cast<idx_t>((thread_num + 1) * size_t(grid.size()) / num_threads);
thread_min[thread_num] = *std::min_element(partitioning+begin,partitioning+end);
}
return *std::min_element(thread_min.begin(), thread_min.end());
};

int min = compute(east - 360.);
Expand Down
40 changes: 19 additions & 21 deletions src/atlas/mesh/actions/BuildHalo.cc
Original file line number Diff line number Diff line change
Expand Up @@ -435,8 +435,6 @@ class Notification {
using Uid2Node = std::unordered_map<uid_t, idx_t>;




void build_lookup_uid2node(Mesh& mesh, Uid2Node& uid2node) {
ATLAS_TRACE();
Notification notes;
Expand All @@ -452,21 +450,21 @@ void build_lookup_uid2node(Mesh& mesh, Uid2Node& uid2node) {
uid_t uid = compute_uid(jnode);
bool inserted = uid2node.insert(std::make_pair(uid, jnode)).second;
if (not inserted) {
int other = uid2node[uid];
std::stringstream msg;
msg << std::setprecision(10) << std::fixed << "Node uid: " << uid << " " << glb_idx(jnode) << " xy("
<< xy(jnode, XX) << "," << xy(jnode, YY) << ")";
if (nodes.has_field("ij")) {
auto ij = array::make_view<idx_t, 2>(nodes.field("ij"));
msg << " ij(" << ij(jnode, XX) << "," << ij(jnode, YY) << ")";
}
msg << " has already been added as node " << glb_idx(other) << " (" << xy(other, XX) << "," << xy(other, YY)
<< ")";
if (nodes.has_field("ij")) {
auto ij = array::make_view<idx_t, 2>(nodes.field("ij"));
msg << " ij(" << ij(other, XX) << "," << ij(other, YY) << ")";
}
notes.add_error(msg.str());
int other = uid2node[uid];
std::stringstream msg;
msg << std::setprecision(10) << std::fixed << "Node uid: " << uid << " " << glb_idx(jnode) << " xy("
<< xy(jnode, XX) << "," << xy(jnode, YY) << ")";
if (nodes.has_field("ij")) {
auto ij = array::make_view<idx_t, 2>(nodes.field("ij"));
msg << " ij(" << ij(jnode, XX) << "," << ij(jnode, YY) << ")";
}
msg << " has already been added as node " << glb_idx(other) << " (" << xy(other, XX) << "," << xy(other, YY)
<< ")";
if (nodes.has_field("ij")) {
auto ij = array::make_view<idx_t, 2>(nodes.field("ij"));
msg << " ij(" << ij(other, XX) << "," << ij(other, YY) << ")";
}
notes.add_error(msg.str());
}
}
if (notes.error()) {
Expand Down Expand Up @@ -495,7 +493,7 @@ void accumulate_elements(const Mesh& mesh, const mpi::BufferView<uid_t>& request

idx_t inode = -1;
// search and get node index for uid
Uid2Node::const_iterator found = uid2node.find(uid);
auto found = uid2node.find(uid);
if (found != uid2node.end()) {
inode = found->second;
}
Expand Down Expand Up @@ -707,7 +705,7 @@ class BuildHaloHelper {
for (it = nodes_uid.begin(); it != nodes_uid.end(); ++it, ++jnode) {
uid_t uid = *it;

Uid2Node::iterator found = uid2node.find(uid);
auto found = uid2node.find(uid);
if (found != uid2node.end()) // Point exists inside domain
{
idx_t node = found->second;
Expand Down Expand Up @@ -773,7 +771,7 @@ class BuildHaloHelper {
for (it = nodes_uid.begin(); it != nodes_uid.end(); ++it, ++jnode) {
uid_t uid = *it;

Uid2Node::iterator found = uid2node.find(uid);
auto found = uid2node.find(uid);
if (found != uid2node.end()) // Point exists inside domain
{
int node = found->second;
Expand Down Expand Up @@ -919,7 +917,7 @@ class BuildHaloHelper {
// make sure new node was not already there
{
uid_t uid = compute_uid(loc_idx);
Uid2Node::iterator found = uid2node.find(uid);
auto found = uid2node.find(uid);
if (found != uid2node.end()) {
int other = found->second;
std::stringstream msg;
Expand Down

0 comments on commit d30709c

Please sign in to comment.