Skip to content

Commit

Permalink
Add flat sparse packs (#877)
Browse files Browse the repository at this point in the history 
* add sparse flat packs

* formatting

* docs

* changelog

* brryan comments

* Update doc/sphinx/src/interface/sparse.rst

Co-authored-by: Philipp Grete <[email protected]>

* incorporate lroberts suggestions

* fix edge case for no vars found in flat case

* oops cellvars are vars now

* testing lower/upper bound on host, not device

* re-enable getlowerobund getupperbound

* fix bug related to lower bounds

* Move allocated check

* format

---------

Co-authored-by: Jonah Maxwell Miller <[email protected]>
Co-authored-by: Philipp Grete <[email protected]>
Co-authored-by: Luke Roberts <[email protected]>
  • Loading branch information
@Yurlungur @jonahm-LANL
@pgrete @lroberts36
4 people authored Jun 1, 2023
1 parent 099c883 commit 74c1d05
Show file tree
Hide file tree
Showing 6 changed files with 217 additions and 62 deletions.
5 changes: 4 additions & 1 deletion CHANGELOG.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,10 @@

## Current develop

### Added (new features/APIs/variables/...)
- [[PR 877]](https://github.com/parthenon-hpc-lab/parthenon/pull/877) Add flat sparse packs
- [[PR 868]](https://github.com/parthenon-hpc-lab/parthenon/pull/868) Add block-local face, edge, and nodal fields and allow for packing

### Changed (changing behavior/API/variables/...)

### Fixed (not changing behavior/API/variables/...)
Expand All @@ -15,7 +19,6 @@
Date: 2023-05-26

### Added (new features/APIs/variables/...)
- [[PR 868]](https://github.com/parthenon-hpc-lab/parthenon/pull/868) Add block-local face, edge, and nodal fields and allow for packing
- [[PR 830]](https://github.com/parthenon-hpc-lab/parthenon/pull/830) Add particle output
- [[PR 840]](https://github.com/parthenon-hpc-lab/parthenon/pull/840) Generalized integrators infrastructure in a backwards compatible way
- [[PR 810]](https://github.com/parthenon-hpc-lab/parthenon/pull/810) Add suport for Ascent in-situ visualization
Expand Down
10 changes: 10 additions & 0 deletions doc/sphinx/src/interface/sparse.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -137,6 +137,16 @@ simulation. The set of possible operations for a sparse field are:
``lroberts36/merge-sparse-with-jdolence-sparse`` that is being used
in ``Riot``. This should probably be brought into ``develop``, since
the “sparse sparse pack” access pattern is probably not desirable.
- *Flatten outer indices:* Most common packs over ``MeshData`` produce
a 5-dimensional data structure, where the slowest two moving indices
are over ``MeshBlock``\ s then tensor/vector components of
``Variable``\ s. However, it is sometimes desirable to construct a
4-dimensional data structure where the slowest moving index ranges
over all variables and all blocks. So for example if we had a system
with 9 blocks and 5 variables per block, the slowest moving index
would be of size 45. This can be enabled ``SparsePack``s with the
``GetFlat`` series of factory functions or by passing the optional
``flat`` boolean into the constructor.

In comparison to a sparse field, a dense field only requires the
operation *Access*.
Expand Down
124 changes: 111 additions & 13 deletions src/interface/sparse_pack.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -125,10 +125,10 @@ class SparsePack : public SparsePackBase {
// The pack will be created and accessible on the device
template <class T>
static SparsePack Get(T *pmd, const std::vector<MetadataFlag> &flags = {},
bool fluxes = false, bool coarse = false) {
bool fluxes = false, bool coarse = false, bool flatten = false) {
const impl::PackDescriptor desc(std::vector<std::string>{Ts::name()...},
std::vector<bool>{Ts::regex()...}, flags, fluxes,
coarse);
coarse, flatten);
return SparsePack(SparsePackBase::GetPack(pmd, desc));
}

Expand All @@ -143,9 +143,9 @@ class SparsePack : public SparsePackBase {
template <class T, class VAR_VEC>
static std::tuple<SparsePack, SparsePackIdxMap>
Get(T *pmd, const VAR_VEC &vars, const std::vector<MetadataFlag> &flags = {},
bool fluxes = false, bool coarse = false) {
bool fluxes = false, bool coarse = false, bool flatten = false) {
static_assert(sizeof...(Ts) == 0, "Cannot create a string/type hybrid pack");
impl::PackDescriptor desc(vars, flags, fluxes, coarse);
impl::PackDescriptor desc(vars, flags, fluxes, coarse, flatten);
return {SparsePack(SparsePackBase::GetPack(pmd, desc)),
SparsePackBase::GetIdxMap(desc)};
}
Expand All @@ -155,7 +155,7 @@ class SparsePack : public SparsePackBase {
static SparsePack GetWithFluxes(T *pmd, const std::vector<MetadataFlag> &flags = {}) {
const bool coarse = false;
const bool fluxes = true;
return Get(pmd, flags, fluxes, coarse);
return Get(pmd, AddFlag_(flags), fluxes, coarse);
}

template <class T, class VAR_VEC>
Expand All @@ -164,7 +164,7 @@ class SparsePack : public SparsePackBase {
const std::vector<MetadataFlag> &flags = {}) {
const bool coarse = false;
const bool fluxes = true;
Get(pmd, vars, flags, fluxes, coarse);
Get(pmd, vars, AddFlag_(flags), fluxes, coarse);
}

template <class T>
Expand All @@ -183,21 +183,79 @@ class SparsePack : public SparsePackBase {
Get(pmd, vars, flags, fluxes, coarse);
}

template <class T>
static SparsePack GetFlat(T *pmd, const std::vector<MetadataFlag> &flags = {}) {
const bool coarse = false;
const bool fluxes = false;
const bool flatten = true;
return Get(pmd, flags, fluxes, coarse, flatten);
}

template <class T, class VAR_VEC>
static std::tuple<SparsePack, SparsePackIdxMap>
GetFlat(T *pmd, const VAR_VEC &vars, const std::vector<MetadataFlag> &flags = {}) {
const bool coarse = false;
const bool fluxes = false;
const bool flatten = true;
return Get(pmd, vars, flags, fluxes, coarse, flatten);
}

template <class T>
static SparsePack GetFlatWithFluxes(T *pmd,
const std::vector<MetadataFlag> &flags = {}) {
const bool coarse = false;
const bool fluxes = true;
const bool flatten = true;
return Get(pmd, AddFlag_(flags), fluxes, coarse, flatten);
}

template <class T, class VAR_VEC>
static std::tuple<SparsePack, SparsePackIdxMap>
GetFlatWithFluxes(T *pmd, const VAR_VEC &vars,
const std::vector<MetadataFlag> &flags = {}) {
const bool coarse = false;
const bool fluxes = true;
const bool flatten = true;
return Get(pmd, vars, AddFlag_(flags), fluxes, coarse, flatten);
}

template <class T>
static SparsePack GetFlatWithCoarse(T *pmd,
const std::vector<MetadataFlag> &flags = {}) {
const bool coarse = true;
const bool fluxes = false;
const bool flatten = true;
return Get(pmd, flags, fluxes, coarse, flatten);
}

template <class T, class VAR_VEC>
static std::tuple<SparsePack, SparsePackIdxMap>
GetFlatWithCoarse(T *pmd, const VAR_VEC &vars,
const std::vector<MetadataFlag> &flags = {}) {
const bool coarse = true;
const bool fluxes = false;
const bool flatten = true;
return Get(pmd, vars, flags, fluxes, coarse, flatten);
}

// Methods for getting parts of the shape of the pack
KOKKOS_FORCEINLINE_FUNCTION
int GetNBlocks() const { return nblocks_; }

KOKKOS_FORCEINLINE_FUNCTION
int GetNDim() const { return ndim_; }
int GetMaxNumberOfVars() const { return pack_.extent_int(2); }

// Returns the total number of vars/components in pack
KOKKOS_FORCEINLINE_FUNCTION
int GetMaxNumberOfVars() const { return pack_.extent_int(2); }
int GetSize() const { return size_; }

KOKKOS_INLINE_FUNCTION
const Coordinates_t &GetCoordinates(const int b) const { return coords_(b)(); }
const Coordinates_t &GetCoordinates(const int b = 0) const { return coords_(b)(); }

// Bound overloads
KOKKOS_INLINE_FUNCTION int GetLowerBound(const int b) const { return 0; }
KOKKOS_INLINE_FUNCTION int GetLowerBound(const int b) const {
return (flat_ && (b > 0)) ? (bounds_(1, b - 1, nvar_) + 1) : 0;
}

KOKKOS_INLINE_FUNCTION int GetUpperBound(const int b) const {
return bounds_(1, b, nvar_);
Expand Down Expand Up @@ -226,7 +284,9 @@ class SparsePack : public SparsePackBase {
}

// Host Bound overloads
KOKKOS_INLINE_FUNCTION int GetLowerBoundHost(const int b) const { return 0; }
KOKKOS_INLINE_FUNCTION int GetLowerBoundHost(const int b) const {
return (flat_ && (b > 0)) ? (bounds_h_(1, b - 1, nvar_) + 1) : 0;
}

KOKKOS_INLINE_FUNCTION int GetUpperBoundHost(const int b) const {
return bounds_h_(1, b, nvar_);
Expand Down Expand Up @@ -268,6 +328,7 @@ class SparsePack : public SparsePackBase {

KOKKOS_INLINE_FUNCTION auto &operator()(const int b, const TE el, PackIdx idx) const {
static_assert(sizeof...(Ts) == 0);
PARTHENON_DEBUG_REQUIRE(!flat_, "Accessor cannot be used for flat packs");
const int n = bounds_(0, b, idx.VariableIdx()) + idx.Offset();
return pack_(static_cast<int>(el) % 3, b, n);
}
Expand All @@ -283,14 +344,22 @@ class SparsePack : public SparsePackBase {

template <class TIn, REQUIRES(IncludesType<TIn, Ts...>::value)>
KOKKOS_INLINE_FUNCTION auto &operator()(const int b, const TIn &t) const {
PARTHENON_DEBUG_REQUIRE(!flat_, "Accessor cannot be used for flat packs");
const int vidx = GetLowerBound(b, t) + t.idx;
return pack_(0, b, vidx);
}

KOKKOS_INLINE_FUNCTION Real &operator()(const int b, const int idx, const int k,
const int j, const int i) const {
KOKKOS_INLINE_FUNCTION
Real &operator()(const int b, const int idx, const int k, const int j,
const int i) const {
PARTHENON_DEBUG_REQUIRE(!flat_, "Accessor cannot be used for flat packs");
return pack_(0, b, idx)(k, j, i);
}
KOKKOS_INLINE_FUNCTION
Real &operator()(int idx, const int k, const int j, const int i) const {
PARTHENON_DEBUG_REQUIRE(flat_, "Accessor valid only for flat packs");
return pack_(0, 0, idx)(k, j, i);
}

KOKKOS_INLINE_FUNCTION Real &operator()(const int b, const TE el, const int idx,
const int k, const int j, const int i) const {
Expand All @@ -300,6 +369,7 @@ class SparsePack : public SparsePackBase {
KOKKOS_INLINE_FUNCTION Real &operator()(const int b, PackIdx idx, const int k,
const int j, const int i) const {
static_assert(sizeof...(Ts) == 0, "Cannot create a string/type hybrid pack");
PARTHENON_DEBUG_REQUIRE(!flat_, "Accessor cannot be used for flat packs");
const int n = bounds_(0, b, idx.VariableIdx()) + idx.Offset();
return pack_(0, b, n)(k, j, i);
}
Expand All @@ -314,6 +384,7 @@ class SparsePack : public SparsePackBase {
template <class TIn, REQUIRES(IncludesType<TIn, Ts...>::value)>
KOKKOS_INLINE_FUNCTION Real &operator()(const int b, const TIn &t, const int k,
const int j, const int i) const {
PARTHENON_DEBUG_REQUIRE(!flat_, "Accessor cannot be used for flat packs");
const int vidx = GetLowerBound(b, t) + t.idx;
return pack_(0, b, vidx)(k, j, i);
}
Expand All @@ -328,21 +399,31 @@ class SparsePack : public SparsePackBase {
// flux() overloads
KOKKOS_INLINE_FUNCTION
auto &flux(const int b, const int dir, const int idx) const {
PARTHENON_DEBUG_REQUIRE(!flat_, "Accessor cannot be used for flat packs");
PARTHENON_DEBUG_REQUIRE(dir > 0 && dir < 4 && with_fluxes_, "Bad input to flux call");
return pack_(dir, b, idx);
}

KOKKOS_INLINE_FUNCTION
Real &flux(const int b, const int dir, const int idx, const int k, const int j,
const int i) const {
PARTHENON_DEBUG_REQUIRE(!flat_, "Accessor cannot be used for flat packs");
PARTHENON_DEBUG_REQUIRE(dir > 0 && dir < 4 && with_fluxes_, "Bad input to flux call");
return pack_(dir, b, idx)(k, j, i);
}

KOKKOS_INLINE_FUNCTION
Real &flux(const int dir, const int idx, const int k, const int j, const int i) const {
PARTHENON_DEBUG_REQUIRE(flat_, "Accessor must only be used for flat packs");
PARTHENON_DEBUG_REQUIRE(dir > 0 && dir < 4 && with_fluxes_, "Bad input to flux call");
return pack_(dir, 0, idx)(k, j, i);
}

KOKKOS_INLINE_FUNCTION
Real &flux(const int b, const int dir, PackIdx idx, const int k, const int j,
const int i) const {
static_assert(sizeof...(Ts) == 0, "Cannot create a string/type hybrid pack");
PARTHENON_DEBUG_REQUIRE(!flat_, "Accessor cannot be used for flat packs");
PARTHENON_DEBUG_REQUIRE(dir > 0 && dir < 4 && with_fluxes_, "Bad input to flux call");
const int n = bounds_(0, b, idx.VariableIdx()) + idx.Offset();
return pack_(dir, b, n)(k, j, i);
Expand All @@ -351,10 +432,27 @@ class SparsePack : public SparsePackBase {
template <class TIn, REQUIRES(IncludesType<TIn, Ts...>::value)>
KOKKOS_INLINE_FUNCTION Real &flux(const int b, const int dir, const TIn &t, const int k,
const int j, const int i) const {
PARTHENON_DEBUG_REQUIRE(!flat_, "Accessor cannot be used for flat packs");
PARTHENON_DEBUG_REQUIRE(dir > 0 && dir < 4 && with_fluxes_, "Bad input to flux call");
const int vidx = GetLowerBound(b, t) + t.idx;
return pack_(dir, b, vidx)(k, j, i);
}

private:
// Must make a copy of the vector, since input vector is const,
// and it may not even be an lvalue.
static std::vector<MetadataFlag> AddFlag_(const std::vector<MetadataFlag> &flags,
MetadataFlag mf = Metadata::WithFluxes) {
if (std::find(flags.begin(), flags.end(), mf) == flags.end()) {
std::vector<MetadataFlag> out;
out.reserve(flags.size() + 1);
out.insert(out.begin(), flags.begin(), flags.end());
out.push_back(mf);
return out;
} else {
return flags;
}
}
};

} // namespace parthenon
Expand Down
Loading

0 comments on commit 74c1d05

Please sign in to comment.