46 #ifndef MUELU_INDEXMANAGER_KOKKOS_DECL_HPP 47 #define MUELU_INDEXMANAGER_KOKKOS_DECL_HPP 52 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 53 #include <KokkosCompat_ClassicNodeAPI_Wrapper.hpp> 54 #include <Kokkos_View.hpp> 56 #include "Teuchos_OrdinalTraits.hpp" 58 #include <Xpetra_Map_fwd.hpp> 59 #include <Xpetra_Vector_fwd.hpp> 60 #include <Xpetra_VectorFactory_fwd.hpp> 82 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
83 class IndexManager_kokkos :
public BaseClass {
84 #undef MUELU_INDEXMANAGER_KOKKOS_SHORT 88 using execution_space =
typename Node::execution_space;
89 using memory_space =
typename Node::memory_space;
90 using device_type = Kokkos::Device<execution_space, memory_space>;
99 int interpolationOrder_;
100 intTupleView coarseRate;
101 intTupleView endRate;
105 LOTupleView lFineNodesPerDir;
109 LOTupleView coarseNodesPerDir;
114 IndexManager_kokkos() =
default;
117 IndexManager_kokkos(
const int NumDimensions,
118 const int interpolationOrder,
120 const ArrayView<const LO> LFineNodesPerDir,
121 const ArrayView<const int> CoarseRate);
123 virtual ~IndexManager_kokkos() {}
126 void setupIM(
const int NumDimensions,
127 const int interpolationOrder,
128 const ArrayView<const int> coarseRate,
129 const ArrayView<const LO> LFineNodesPerDir);
133 void computeMeshParameters();
135 int getNumDimensions()
const {
return numDimensions;}
137 int getInterpolationOrder()
const {
return interpolationOrder_;}
139 LO getNumLocalFineNodes()
const {
return lNumFineNodes;}
141 LO getNumCoarseNodes()
const {
return numCoarseNodes;}
143 KOKKOS_INLINE_FUNCTION
144 intTupleView getCoarseningRates()
const {
return coarseRate;}
146 KOKKOS_INLINE_FUNCTION
147 intTupleView getCoarseningEndRates()
const {
return endRate;}
149 KOKKOS_INLINE_FUNCTION
150 LOTupleView getLocalFineNodesPerDir()
const {
return lFineNodesPerDir;}
152 KOKKOS_INLINE_FUNCTION
153 LOTupleView getCoarseNodesPerDir()
const {
return coarseNodesPerDir;}
155 Array<LO> getCoarseNodesPerDirArray()
const;
157 KOKKOS_INLINE_FUNCTION
158 void getFineLID2FineTuple(
const LO myLID, LO (&tuple)[3])
const {
160 tuple[2] = myLID / (lFineNodesPerDir(1)*lFineNodesPerDir(0));
161 tmp = myLID % (lFineNodesPerDir(1)*lFineNodesPerDir(0));
162 tuple[1] = tmp / lFineNodesPerDir(0);
163 tuple[0] = tmp % lFineNodesPerDir(0);
166 KOKKOS_INLINE_FUNCTION
167 void getFineTuple2FineLID(
const LO tuple[3], LO& myLID)
const {
168 myLID = tuple[2]*lNumFineNodes10 + tuple[1]*lFineNodesPerDir[0] + tuple[0];
171 KOKKOS_INLINE_FUNCTION
172 void getCoarseLID2CoarseTuple(
const LO myLID, LO (&tuple)[3])
const {
174 tuple[2] = myLID / numCoarseNodes10;
175 tmp = myLID % numCoarseNodes10;
176 tuple[1] = tmp / coarseNodesPerDir[0];
177 tuple[0] = tmp % coarseNodesPerDir[0];
180 KOKKOS_INLINE_FUNCTION
181 void getCoarseTuple2CoarseLID(
const LO i,
const LO j,
const LO k, LO& myLID)
const {
182 myLID = k*numCoarseNodes10 + j*coarseNodesPerDir[0] + i;
189 #define MUELU_INDEXMANAGER_KOKKOS_SHORT 190 #endif // HAVE_MUELU_KOKKOS_REFACTOR 191 #endif // MUELU_INDEXMANAGER_KOKKOS_DECL_HPP
Namespace for MueLu classes and methods.