46 #ifndef MUELU_STRUCTUREDAGGREGATIONFACTORY_KOKKOS_DEF_HPP
47 #define MUELU_STRUCTUREDAGGREGATIONFACTORY_KOKKOS_DEF_HPP
49 #ifdef HAVE_MUELU_KOKKOS_REFACTOR
59 #include "MueLu_Utilities.hpp"
62 #include "MueLu_LWGraph_kokkos.hpp"
63 #include "MueLu_Aggregates_kokkos.hpp"
64 #include "MueLu_IndexManager_kokkos.hpp"
65 #include "MueLu_AggregationStructuredAlgorithm_kokkos.hpp"
71 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
72 StructuredAggregationFactory_kokkos<LocalOrdinal, GlobalOrdinal, Node>::
73 StructuredAggregationFactory_kokkos() : bDefinitionPhase_(true) { }
75 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
76 RCP<const ParameterList> StructuredAggregationFactory_kokkos<LocalOrdinal, GlobalOrdinal, Node>::
77 GetValidParameterList()
const {
78 RCP<ParameterList> validParamList =
rcp(
new ParameterList());
80 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name))
83 SET_VALID_ENTRY(
"aggregation: error on nodes with no on-rank neighbors");
85 #undef SET_VALID_ENTRY
88 validParamList->set<std::string> (
"aggregation: output type",
"Aggregates",
89 "Type of object holding the aggregation data: Aggregtes or CrsGraph");
90 validParamList->set<std::string> (
"aggregation: coarsening rate",
"{3}",
91 "Coarsening rate per spatial dimensions");
92 validParamList->set<
int> (
"aggregation: coarsening order", 0,
93 "The interpolation order used to construct grid transfer operators based off these aggregates.");
95 validParamList->set<RCP<const FactoryBase> >(
"Graph", Teuchos::null,
96 "Graph of the matrix after amalgamation but without dropping.");
97 validParamList->set<RCP<const FactoryBase> >(
"DofsPerNode", Teuchos::null,
98 "Number of degrees of freedom per mesh node, provided by the coalsce drop factory.");
99 validParamList->set<RCP<const FactoryBase> >(
"numDimensions", Teuchos::null,
100 "Number of spatial dimension provided by CoordinatesTransferFactory.");
101 validParamList->set<RCP<const FactoryBase> >(
"lNodesPerDim", Teuchos::null,
102 "Number of nodes per spatial dimmension provided by CoordinatesTransferFactory.");
104 return validParamList;
107 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
108 void StructuredAggregationFactory_kokkos<LocalOrdinal, GlobalOrdinal, Node>::
109 DeclareInput(Level& currentLevel)
const {
110 Input(currentLevel,
"Graph");
111 Input(currentLevel,
"DofsPerNode");
114 if(currentLevel.GetLevelID() == 0) {
115 if(currentLevel.IsAvailable(
"numDimensions", NoFactory::get())) {
116 currentLevel.DeclareInput(
"numDimensions", NoFactory::get(),
this);
119 Exceptions::RuntimeError,
120 "numDimensions was not provided by the user on level0!");
122 if(currentLevel.IsAvailable(
"lNodesPerDim", NoFactory::get())) {
123 currentLevel.DeclareInput(
"lNodesPerDim", NoFactory::get(),
this);
126 Exceptions::RuntimeError,
127 "lNodesPerDim was not provided by the user on level0!");
130 Input(currentLevel,
"lNodesPerDim");
131 Input(currentLevel,
"numDimensions");
135 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
136 void StructuredAggregationFactory_kokkos<LocalOrdinal, GlobalOrdinal, Node>::
137 Build(Level ¤tLevel)
const {
138 FactoryMonitor m(*
this,
"Build", currentLevel);
140 RCP<Teuchos::FancyOStream> out;
141 if(
const char* dbg = std::getenv(
"MUELU_STRUCTUREDAGGREGATION_DEBUG")) {
142 out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
143 out->setShowAllFrontMatter(
false).setShowProcRank(
true);
148 typedef typename LWGraph_kokkos::local_graph_type::device_type::execution_space execution_space;
149 typedef typename LWGraph_kokkos::local_graph_type::device_type::memory_space memory_space;
151 *out <<
"Entering structured aggregation" << std::endl;
153 ParameterList pL = GetParameterList();
154 bDefinitionPhase_ =
false;
157 RCP<const LWGraph_kokkos> graph = Get<RCP<LWGraph_kokkos> >(currentLevel,
"Graph");
158 RCP<const Map> fineMap = graph->GetDomainMap();
159 const int myRank = fineMap->getComm()->getRank();
160 const LO dofsPerNode = Get<LO>(currentLevel,
"DofsPerNode");
164 const int interpolationOrder = pL.get<
int>(
"aggregation: coarsening order");
165 std::string outputType = pL.get<std::string>(
"aggregation: output type");
166 const bool outputAggregates = (outputType ==
"Aggregates" ?
true :
false);
167 Array<LO> lFineNodesPerDir(3);
169 if(currentLevel.GetLevelID() == 0) {
171 lFineNodesPerDir = currentLevel.Get<Array<LO> >(
"lNodesPerDim", NoFactory::get());
172 numDimensions = currentLevel.Get<
int>(
"numDimensions", NoFactory::get());
175 lFineNodesPerDir = Get<Array<LO> >(currentLevel,
"lNodesPerDim");
176 numDimensions = Get<int>(currentLevel,
"numDimensions");
181 for(
int dim = 0; dim < 3; ++dim) {
182 if(dim >= numDimensions) {
183 lFineNodesPerDir[dim] = 1;
188 std::string coarseningRate = pL.get<std::string>(
"aggregation: coarsening rate");
191 coarseRate = Teuchos::fromStringToArray<LO>(coarseningRate);
193 GetOStream(
Errors,-1) <<
" *** \"aggregation: coarsening rate\" must be a string convertible into an array! *** "
198 Exceptions::RuntimeError,
199 "\"aggregation: coarsening rate\" must have at least as many"
200 " components as the number of spatial dimensions in the problem.");
203 RCP<IndexManager_kokkos> geoData =
rcp(
new IndexManager_kokkos(numDimensions,
204 interpolationOrder, myRank,
208 *out <<
"The index manager has now been built" << std::endl;
210 !=
static_cast<size_t>(geoData->getNumLocalFineNodes()),
211 Exceptions::RuntimeError,
212 "The local number of elements in the graph's map is not equal to "
213 "the number of nodes given by: lNodesPerDim!");
217 RCP<AggregationStructuredAlgorithm_kokkos> myStructuredAlgorithm
218 =
rcp(
new AggregationStructuredAlgorithm_kokkos());
220 if(interpolationOrder == 0 && outputAggregates){
221 RCP<Aggregates_kokkos> aggregates =
rcp(
new Aggregates_kokkos(graph->GetDomainMap()));
222 aggregates->setObjectLabel(
"ST");
223 aggregates->SetIndexManager(geoData);
224 aggregates->AggregatesCrossProcessors(
false);
225 aggregates->SetNumAggregates(geoData->getNumCoarseNodes());
227 LO numNonAggregatedNodes = geoData->getNumLocalFineNodes();
231 KOKKOS_LAMBDA(
const LO nodeIdx) {aggStat(nodeIdx) =
READY;});
233 myStructuredAlgorithm->BuildAggregates(pL, *graph, *aggregates, aggStat,
234 numNonAggregatedNodes);
236 *out <<
"numNonAggregatedNodes: " << numNonAggregatedNodes << std::endl;
239 "MueLu::StructuredAggregationFactory::Build: Leftover nodes found! Error!");
240 aggregates->ComputeAggregateSizes(
true);
241 GetOStream(
Statistics1) << aggregates->description() << std::endl;
242 Set(currentLevel,
"Aggregates", aggregates);
246 RCP<CrsGraph> myGraph;
247 myStructuredAlgorithm->BuildGraph(*graph, geoData, dofsPerNode, myGraph);
248 Set(currentLevel,
"prolongatorGraph", myGraph);
251 Set(currentLevel,
"lCoarseNodesPerDim", geoData->getCoarseNodesPerDirArray());
252 Set(currentLevel,
"indexManager", geoData);
253 Set(currentLevel,
"interpolationOrder", interpolationOrder);
254 Set(currentLevel,
"numDimensions", numDimensions);
260 #endif // HAVE_MUELU_KOKKOS_REFACTOR
void parallel_for(const ExecPolicy &policy, const FunctorType &functor, const std::string &str="", typename Impl::enable_if< Kokkos::Impl::is_execution_policy< ExecPolicy >::value >::type *=0)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
#define SET_VALID_ENTRY(name)