decomposition.hpp

Go to the documentation of this file.

// SPDX-FileCopyrightText: 2025 VTT Technical Research Centre of Finland Ltd
// SPDX-License-Identifier: AGPL-3.0-or-later
 
#pragma once
 
#include "openpfc/core/world.hpp"
#include <array>
#include <heffte.h>
#include <ostream>
#include <stdexcept>
#include <vector>
 
namespace pfc {
 
namespace world {
 
template <typename T>
inline auto create(const World<T> &world, const heffte::box3d<int> &box) {
  return World(box.low, box.high, get_coordinate_system(world));
}
 
template <typename T> inline auto to_indices(const World<T> &world) {
  std::array<int, 3> lower = get_lower(world);
  std::array<int, 3> upper = get_upper(world);
  return heffte::box3d<int>(lower, upper);
}
 
template <typename T>
inline auto split_world(const World<T> &world, const Int3 &grid) {
  std::vector<World<T>> sub_worlds;
  for (const auto &box : heffte::split_world(to_indices(world), grid)) {
    sub_worlds.push_back(create(world, box));
  }
  return sub_worlds;
}
 
} // namespace world
 
namespace decomposition {
 
using pfc::csys::CartesianTag;
using pfc::types::Bool3;
using pfc::types::Int3;
using pfc::types::Real3;
 
using World = pfc::world::World<CartesianTag>;
using Int3 = pfc::types::Int3;
 
using heffte::proc_setup_min_surface;
 
struct Decomposition {
 
  const pfc::World &m_global_world; 
  const std::array<int, 3> m_grid;  
  const std::vector<pfc::World> m_subworlds; 
 
  Decomposition(const World &world, const Int3 grid)
      : m_global_world(world), m_grid(grid), m_subworlds(split_world(world, grid)) {}
 
  friend std::ostream &operator<<(std::ostream &os, const Decomposition &d) {
    os << "Decomposition:\n";
    os << "  Global World: " << d.m_global_world << "\n";
    os << "  Grid: [" << d.m_grid.at(0) << ", " << d.m_grid.at(1) << ", "
       << d.m_grid.at(2) << "]\n";
    os << "  Sub-worlds:\n";
    for (size_t i = 0; i < d.m_subworlds.size(); ++i) {
      os << "    Sub-world " << i << ": " << d.m_subworlds[i] << "\n";
    }
    return os;
  }
};
 
inline const auto &get_global_world(const Decomposition &decomposition) noexcept {
  return decomposition.m_global_world;
}
 
inline const auto &get_world(const Decomposition &decomposition) noexcept {
  return get_global_world(decomposition);
}
 
inline const auto &get_grid(const Decomposition &decomposition) noexcept {
  return decomposition.m_grid;
}
 
inline const auto &get_subworlds(const Decomposition &decomposition) noexcept {
  return decomposition.m_subworlds;
}
 
inline const auto &get_subworld(const Decomposition &decomposition, int i) {
  return get_subworlds(decomposition).at(i);
}
 
inline auto create(const World &world, const Int3 &grid) noexcept {
  return Decomposition(world, grid);
};
 
inline auto create(const World &world, const int &nparts) noexcept {
  auto indices = to_indices(world);
  auto grid = proc_setup_min_surface(indices, nparts);
  return create(world, grid);
}
 
inline int get_num_domains(const Decomposition &decomposition) noexcept {
  return get_subworlds(decomposition).size();
}
 
} // namespace decomposition
 
using Decomposition = decomposition::Decomposition;
 
} // namespace pfc