Store block gid and neighbor refinement levels in sparse packs

CHANGELOG.md

@@ -25,6 +25,7 @@
 Date: 2024-08-30
 
 ### Added (new features/APIs/variables/...)
+- [[PR 1167]](https://github.com/parthenon-hpc-lab/parthenon/pull/1167) Store block gid and neighbor refinement levels in sparse packs
 - [[PR 1151]](https://github.com/parthenon-hpc-lab/parthenon/pull/1151) Add time offset `c` to LowStorageIntegrator
 - [[PR 1147]](https://github.com/parthenon-hpc-lab/parthenon/pull/1147) Add `par_reduce_inner` functions
 - [[PR 1159]](https://github.com/parthenon-hpc-lab/parthenon/pull/1159) Add additional timestep controllers in parthenon/time.

src/interface/sparse_pack.hpp

@@ -189,6 +189,19 @@ class SparsePack : public SparsePackBase {
     return bounds_h_(1, b, vidx);
   }
 
+  KOKKOS_INLINE_FUNCTION int GetLevel(const int b, const int off3, const int off2,
+                                      const int off1) const {
+    return block_props_(b, (off3 + 1) + 3 * ((off2 + 1) + 3 * (off1 + 1)));
+  }
+
+  KOKKOS_INLINE_FUNCTION int GetGID(const int b) const { return block_props_(b, 27); }
+
+  int GetLevelHost(const int b, const int off3, const int off2, const int off1) const {
+    return block_props_h_(b, (off3 + 1) + 3 * ((off2 + 1) + 3 * (off1 + 1)));
+  }
+
+  int GetGIDHost(const int b) const { return block_props_h_(b, 27); }
+
   // Number of components of a variable on a block
   template <typename T>
   KOKKOS_INLINE_FUNCTION int GetSize(const int b, const T &t) const {

src/interface/sparse_pack_base.cpp

@@ -162,6 +162,10 @@ SparsePackBase SparsePackBase::Build(T *pmd, const PackDescriptor &desc,
   pack.bounds_ = bounds_t("bounds", 2, nblocks, nvar + 1);
   pack.bounds_h_ = Kokkos::create_mirror_view(pack.bounds_);
 
+  // This array stores refinement levels of current block and all neighboring blocks.
+  pack.block_props_ = block_props_t("block_props", nblocks, 27 + 1);
+  pack.block_props_h_ = Kokkos::create_mirror_view(pack.block_props_);
+
   pack.coords_ = coords_t("coords", desc.flat ? max_size : nblocks);
   auto coords_h = Kokkos::create_mirror_view(pack.coords_);
 
@@ -180,6 +184,24 @@ SparsePackBase SparsePackBase::Build(T *pmd, const PackDescriptor &desc,
       coords_h(b) = pmbd->GetBlockPointer()->coords_device;
     }
 
+    // Initialize block refinement levels to current block level to provide default if
+    // neighbors not present
+    for (int n = 0; n < 27; n++) {
+      pack.block_props_h_(blidx, (1 + 3 * (1 + 3 * 1))) =
+          pmbd->GetBlockPointer()->loc.level();
+    }
+    // This block's gid stored in central (1, 1, 1, 1) element
+    pack.block_props_h_(blidx, 27) = pmbd->GetBlockPointer()->gid;
+    for (auto &neighbor : pmbd->GetBlockPointer()->neighbors) {
+      // Multiple refined neighbors may write to the same index but they will always have
+      // the same refinement level.
+      pack.block_props_h_(
+          blidx, (neighbor.offsets[2] + 1) +
+                     3 * ((neighbor.offsets[1] + 1) + 3 * (neighbor.offsets[0] + 1))) =
+          neighbor.loc.level();
+      // Currently not storing neighbor gids
+    }
+
     for (int i = 0; i < nvar; ++i) {
       pack.bounds_h_(0, blidx, i) = idx;
       for (const auto &[var_name, uid] : desc.var_groups[i]) {
@@ -278,6 +300,7 @@ SparsePackBase SparsePackBase::Build(T *pmd, const PackDescriptor &desc,
   });
   Kokkos::deep_copy(pack.pack_, pack.pack_h_);
   Kokkos::deep_copy(pack.bounds_, pack.bounds_h_);
+  Kokkos::deep_copy(pack.block_props_, pack.block_props_h_);
   Kokkos::deep_copy(pack.coords_, coords_h);
 
   return pack;

src/interface/sparse_pack_base.hpp

@@ -59,6 +59,8 @@ class SparsePackBase {
   using pack_h_t = typename pack_t::HostMirror;
   using bounds_t = ParArray3D<int>;
   using bounds_h_t = typename bounds_t::HostMirror;
+  using block_props_t = ParArray2D<int>;
+  using block_props_h_t = typename block_props_t::HostMirror;
   using coords_t = ParArray1D<ParArray0D<Coordinates_t>>;
 
   // Returns a SparsePackBase object that is either newly created or taken
@@ -90,6 +92,8 @@ class SparsePackBase {
   pack_h_t pack_h_;
   bounds_t bounds_;
   bounds_h_t bounds_h_;
+  block_props_t block_props_;
+  block_props_h_t block_props_h_;
   coords_t coords_;
 
   int flx_idx_;

src/outputs/parthenon_xdmf.cpp

@@ -211,7 +211,6 @@ void genXDMF(std::string hdfFile, Mesh *pm, SimTime *tm, IndexDomain domain, int
   if (swarm_xdmf && all_swarm_info.all_info.size() > 0) {
     std::string sfilename_aux = hdfFile + ".swarm.xdmf";
     std::ofstream pxdmf;
-    hsize_t dims[H5_NDIM] = {0}; // zero-initialized
 
     // open file
     pxdmf = std::ofstream(sfilename_aux.c_str(), std::ofstream::trunc);