46 #ifndef MUELU_PRESERVEDIRICHLETAGGREGATIONALGORITHM_DEF_HPP_
47 #define MUELU_PRESERVEDIRICHLETAGGREGATIONALGORITHM_DEF_HPP_
49 #include <Teuchos_Comm.hpp>
50 #include <Teuchos_CommHelpers.hpp>
56 #include "MueLu_LWGraph.hpp"
57 #include "MueLu_Aggregates.hpp"
63 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
65 Monitor m(*
this,
"BuildAggregatesNonKokkos");
67 bool preserve = params.
get<
bool>(
"aggregation: preserve Dirichlet points");
70 const int myRank = graph.
GetComm()->getRank();
77 for (
LO i = 0; i < numRows; i++)
80 numNonAggregatedNodes--;
85 vertex2AggId[i] = numLocalAggregates++;
86 procWinner[i] = myRank;
94 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
100 LO& numNonAggregatedNodes)
const {
104 Monitor m(*
this,
"BuildAggregates");
109 const bool preserve = params.
get<
bool>(
"aggregation: preserve Dirichlet points");
113 const int myRank = graph.
GetComm()->getRank();
116 auto vertex2AggId = aggregates.
GetVertex2AggId()->getDeviceLocalView(Xpetra::Access::ReadWrite);
117 auto procWinner = aggregates.
GetProcWinner()->getDeviceLocalView(Xpetra::Access::ReadWrite);
120 Kokkos::View<LO, device_type> aggCount(
"aggCount");
124 Kokkos::parallel_for(
125 "MueLu - PreserveDirichlet: tagging ignored nodes",
126 Kokkos::RangePolicy<local_ordinal_type, execution_space>(0, numNodes),
127 KOKKOS_LAMBDA(
const local_ordinal_type nodeIdx) {
130 const LO aggIdx = Kokkos::atomic_fetch_add(&aggCount(), 1);
135 vertex2AggId(nodeIdx, 0) = aggIdx;
136 procWinner(nodeIdx, 0) = myRank;
140 typename Kokkos::View<LO, device_type>::HostMirror aggCount_h = Kokkos::create_mirror_view(aggCount);
141 Kokkos::deep_copy(aggCount_h, aggCount);
Kokkos::View< unsigned *, typename LWGraphHostType::device_type > AggStatHostType
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
T & get(const std::string &name, T def_value)
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
void BuildAggregates(const Teuchos::ParameterList ¶ms, const LWGraph_kokkos &graph, Aggregates &aggregates, typename AggregationAlgorithmBase< LocalOrdinal, GlobalOrdinal, Node >::AggStatType &aggStat, LO &numNonAggregatedNodes) const
LocalOrdinal local_ordinal_type
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.
void BuildAggregatesNonKokkos(const Teuchos::ParameterList ¶ms, const LWGraph &graph, Aggregates &aggregates, typename AggregationAlgorithmBase< LocalOrdinal, GlobalOrdinal, Node >::AggStatHostType &aggStat, LO &numNonAggregatedNodes) const
Local aggregation.
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.