46 #ifndef MUELU_AGGREGATIONSTRUCTUREDALGORITHM_DEF_HPP_
47 #define MUELU_AGGREGATIONSTRUCTUREDALGORITHM_DEF_HPP_
49 #include <Teuchos_Comm.hpp>
50 #include <Teuchos_CommHelpers.hpp>
52 #include <Xpetra_MapFactory.hpp>
53 #include <Xpetra_Map.hpp>
59 #include "MueLu_LWGraph.hpp"
60 #include "MueLu_LWGraph_kokkos.hpp"
61 #include "MueLu_Aggregates.hpp"
62 #include "MueLu_IndexManager.hpp"
65 #include "MueLu_IndexManager_kokkos.hpp"
69 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
74 LO& numNonAggregatedNodes)
const {
75 Monitor m(*
this,
"BuildAggregates");
78 if (
const char* dbg = std::getenv(
"MUELU_STRUCTUREDALGORITHM_DEBUG")) {
79 out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
94 *out <<
"Extract data for ghosted nodes" << std::endl;
96 ghostedCoarseNodeCoarsePIDs, ghostedCoarseNodeCoarseGIDs);
100 LO ghostedCoarseNodeCoarseLID, aggId;
101 *out <<
"Loop over fine nodes and assign them to an aggregate and a rank" << std::endl;
106 for (
int dim = 0; dim < 3; ++dim) {
117 if (rem > (rate / 2)) {
127 ghostedCoarseNodeCoarseLID);
129 aggId = ghostedCoarseNodeCoarseLIDs[ghostedCoarseNodeCoarseLID];
130 vertex2AggId[nodeIdx] = aggId;
131 procWinner[nodeIdx] = ghostedCoarseNodeCoarsePIDs[ghostedCoarseNodeCoarseLID];
133 --numNonAggregatedNodes;
138 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
143 Monitor m(*
this,
"BuildGraphP");
146 if (
const char* dbg = std::getenv(
"MUELU_STRUCTUREDALGORITHM_DEBUG")) {
147 out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
156 int numInterpolationPoints = 0;
158 numInterpolationPoints = 1;
163 *out <<
"numInterpolationPoints=" << numInterpolationPoints << std::endl;
172 *out <<
"Compute prolongatorGraph data" << std::endl;
174 ComputeGraphDataConstant(graph, geoData, dofsPerNode, numInterpolationPoints,
175 nnzOnRow, rowPtr, colIndex);
177 ComputeGraphDataLinear(graph, geoData, dofsPerNode, numInterpolationPoints,
178 nnzOnRow, rowPtr, colIndex);
184 *out <<
"Compute domain and column maps of the CrsGraph" << std::endl;
186 *out <<
"Extract data for ghosted nodes" << std::endl;
191 ghostedCoarseNodeCoarsePIDs, ghostedCoarseNodeCoarseGIDs);
197 ghostedCoarseNodeCoarseGIDs(),
201 LO coarseNodeIdx = 0;
202 Array<GO> coarseNodeCoarseGIDs, coarseNodeFineGIDs;
204 for (
LO nodeIdx = 0; nodeIdx < ghostedCoarseNodeCoarseGIDs.
size(); ++nodeIdx) {
205 if (ghostedCoarseNodeCoarsePIDs[nodeIdx] == colMap->getComm()->getRank()) {
206 coarseNodeCoarseGIDs[coarseNodeIdx] = ghostedCoarseNodeCoarseGIDs[nodeIdx];
210 domainMap = MapFactory::Build(graph.
GetDomainMap()->lib(),
212 coarseNodeCoarseGIDs(),
215 coarseCoordinatesMap = MapFactory::Build(graph.
GetDomainMap()->lib(),
217 coarseNodeCoarseGIDs(),
220 coarseCoordinatesFineMap = MapFactory::Build(graph.
GetDomainMap()->lib(),
222 coarseNodeFineGIDs(),
238 coarseCoordinatesMap = MapFactory::Build(graph.
GetDomainMap()->lib(),
243 coarseCoordinatesFineMap = MapFactory::Build(graph.
GetDomainMap()->lib(),
245 coarseNodeFineGIDs(),
250 *out <<
"Call constructor of CrsGraph" << std::endl;
251 myGraph = CrsGraphFactory::Build(rowMap,
255 *out <<
"Fill CrsGraph" << std::endl;
258 for (
LO dof = 0; dof < dofsPerNode; ++dof) {
259 rowIdx = nodeIdx * dofsPerNode + dof;
260 myGraph->insertLocalIndices(rowIdx, colIndex(rowPtr[rowIdx], nnzOnRow[rowIdx]));
264 *out <<
"Call fillComplete on CrsGraph" << std::endl;
265 myGraph->fillComplete(domainMap, rowMap);
266 *out <<
"Prolongator CrsGraph computed" << std::endl;
270 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
273 const LO dofsPerNode,
const int ,
277 if (
const char* dbg = std::getenv(
"MUELU_STRUCTUREDALGORITHM_DEBUG")) {
278 out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
288 ghostedCoarseNodeCoarsePIDs, ghostedCoarseNodeCoarseGIDs);
290 LO ghostedCoarseNodeCoarseLID, rem, rate;
296 for (
int dim = 0; dim < 3; ++dim) {
307 if (rem > (rate / 2)) {
317 ghostedCoarseNodeCoarseLID);
319 for (
LO dof = 0; dof < dofsPerNode; ++dof) {
320 nnzOnRow[nodeIdx * dofsPerNode + dof] = 1;
321 rowPtr[nodeIdx * dofsPerNode + dof + 1] = rowPtr[nodeIdx * dofsPerNode + dof] + 1;
322 colIndex[rowPtr[nodeIdx * dofsPerNode + dof]] =
323 ghostedCoarseNodeCoarseLIDs[ghostedCoarseNodeCoarseLID] * dofsPerNode + dof;
329 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
332 const LO dofsPerNode,
const int numInterpolationPoints,
336 if (
const char* dbg = std::getenv(
"MUELU_STRUCTUREDALGORITHM_DEBUG")) {
337 out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
348 const LO coarsePointOffset[8][3] = {{0, 0, 0}, {1, 0, 0}, {0, 1, 0}, {1, 1, 0}, {0, 0, 1}, {1, 0, 1}, {0, 1, 1}, {1, 1, 1}};
353 for (
int dim = 0; dim < numDimensions; dim++) {
369 bool allCoarse =
true;
371 for (
int dim = 0; dim < numDimensions; ++dim) {
372 isCoarse[dim] =
false;
373 if (ijkRem[dim] == 0)
374 isCoarse[dim] =
true;
379 isCoarse[dim] =
true;
382 isCoarse[dim] =
true;
389 LO rowIdx = 0, colIdx = 0;
391 for (
LO dof = 0; dof < dofsPerNode; ++dof) {
392 rowIdx = nodeIdx * dofsPerNode + dof;
393 nnzOnRow[rowIdx] = 1;
394 rowPtr[rowIdx + 1] = rowPtr[rowIdx] + 1;
398 colIndex[rowPtr[rowIdx]] = colIdx * dofsPerNode + dof;
402 for (
int dim = 0; dim < numDimensions; ++dim) {
407 for (
LO dof = 0; dof < dofsPerNode; ++dof) {
409 rowIdx = nodeIdx * dofsPerNode + dof;
410 nnzOnRow[rowIdx] = Teuchos::as<size_t>(numInterpolationPoints);
411 rowPtr[rowIdx + 1] = rowPtr[rowIdx] + Teuchos::as<LO>(numInterpolationPoints);
413 for (
LO interpIdx = 0; interpIdx < numInterpolationPoints; ++interpIdx) {
415 coarseIdx[1] + coarsePointOffset[interpIdx][1],
416 coarseIdx[2] + coarsePointOffset[interpIdx][2],
418 colIndex[rowPtr[rowIdx] + interpIdx] = colIdx * dofsPerNode + dof;
425 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
430 LO& numNonAggregatedNodes)
const {
431 Monitor m(*
this,
"BuildAggregates");
434 if (
const char* dbg = std::getenv(
"MUELU_STRUCTUREDALGORITHM_DEBUG")) {
435 out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
447 *out <<
"Loop over fine nodes and assign them to an aggregate and a rank" << std::endl;
448 LO numAggregatedNodes;
454 Kokkos::parallel_reduce(
"StructuredAggregation: fill aggregates data",
455 Kokkos::RangePolicy<execution_space>(0, numLocalFineNodes),
459 *out <<
"numCoarseNodes= " << numCoarseNodes
460 <<
", numAggregatedNodes= " << numAggregatedNodes << std::endl;
461 numNonAggregatedNodes = numNonAggregatedNodes - numAggregatedNodes;
465 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
469 Monitor m(*
this,
"BuildGraphP");
472 if (
const char* dbg = std::getenv(
"MUELU_STRUCTUREDALGORITHM_DEBUG")) {
473 out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
480 int numInterpolationPoints = 0;
482 numInterpolationPoints = 1;
487 *out <<
"numInterpolationPoints=" << numInterpolationPoints << std::endl;
491 const LO numNnzEntries = dofsPerNode * (numCoarseNodes + numInterpolationPoints * (numLocalFineNodes - numCoarseNodes));
494 entries_type colIndex(
"Prolongator graph, colIndices", numNnzEntries);
496 *out <<
"Compute prolongatorGraph data" << std::endl;
506 Kokkos::parallel_for(
"Structured Aggregation: compute loca graph data",
507 Kokkos::RangePolicy<execution_space>(0, numLocalFineNodes),
516 numInterpolationPoints,
521 Kokkos::parallel_scan(
"Structured Aggregation: compute rowPtr for prolongator graph",
522 Kokkos::RangePolicy<execution_space>(0, numLocalFineNodes + 1),
529 numInterpolationPoints,
536 Kokkos::parallel_for(
"Structured Aggregation: compute loca graph data",
537 Kokkos::RangePolicy<execution_space>(0, numLocalFineNodes),
545 *out <<
"Compute domain and column maps of the CrsGraph" << std::endl;
553 myGraph = CrsGraphFactory::Build(myLocalGraph, graph.
GetDomainMap(), colMap,
558 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
562 Kokkos::View<unsigned*, device_type> aggStat,
568 , vertex2AggID_(vertex2AggID)
569 , procWinner_(procWinner) {}
571 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
576 LO coarseNodeCoarseLID;
577 LO nodeFineTuple[3], coarseIdx[3];
578 auto coarseRate = geoData_.getCoarseningRates();
579 auto endRate = geoData_.getCoarseningEndRates();
580 auto lFineNodesPerDir = geoData_.getLocalFineNodesPerDir();
582 geoData_.getFineLID2FineTuple(nodeIdx, nodeFineTuple);
584 for (
int dim = 0; dim < 3; ++dim) {
585 coarseIdx[dim] = nodeFineTuple[dim] / coarseRate(dim);
586 rem = nodeFineTuple[dim] % coarseRate(dim);
587 rate = (nodeFineTuple[dim] < lFineNodesPerDir(dim) - endRate(dim)) ? coarseRate(dim) : endRate(dim);
588 if (rem > (rate / 2)) {
593 geoData_.getCoarseTuple2CoarseLID(coarseIdx[0], coarseIdx[1], coarseIdx[2],
594 coarseNodeCoarseLID);
596 vertex2AggID_(nodeIdx, 0) = coarseNodeCoarseLID;
597 procWinner_(nodeIdx, 0) = myRank_;
599 ++lNumAggregatedNodes;
603 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
607 const LO NumGhostedNodes,
608 const LO dofsPerNode,
615 , numGhostedNodes_(NumGhostedNodes)
616 , dofsPerNode_(dofsPerNode)
617 , coarseRate_(coarseRate)
619 , lFineNodesPerDir_(lFineNodesPerDir)
621 , colIndex_(colIndex) {
624 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
627 LO nodeFineTuple[3] = {0, 0, 0};
628 LO nodeCoarseTuple[3] = {0, 0, 0};
631 geoData_.getFineLID2FineTuple(nodeIdx, nodeFineTuple);
635 LO rem, rate, coarseNodeCoarseLID;
636 for (
int dim = 0; dim < 3; ++dim) {
637 nodeCoarseTuple[dim] = nodeFineTuple[dim] / coarseRate_(dim);
638 rem = nodeFineTuple[dim] % coarseRate_(dim);
639 if (nodeFineTuple[dim] < (lFineNodesPerDir_(dim) - endRate_(dim))) {
640 rate = coarseRate_(dim);
642 rate = endRate_(dim);
644 if (rem > (rate / 2)) {
645 ++nodeCoarseTuple[dim];
650 geoData_.getCoarseTuple2CoarseLID(nodeCoarseTuple[0], nodeCoarseTuple[1], nodeCoarseTuple[2],
651 coarseNodeCoarseLID);
654 for (
LO dof = 0; dof < dofsPerNode_; ++dof) {
655 rowPtr_(nodeIdx * dofsPerNode_ + dof + 1) = nodeIdx * dofsPerNode_ + dof + 1;
656 colIndex_(nodeIdx * dofsPerNode_ + dof) = coarseNodeCoarseLID * dofsPerNode_ + dof;
661 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
664 const LO dofsPerNode,
665 const int numInterpolationPoints,
666 const LO numLocalRows,
671 , dofsPerNode_(dofsPerNode)
672 , numInterpolationPoints_(numInterpolationPoints)
673 , numLocalRows_(numLocalRows)
674 , coarseRate_(coarseRate)
675 , lFineNodesPerDir_(lFineNodesPerDir)
678 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
687 rowPtr_(rowIdx) = update;
689 if (rowIdx < numLocalRows_) {
690 LO nodeIdx = rowIdx / dofsPerNode_;
691 bool allCoarse =
true;
692 LO nodeFineTuple[3] = {0, 0, 0};
693 geoData_.getFineLID2FineTuple(nodeIdx, nodeFineTuple);
694 for (
int dim = 0; dim < 3; ++dim) {
695 const LO rem = nodeFineTuple[dim] % coarseRate_(dim);
698 allCoarse = (allCoarse && ((rem == 0) || (nodeFineTuple[dim] == lFineNodesPerDir_(dim) - 1)));
700 update += (allCoarse ? 1 : numInterpolationPoints_);
704 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
707 const int numDimensions,
708 const LO numGhostedNodes,
709 const LO dofsPerNode,
710 const int numInterpolationPoints,
718 , numDimensions_(numDimensions)
719 , numGhostedNodes_(numGhostedNodes)
720 , dofsPerNode_(dofsPerNode)
721 , numInterpolationPoints_(numInterpolationPoints)
722 , coarseRate_(coarseRate)
724 , lFineNodesPerDir_(lFineNodesPerDir)
725 , ghostedNodesPerDir_(ghostedNodesPerDir)
727 , colIndex_(colIndex) {
730 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
733 LO nodeFineTuple[3] = {0, 0, 0};
734 LO nodeCoarseTuple[3] = {0, 0, 0};
737 geoData_.getFineLID2FineTuple(nodeIdx, nodeFineTuple);
739 LO coarseNodeCoarseLID;
740 bool allCoarse =
false;
741 for (
int dim = 0; dim < 3; ++dim) {
742 nodeCoarseTuple[dim] = nodeFineTuple[dim] / coarseRate_(dim);
744 if (rowPtr_(nodeIdx + 1) == rowPtr_(nodeIdx) + 1) {
748 geoData_.getCoarseTuple2CoarseLID(nodeCoarseTuple[0], nodeCoarseTuple[1], nodeCoarseTuple[2],
749 coarseNodeCoarseLID);
753 for (
LO dof = 0; dof < dofsPerNode_; ++dof) {
754 colIndex_(rowPtr_(nodeIdx * dofsPerNode_ + dof)) = coarseNodeCoarseLID * dofsPerNode_ + dof;
757 for (
int dim = 0; dim < numDimensions_; ++dim) {
758 if (nodeCoarseTuple[dim] == ghostedNodesPerDir_(dim) - 1) {
759 --nodeCoarseTuple[dim];
765 for (
LO dof = 0; dof < dofsPerNode_; ++dof) {
766 geoData_.getCoarseTuple2CoarseLID(nodeCoarseTuple[0], nodeCoarseTuple[1], nodeCoarseTuple[2], colIndex_(rowPtr_(nodeIdx * dofsPerNode_ + dof) + 0));
767 geoData_.getCoarseTuple2CoarseLID(nodeCoarseTuple[0] + 1, nodeCoarseTuple[1], nodeCoarseTuple[2], colIndex_(rowPtr_(nodeIdx * dofsPerNode_ + dof) + 1));
768 if (numDimensions_ > 1) {
769 geoData_.getCoarseTuple2CoarseLID(nodeCoarseTuple[0], nodeCoarseTuple[1] + 1, nodeCoarseTuple[2], colIndex_(rowPtr_(nodeIdx * dofsPerNode_ + dof) + 2));
770 geoData_.getCoarseTuple2CoarseLID(nodeCoarseTuple[0] + 1, nodeCoarseTuple[1] + 1, nodeCoarseTuple[2], colIndex_(rowPtr_(nodeIdx * dofsPerNode_ + dof) + 3));
771 if (numDimensions_ > 2) {
772 geoData_.getCoarseTuple2CoarseLID(nodeCoarseTuple[0], nodeCoarseTuple[1], nodeCoarseTuple[2] + 1, colIndex_(rowPtr_(nodeIdx * dofsPerNode_ + dof) + 4));
773 geoData_.getCoarseTuple2CoarseLID(nodeCoarseTuple[0] + 1, nodeCoarseTuple[1], nodeCoarseTuple[2] + 1, colIndex_(rowPtr_(nodeIdx * dofsPerNode_ + dof) + 5));
774 geoData_.getCoarseTuple2CoarseLID(nodeCoarseTuple[0], nodeCoarseTuple[1] + 1, nodeCoarseTuple[2] + 1, colIndex_(rowPtr_(nodeIdx * dofsPerNode_ + dof) + 6));
775 geoData_.getCoarseTuple2CoarseLID(nodeCoarseTuple[0] + 1, nodeCoarseTuple[1] + 1, nodeCoarseTuple[2] + 1, colIndex_(rowPtr_(nodeIdx * dofsPerNode_ + dof) + 7));
Kokkos::View< unsigned *, typename LWGraphHostType::device_type > AggStatHostType
GO getStartIndex(int const dim) const
LO getLocalCoarseNodesInDir(const int dim) const
virtual void getFineNodeGhostedTuple(const LO myLID, LO &i, LO &j, LO &k) const =0
bool getMeshEdge(const int dir) const
Lightweight MueLu representation of a compressed row storage graph.
const RCP< LOVector > & GetProcWinner() const
Returns constant vector that maps local node IDs to owning processor IDs.
void ComputeGraphDataConstant(const LWGraph &graph, RCP< IndexManager > &geoData, const LO dofsPerNode, const int numInterpolationPoints, ArrayRCP< size_t > &nnzOnRow, Array< size_t > &rowPtr, Array< LO > &colIndex) const
Container class for aggregation information.
RCP< IndexManager > & GetIndexManager()
Get the index manager used by various aggregation algorithms. This has to be done by the aggregation ...
int getInterpolationOrder() const
LO getGhostedNodesInDir(const int dim) const
fillAggregatesFunctor(RCP< IndexManager_kokkos > geoData, const int myRank, Kokkos::View< unsigned *, device_type > aggStat, LOVectorView vertex2AggID, LOVectorView procWinner)
basic_FancyOStream & setShowProcRank(const bool showProcRank)
KOKKOS_INLINE_FUNCTION intTupleView getCoarseningEndRates() const
KOKKOS_INLINE_FUNCTION void operator()(const LO nodeIdx) const
KOKKOS_INLINE_FUNCTION intTupleView getCoarseningRates() const
decltype(std::declval< LOVector >().getDeviceLocalView(Xpetra::Access::ReadWrite)) LOVectorView
int getInterpolationOrder() const
GO getStartGhostedCoarseNode(int const dim) const
int getCoarseningRate(const int dim) const
LO getNumCoarseNodes() const
KOKKOS_INLINE_FUNCTION void operator()(const LO nodeIdx) const
typename Kokkos::View< const LO[3], device_type > constLOTupleView
virtual void getGhostedNodesData(const RCP< const Map > fineMap, Array< LO > &ghostedNodeCoarseLIDs, Array< int > &ghostedNodeCoarsePIDs, Array< GO > &ghostedNodeCoarseGIDs) const =0
computeGraphRowPtrFunctor(RCP< IndexManager_kokkos > geoData, const LO dofsPerNode, const int numInterpolationPoints, const LO numLocalRows, constIntTupleView coarseRate, constLOTupleView lFineNodesPerDir, non_const_row_map_type rowPtr)
int getNumDimensions() const
LO getNumLocalFineNodes() const
computeGraphDataLinearFunctor(RCP< IndexManager_kokkos > geoData, const int numDimensions, const LO numGhostedNodes, const LO dofsPerNode, const int numInterpolationPoints, constIntTupleView coarseRate, constIntTupleView endRate, constLOTupleView lFineNodesPerDir, constLOTupleView ghostedNodesPerDir, non_const_row_map_type rowPtr, entries_type colIndex)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
LO getNumLocalFineNodes() const
void BuildGraphOnHost(const LWGraph &graph, RCP< IndexManager > &geoData, const LO dofsPerNode, RCP< CrsGraph > &myGraph, RCP< const Map > &coarseCoordinatesFineMap, RCP< const Map > &coarseCoordinatesMap) const
Local aggregation.
GO getNumGlobalCoarseNodes() const
virtual void getCoarseNodeGhostedLID(const LO i, const LO j, const LO k, LO &myLID) const =0
virtual void getCoarseNodesData(const RCP< const Map > fineCoordinatesMap, Array< GO > &coarseNodeCoarseGIDs, Array< GO > &coarseNodeFineGIDs) const =0
void BuildAggregates(const Teuchos::ParameterList ¶ms, const LWGraph_kokkos &graph, Aggregates &aggregates, typename AggregationAlgorithmBase< LocalOrdinal, GlobalOrdinal, Node >::AggStatType &aggStat, LO &numNonAggregatedNodes) const
Build aggregates object.
LO getNumLocalCoarseNodes() const
typename local_graph_type::row_map_type::non_const_type non_const_row_map_type
basic_FancyOStream & setShowAllFrontMatter(const bool showAllFrontMatter)
const RCP< LOMultiVector > & GetVertex2AggId() const
Returns constant vector that maps local node IDs to local aggregates IDs.
KOKKOS_INLINE_FUNCTION void operator()(const LO rowIdx, GO &update, const bool final) const
typename LWGraph_kokkos::local_graph_type local_graph_type
LO getLocalFineNodesInDir(const int dim) const
KOKKOS_INLINE_FUNCTION void operator()(const LO nodeIdx, LO &lNumAggregatedNodes) const
bool isAggregationCoupled() const
void BuildAggregatesNonKokkos(const Teuchos::ParameterList ¶ms, const LWGraph &graph, Aggregates &aggregates, typename AggregationAlgorithmBase< LocalOrdinal, GlobalOrdinal, Node >::AggStatHostType &aggStat, LO &numNonAggregatedNodes) const
Local aggregation.
LO getOffset(int const dim) const
Timer to be used in non-factories.
typename local_graph_type::entries_type entries_type
computeGraphDataConstantFunctor(RCP< IndexManager_kokkos > geoData, const LO numGhostedNodes, const LO dofsPerNode, constIntTupleView coarseRate, constIntTupleView endRate, constLOTupleView lFineNodesPerDir, non_const_row_map_type rowPtr, entries_type colIndex)
bool isSingleCoarsePoint() const
int getCoarseningEndRate(const int dim) const
const RCP< const Map > GetDomainMap() const
Lightweight MueLu representation of a compressed row storage graph.
void ComputeGraphDataLinear(const LWGraph &graph, RCP< IndexManager > &geoData, const LO dofsPerNode, const int numInterpolationPoints, ArrayRCP< size_t > &nnzOnRow, Array< size_t > &rowPtr, Array< LO > &colIndex) const
KOKKOS_INLINE_FUNCTION LOTupleView getCoarseNodesPerDir() const
KOKKOS_INLINE_FUNCTION LOTupleView getLocalFineNodesPerDir() const
void BuildGraph(const LWGraph_kokkos &graph, RCP< IndexManager_kokkos > &geoData, const LO dofsPerNode, RCP< CrsGraph > &myGraph) const
Build a CrsGraph instead of aggregates.
const RCP< const Teuchos::Comm< int > > GetComm() const
typename Kokkos::View< const int[3], device_type > constIntTupleView
Kokkos::View< unsigned *, typename LWGraphType::device_type > AggStatType
int getNumDimensions() const
RCP< IndexManager_kokkos > & GetIndexManagerKokkos()
Get the index manager used by structured aggregation algorithms. This has to be done by the aggregati...