test_petsc_numbering_empty_rank.cc

Go to the documentation of this file.

// SPDX-FileCopyrightText: Copyright (c) NumeriCor GmbH
// SPDX-License-Identifier: Apache-2.0
 
#include "SF_base.h"
#include "basics.h"
 
#include <mpi.h>
 
#include <iostream>
#include <string>
#include <vector>
 
namespace {
 
using SF::vector;
using mesh_t = SF::meshdata<SF_int, SF_real>;
 
constexpr int NX = 4;
constexpr int NY = 4;
constexpr int NUM_PTS = NX * NY;
 
int failures = 0;
 
void check(bool ok, const std::string& what, int rank)
{
  if (!ok) {
    std::cerr << "FAILED [rank " << rank << "]: " << what << '\n';
    failures++;
  }
}
 
void global_indices(int rank, vector<SF_int>& gidx)
{
  gidx.resize(rank == 0 ? NUM_PTS : 0);
  for (size_t i = 0; i < gidx.size(); i++) gidx[i] = SF_int(i);
}
 
void build_mesh(int rank, MPI_Comm comm, mesh_t& mesh)
{
  mesh.comm = comm;
  mesh.g_numpts = NUM_PTS;
 
  vector<SF_int> cnt;
  if (rank == 0) {
    for (int y = 0; y + 1 < NY; y++) {
      for (int x = 0; x + 1 < NX; x++) {
        const SF_int n = SF_int(y * NX + x);
        mesh.con.push_back(n);
        mesh.con.push_back(n + 1);
        mesh.con.push_back(n + NX + 1);
        mesh.con.push_back(n + NX);
        cnt.push_back(4);
      }
    }
  }
 
  mesh.l_numelem = cnt.size();
  SF::dsp_from_cnt(cnt, mesh.dsp);
 
  vector<SF_int> gidx;
  global_indices(rank, gidx);
  mesh.pl.assign(gidx, comm);
  mesh.l_numpts = mesh.pl.num_local_idx();
}
 
void run_case(bool use_rcm, MPI_Comm comm, int rank)
{
  const std::string what = use_rcm ? "renumbered" : "ascending";
 
  mesh_t mesh;
  build_mesh(rank, comm, mesh);
 
  const vector<SF_int>& alg_nod = mesh.pl.algebraic_nodes();
  if (rank == 0)
    check(alg_nod.size() == NUM_PTS, what + ": rank 0 owns the whole grid", rank);
  else
    check(alg_nod.size() == 0, what + ": rank owns no algebraic node", rank);
 
  SF::petsc_numbering<SF_int, SF_real> numbering(mesh.pl, use_rcm);
  numbering(mesh);
 
  const vector<SF_int>& petsc_idx = mesh.get_numbering(SF::NBR_PETSC);
  const vector<SF_int>& alg_dsp = mesh.pl.algebraic_layout();
  const SF_int alg_start = alg_dsp[rank];
 
  // The invariant that makes a renumbering invisible to everything reading the local
  // PETSc arrays by algebraic position.
  bool positional = true;
  for (size_t i = 0; i < alg_nod.size(); i++)
    if (petsc_idx[alg_nod[i]] != alg_start + SF_int(i)) positional = false;
  check(positional, what + ": algebraic position matches the PETSc slot", rank);
 
  // Renumbering permutes the algebraic node set, it must not change which nodes are in it.
  std::vector<bool> seen(mesh.pl.num_local_idx(), false);
  bool set_intact = true;
  for (size_t i = 0; i < alg_nod.size(); i++) {
    const SF_int n = alg_nod[i];
    if (n < 0 || size_t(n) >= seen.size() || seen[n]) set_intact = false;
    else seen[n] = true;
  }
  check(set_intact, what + ": algebraic node set is free of gaps and duplicates", rank);
}
 
}  // namespace
 
int main(int argc, char** argv)
{
  MPI_Init(&argc, &argv);
 
  MPI_Comm comm = MPI_COMM_WORLD;
  int size, rank;
  MPI_Comm_size(comm, &size);
  MPI_Comm_rank(comm, &rank);
 
  if (size < 2) {
    if (!rank) std::cerr << "FAILED: this test needs at least 2 ranks to leave one empty\n";
    MPI_Finalize();
    return 1;
  }
 
  run_case(false, comm, rank);
  run_case(true, comm, rank);
 
  int total = 0;
  MPI_Allreduce(&failures, &total, 1, MPI_INT, MPI_SUM, comm);
 
  MPI_Finalize();
  return total ? 1 : 0;
}