leaf_field.hpp

Go to the documentation of this file.

// core — leaf-indexed fields and world geometry for a BlockOctree.
//
// A BlockOctree addresses leaves in integer *fine units* relative to the block
// origin. AmrGeometry maps those to world coordinates (origin + h0 * fine), and
// LeafField<T> is a value-per-leaf array in the octree's leaf (Z-order) slot
// order — the simulation field the device kernels and the VTU writer operate on.
// The on-device counterpart is simply a peclet::core::View<T> of length numLeaves(); this
// host container is the I/O / setup form.
//
// Header-only, guarded by PECLET_CORE_HAVE_MORTON (so it tracks block_octree.hpp).
#ifndef PECLET_CORE_AMR_LEAF_FIELD_HPP
#define PECLET_CORE_AMR_LEAF_FIELD_HPP
 
#ifdef PECLET_CORE_HAVE_MORTON
 
#include <array>
#include <cstddef>
#include <vector>
 
#include "peclet/core/amr/block_octree.hpp"
#include "peclet/core/common/types.hpp"
 
namespace peclet::core::amr {
 
template <int Dim>
struct AmrGeometry {
  Vec<Dim> origin{};
  Real h0 = 1.0;
 
  Real leafSize(unsigned level) const { return h0 * static_cast<Real>(Index(1) << level); }
 
  template <class Coord>
  Vec<Dim> lowerCorner(const std::array<Coord, Dim>& lo) const {
    Vec<Dim> p{};
    for (int d = 0; d < Dim; ++d)
      p[d] = origin[d] + static_cast<Real>(lo[d]) * h0;
    return p;
  }
 
  template <class Coord>
  Vec<Dim> center(const std::array<std::array<Coord, Dim>, 2>& b) const {
    Vec<Dim> p{};
    for (int d = 0; d < Dim; ++d)
      p[d] = origin[d] +
             (static_cast<Real>(b[0][d]) + 0.5 * (static_cast<Real>(b[1][d] - b[0][d]) + 1.0)) * h0;
    return p;
  }
};
 
template <class T>
struct LeafField {
  std::vector<T> values;
 
  LeafField() = default;
  explicit LeafField(Index n, T init = T{}) : values(static_cast<std::size_t>(n), init) {}
 
  template <int Dim, unsigned Bits>
  explicit LeafField(const BlockOctree<Dim, Bits>& t, T init = T{})
      : values(static_cast<std::size_t>(t.numLeaves()), init) {}
 
  Index size() const { return static_cast<Index>(values.size()); }
  T& operator[](Index i) { return values[static_cast<std::size_t>(i)]; }
  const T& operator[](Index i) const { return values[static_cast<std::size_t>(i)]; }
  T* data() { return values.data(); }
  const T* data() const { return values.data(); }
};
 
}  // namespace peclet::core::amr
 
#endif  // PECLET_CORE_HAVE_MORTON
#endif  // PECLET_CORE_AMR_LEAF_FIELD_HPP