MueLu Version of the Day
Loading...
Searching...
No Matches
MueLu_IndexManager_kokkos_decl.hpp
Go to the documentation of this file.
1// @HEADER
2//
3// ***********************************************************************
4//
5// MueLu: A package for multigrid based preconditioning
6// Copyright 2012 Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact
39// Jonathan Hu (jhu@sandia.gov)
40// Ray Tuminaro (rstumin@sandia.gov)
41// Luc Berger-Vergiat (lberge@sandia.gov)
42//
43// ***********************************************************************
44//
45// @HEADER
46#ifndef MUELU_INDEXMANAGER_KOKKOS_DECL_HPP
47#define MUELU_INDEXMANAGER_KOKKOS_DECL_HPP
48
49#include "MueLu_ConfigDefs.hpp"
50#include "MueLu_Types.hpp"
51
52#include <KokkosCompat_ClassicNodeAPI_Wrapper.hpp>
53
54#include "Teuchos_OrdinalTraits.hpp"
55
56#include <Xpetra_Map_fwd.hpp>
57#include <Xpetra_Vector_fwd.hpp>
58#include <Xpetra_VectorFactory_fwd.hpp>
59
60#include "MueLu_BaseClass.hpp"
62
63/*****************************************************************************
64
65****************************************************************************/
66
67namespace MueLu {
68
80 template <class LocalOrdinal, class GlobalOrdinal, class Node>
82#undef MUELU_INDEXMANAGER_KOKKOS_SHORT
84
85 public:
86 using execution_space = typename Node::execution_space;
87 using memory_space = typename Node::memory_space;
88 using device_type = Kokkos::Device<execution_space, memory_space>;
89 using intTupleView = typename Kokkos::View<int[3], device_type>;
90 using LOTupleView = typename Kokkos::View<LO[3], device_type>;
91
92 private:
93
94 const int meshLayout = UNCOUPLED;
95 int myRank = -1;
100
104
108
109 public:
110
113
115 IndexManager_kokkos(const int NumDimensions,
116 const int interpolationOrder,
117 const int MyRank,
118 const ArrayView<const LO> LFineNodesPerDir,
119 const ArrayView<const int> CoarseRate);
120
122
124 void setupIM(const int NumDimensions,
125 const int interpolationOrder,
126 const ArrayView<const int> coarseRate,
127 const ArrayView<const LO> LFineNodesPerDir);
128
132
133 int getNumDimensions() const {return numDimensions;}
134
136
138
140
141 KOKKOS_INLINE_FUNCTION
143
144 KOKKOS_INLINE_FUNCTION
146
147 KOKKOS_INLINE_FUNCTION
149
150 KOKKOS_INLINE_FUNCTION
152
153 Array<LO> getCoarseNodesPerDirArray() const;
154
155 KOKKOS_INLINE_FUNCTION
156 void getFineLID2FineTuple(const LO myLID, LO (&tuple)[3]) const {
157 LO tmp;
158 tuple[2] = myLID / (lFineNodesPerDir(1)*lFineNodesPerDir(0));
159 tmp = myLID % (lFineNodesPerDir(1)*lFineNodesPerDir(0));
160 tuple[1] = tmp / lFineNodesPerDir(0);
161 tuple[0] = tmp % lFineNodesPerDir(0);
162 } // getFineNodeLocalTuple
163
164 KOKKOS_INLINE_FUNCTION
165 void getFineTuple2FineLID(const LO tuple[3], LO& myLID) const {
166 myLID = tuple[2]*lNumFineNodes10 + tuple[1]*lFineNodesPerDir[0] + tuple[0];
167 } // getFineNodeLID
168
169 KOKKOS_INLINE_FUNCTION
170 void getCoarseLID2CoarseTuple(const LO myLID, LO (&tuple)[3]) const {
171 LO tmp;
172 tuple[2] = myLID / numCoarseNodes10;
173 tmp = myLID % numCoarseNodes10;
174 tuple[1] = tmp / coarseNodesPerDir[0];
175 tuple[0] = tmp % coarseNodesPerDir[0];
176 } // getCoarseNodeLocalTuple
177
178 KOKKOS_INLINE_FUNCTION
179 void getCoarseTuple2CoarseLID(const LO i, const LO j, const LO k, LO& myLID) const {
180 myLID = k*numCoarseNodes10 + j*coarseNodesPerDir[0] + i;
181 } // getCoarseNodeLID
182
183 };
184
185} //namespace MueLu
186
187#define MUELU_INDEXMANAGER_KOKKOS_SHORT
188#endif // MUELU_INDEXMANAGER_KOKKOS_DECL_HPP
Base class for MueLu classes.
Container class for mesh layout and indices calculation.
LO numCoarseNodes10
local number of nodes per 0-1 slice remaining after coarsening.
LOTupleView coarseNodesPerDir
local number of nodes per direction remaing after coarsening.
KOKKOS_INLINE_FUNCTION LOTupleView getCoarseNodesPerDir() const
intTupleView endRate
adapted coarsening rate at the edge of the mesh in each direction.
LO lNumFineNodes10
local number of nodes per 0-1 slice.
Kokkos::Device< execution_space, memory_space > device_type
KOKKOS_INLINE_FUNCTION void getCoarseLID2CoarseTuple(const LO myLID, LO(&tuple)[3]) const
KOKKOS_INLINE_FUNCTION void getCoarseTuple2CoarseLID(const LO i, const LO j, const LO k, LO &myLID) const
KOKKOS_INLINE_FUNCTION LOTupleView getLocalFineNodesPerDir() const
typename Node::execution_space execution_space
KOKKOS_INLINE_FUNCTION void getFineLID2FineTuple(const LO myLID, LO(&tuple)[3]) const
LOTupleView lFineNodesPerDir
local number of nodes per direction.
typename Kokkos::View< LO[3], device_type > LOTupleView
LO numCoarseNodes
local number of nodes remaining after coarsening.
KOKKOS_INLINE_FUNCTION intTupleView getCoarseningRates() const
IndexManager_kokkos()=default
Default constructor, return empty object.
intTupleView coarseRate
coarsening rate in each direction
KOKKOS_INLINE_FUNCTION intTupleView getCoarseningEndRates() const
int numDimensions
Number of spacial dimensions in the problem.
KOKKOS_INLINE_FUNCTION void getFineTuple2FineLID(const LO tuple[3], LO &myLID) const
int interpolationOrder_
Interpolation order used by grid transfer operators using these aggregates.
void setupIM(const int NumDimensions, const int interpolationOrder, const ArrayView< const int > coarseRate, const ArrayView< const LO > LFineNodesPerDir)
Common setup pattern used for all the different types of undelying mesh.
typename Kokkos::View< int[3], device_type > intTupleView
Namespace for MueLu classes and methods.