53 #ifndef MUELU_ONEPTAGGREGATIONALGORITHM_DEF_HPP_
54 #define MUELU_ONEPTAGGREGATIONALGORITHM_DEF_HPP_
56 #include <Teuchos_Comm.hpp>
57 #include <Teuchos_CommHelpers.hpp>
63 #include "MueLu_LWGraph.hpp"
64 #include "MueLu_Aggregates.hpp"
70 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
74 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
76 Monitor m(*
this,
"BuildAggregatesNonKokkos");
79 const int myRank = graph.
GetComm()->getRank();
90 while (iNode1 < nRows) {
91 if (aggStat[iNode1] ==
ONEPT) {
93 std::vector<int> aggList;
94 aggList.push_back(iNode1);
95 int aggIndex = nLocalAggregates++;
97 for (
size_t k = 0; k < aggList.size(); k++) {
99 vertex2AggId[aggList[k]] = aggIndex;
100 procWinner[aggList[k]] = myRank;
102 numNonAggregatedNodes -= aggList.
size();
112 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
118 LO& numNonAggregatedNodes)
const {
122 Monitor m(*
this,
"BuildAggregates");
124 typename Kokkos::View<unsigned*, device_type>::HostMirror aggstatHost = Kokkos::create_mirror(aggstat);
125 Kokkos::deep_copy(aggstatHost, aggstat);
126 std::vector<unsigned> aggStat;
127 aggStat.resize(aggstatHost.extent(0));
128 for (
size_t idx = 0; idx < aggstatHost.extent(0); ++idx) {
129 aggStat[idx] = aggstatHost(idx);
133 const int myRank = graph.
GetComm()->getRank();
144 while (iNode1 < nRows) {
145 if (aggStat[iNode1] ==
ONEPT) {
147 std::vector<int> aggList;
148 aggList.push_back(iNode1);
149 int aggIndex = nLocalAggregates++;
152 for (
size_t k = 0; k < aggList.size(); k++) {
154 vertex2AggId[aggList[k]] = aggIndex;
155 procWinner[aggList[k]] = myRank;
157 numNonAggregatedNodes -= aggList.
size();
163 for (
size_t idx = 0; idx < aggstatHost.extent(0); ++idx) {
164 aggstatHost(idx) = aggStat[idx];
166 Kokkos::deep_copy(aggstat, aggstatHost);
Kokkos::View< unsigned *, typename LWGraphHostType::device_type > AggStatHostType
void BuildAggregatesNonKokkos(Teuchos::ParameterList const ¶ms, LWGraph const &graph, Aggregates &aggregates, typename AggregationAlgorithmBase< LocalOrdinal, GlobalOrdinal, Node >::AggStatHostType &aggStat, LO &numNonAggregatedNodes) const
Local aggregation.
MueLu::DefaultLocalOrdinal LocalOrdinal
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.
Container class for aggregation information.
KOKKOS_INLINE_FUNCTION LO GetNumAggregates() const
KOKKOS_INLINE_FUNCTION size_type GetNodeNumVertices() const
Return number of graph vertices.
void SetIsRoot(LO i, bool value=true)
Set root node information.
typename device_type::execution_space execution_space
OnePtAggregationAlgorithm(RCP< const FactoryBase > const &graphFact=Teuchos::null)
Constructor.
void BuildAggregates(Teuchos::ParameterList const ¶ms, LWGraph_kokkos const &graph, Aggregates &aggregates, typename AggregationAlgorithmBase< LocalOrdinal, GlobalOrdinal, Node >::AggStatType &aggStat, LO &numNonAggregatedNodes) const
const RCP< LOMultiVector > & GetVertex2AggId() const
Returns constant vector that maps local node IDs to local aggregates IDs.
Timer to be used in non-factories.
Lightweight MueLu representation of a compressed row storage graph.
typename std::conditional< OnHost, Kokkos::Device< Kokkos::Serial, Kokkos::HostSpace >, typename Node::device_type >::type device_type
const RCP< const Teuchos::Comm< int > > GetComm() const
Kokkos::View< unsigned *, typename LWGraphType::device_type > AggStatType
void SetNumAggregates(LO nAggregates)
Set number of local aggregates on current processor.