93 Kokkos::View<unsigned*, typename LWGraph_kokkos::device_type>& aggStat,
94 LO& numNonAggregatedNodes)
const {
96 bool error_on_isolated = params.get<
bool>(
"aggregation: error on nodes with no on-rank neighbors");
97 bool makeNonAdjAggs = params.get<
bool>(
"aggregation: phase3 avoid singletons");
99 const LO numRows = graph.GetNodeNumVertices();
100 const int myRank = graph.GetComm()->getRank();
102 auto vertex2AggId = aggregates.GetVertex2AggId()->getDeviceLocalView(Xpetra::Access::ReadWrite);
103 auto procWinner = aggregates.GetProcWinner() ->getDeviceLocalView(Xpetra::Access::ReadWrite);
104 auto colors = aggregates.GetGraphColors();
105 const LO numColors = aggregates.GetGraphNumColors();
107 auto lclLWGraph = graph.getLocalLWGraph();
109 Kokkos::View<LO, device_type> numAggregates(
"numAggregates");
110 Kokkos::deep_copy(numAggregates, aggregates.GetNumAggregates());
112 Kokkos::View<unsigned*, device_type> aggStatOld(
"Initial aggregation status", aggStat.extent(0));
113 Kokkos::deep_copy(aggStatOld, aggStat);
114 Kokkos::View<LO, device_type> numNonAggregated(
"numNonAggregated");
115 Kokkos::deep_copy(numNonAggregated, numNonAggregatedNodes);
116 for(
int color = 1; color < numColors + 1; ++color) {
117 Kokkos::parallel_for(
"Aggregation Phase 3: aggregates clean-up",
118 Kokkos::RangePolicy<execution_space>(0, numRows),
119 KOKKOS_LAMBDA(
const LO nodeIdx) {
121 if( (colors(nodeIdx) != color) ||
123 (aggStatOld(nodeIdx) ==
IGNORED) ){
return; }
126 auto neighbors = lclLWGraph.getNeighborVertices(nodeIdx);
131 bool isNewAggregate =
false;
132 for(
int neigh = 0; neigh < neighbors.length; ++neigh) {
133 neighIdx = neighbors(neigh);
135 if((neighIdx != nodeIdx) &&
136 lclLWGraph.isLocalNeighborVertex(neighIdx) &&
137 (aggStatOld(neighIdx) ==
READY)) {
138 isNewAggregate =
true;
147 const LO aggId = Kokkos::atomic_fetch_add(&numAggregates(), 1);
149 procWinner(nodeIdx, 0) = myRank;
150 vertex2AggId(nodeIdx, 0) = aggId;
152 Kokkos::atomic_decrement(&numNonAggregated());
153 for(
int neigh = 0; neigh < neighbors.length; ++neigh) {
154 neighIdx = neighbors(neigh);
155 if((neighIdx != nodeIdx) &&
156 lclLWGraph.isLocalNeighborVertex(neighIdx) &&
157 (aggStatOld(neighIdx) ==
READY)) {
159 procWinner(neighIdx, 0) = myRank;
160 vertex2AggId(neighIdx, 0) = aggId;
161 Kokkos::atomic_decrement(&numNonAggregated());
169 for(
int neigh = 0; neigh < neighbors.length; ++neigh) {
170 neighIdx = neighbors(neigh);
171 if (lclLWGraph.isLocalNeighborVertex(neighIdx) &&
174 procWinner(nodeIdx, 0) = myRank;
175 vertex2AggId(nodeIdx, 0) = vertex2AggId(neighIdx, 0);
176 Kokkos::atomic_decrement(&numNonAggregated());
184 for(LO otherNodeIdx = 0; otherNodeIdx < numRows; ++otherNodeIdx) {
185 if((otherNodeIdx != nodeIdx) &&
188 procWinner(nodeIdx, 0) = myRank;
189 vertex2AggId(nodeIdx, 0) = vertex2AggId(otherNodeIdx, 0);
190 Kokkos::atomic_decrement(&numNonAggregated());
198 if(!error_on_isolated) {
199 const LO aggId = Kokkos::atomic_fetch_add(&numAggregates(), 1);
201 procWinner(nodeIdx, 0) = myRank;
202 vertex2AggId(nodeIdx, 0) = aggId;
203 Kokkos::atomic_decrement(&numNonAggregated());
209 Kokkos::deep_copy(aggStatOld, aggStat);
212 auto numNonAggregated_h = Kokkos::create_mirror_view(numNonAggregated);
213 Kokkos::deep_copy(numNonAggregated_h, numNonAggregated);
214 numNonAggregatedNodes = numNonAggregated_h();
215 if( (error_on_isolated) && (numNonAggregatedNodes > 0) ) {
217 std::ostringstream oss;
218 oss<<
"MueLu::AggregationPhase3Algorithm::BuildAggregates: MueLu has detected a non-Dirichlet node that has no on-rank neighbors and is terminating (by user request). "<<std::endl;
219 oss<<
"If this error is being generated at level 0, this is due to an initial partitioning problem in your matrix."<<std::endl;
220 oss<<
"If this error is being generated at any other level, try turning on repartitioning, which may fix this problem."<<std::endl;
225 auto numAggregates_h = Kokkos::create_mirror_view(numAggregates);
226 Kokkos::deep_copy(numAggregates_h, numAggregates);
227 aggregates.SetNumAggregates(numAggregates_h());