47#ifndef MUELU_CLASSICALMAPFACTORY_DEF_HPP_
48#define MUELU_CLASSICALMAPFACTORY_DEF_HPP_
52#include <Teuchos_Array.hpp>
53#include <Teuchos_ArrayRCP.hpp>
57#include <Teuchos_DefaultMpiComm.hpp>
60#include <Xpetra_Vector.hpp>
61#include <Xpetra_StridedMapFactory.hpp>
62#include <Xpetra_VectorFactory.hpp>
63#include <Xpetra_Import.hpp>
64#include <Xpetra_IO.hpp>
72#include "MueLu_Graph.hpp"
73#include "MueLu_LWGraph.hpp"
75#ifdef HAVE_MUELU_ZOLTAN2
77#include <Zoltan2_XpetraCrsGraphAdapter.hpp>
78#include <Zoltan2_ColoringProblem.hpp>
79#include <Zoltan2_ColoringSolution.hpp>
83#include "MueLu_LWGraph_kokkos.hpp"
84#include <KokkosGraph_Distance1ColorHandle.hpp>
85#include <KokkosGraph_Distance1Color.hpp>
89 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
92 RCP<ParameterList> validParamList = rcp(
new ParameterList());
93#define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name))
98 validParamList->set< RCP<const FactoryBase> >(
"A", Teuchos::null,
"Generating factory of the matrix A");
99 validParamList->set< RCP<const FactoryBase> >(
"UnAmalgamationInfo", Teuchos::null,
"Generating factory of UnAmalgamationInfo");
100 validParamList->set< RCP<const FactoryBase> >(
"Graph", null,
"Generating factory of the graph");
101 validParamList->set< RCP<const FactoryBase> >(
"Coloring Graph", null,
"Generating factory of the graph");
102 validParamList->set< RCP<const FactoryBase> >(
"DofsPerNode", null,
"Generating factory for variable \'DofsPerNode\', usually the same as for \'Graph\'");
104 return validParamList;
107 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
110 Input(currentLevel,
"A");
111 Input(currentLevel,
"UnAmalgamationInfo");
112 Input(currentLevel,
"Graph");
114 const ParameterList& pL = GetParameterList();
115 bool use_color_graph = pL.get<
bool>(
"aggregation: coloring: use color graph");
117 Input(currentLevel,
"Coloring Graph");
120 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
124 const ParameterList& pL = GetParameterList();
125 RCP<const Matrix> A = Get<RCP<Matrix> >(currentLevel,
"A");
127 RCP<const GraphBase> graph;
128 bool use_color_graph = pL.get<
bool>(
"aggregation: coloring: use color graph");
129 if(use_color_graph) graph = Get<RCP<GraphBase> >(currentLevel,
"Coloring Graph");
130 else graph = Get<RCP<GraphBase> >(currentLevel,
"Graph");
137 ArrayRCP<LO> myColors;
140 RCP<LocalOrdinalVector> fc_splitting;
141 std::string coloringAlgo = pL.get<std::string>(
"aggregation: coloring algorithm");
144#ifdef HAVE_MUELU_ZOLTAN2
145 int numProcs = A->getRowMap()->getComm()->getSize();
146 if(coloringAlgo!=
"file" && numProcs>1 && graph->GetDomainMap()->lib() == Xpetra::UseTpetra)
147 coloringAlgo=
"Zoltan2";
153 int rank = graph->GetDomainMap()->getComm()->getRank();
155 printf(
"[%d,%d] graph local size = %dx%d\n",rank,currentLevel.
GetLevelID(),(
int)graph->GetDomainMap()->getLocalNumElements(),(
int)graph->GetImportMap()->getLocalNumElements());
157 std::ofstream ofs(std::string(
"m_dropped_graph_") + std::to_string(currentLevel.
GetLevelID())+std::string(
"_") + std::to_string(rank) + std::string(
".dat"),std::ofstream::out);
158 RCP<Teuchos::FancyOStream> fancy = Teuchos::fancyOStream(Teuchos::rcpFromRef(ofs));
159 graph->print(*fancy,
Debug);
162 A->getRowMap()->getComm()->barrier();
170 if(coloringAlgo!=
"file" && graph->GetDomainMap()->lib() == Xpetra::UseEpetra)
174 if(coloringAlgo ==
"file") {
177 std::string map_file = std::string(
"map_fcsplitting_") + std::to_string(currentLevel.
GetLevelID()) + std::string(
".m");
178 std::string color_file = std::string(
"fcsplitting_") + std::to_string(currentLevel.
GetLevelID()) + std::string(
".m");
180 FILE * mapfile = fopen(map_file.c_str(),
"r");
181 using real_type =
typename Teuchos::ScalarTraits<SC>::magnitudeType;
182 using RealValuedMultiVector =
typename Xpetra::MultiVector<real_type,LO,GO,NO>;
183 RCP<RealValuedMultiVector> mv;
186 GetOStream(
Statistics1) <<
"Reading FC splitting from " << color_file <<
", using map file " << map_file <<
". On rank " << A->getRowMap()->getComm()->getRank() <<
" local size is " << A->getRowMap()->getLocalNumElements() << std::endl;
189 RCP<const Map> colorMap = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::ReadMap(map_file, A->getRowMap()->lib(), A->getRowMap()->getComm());
190 TEUCHOS_TEST_FOR_EXCEPTION(!colorMap->isCompatible(*A->getRowMap()),std::invalid_argument,
"Coloring on disk has incompatible map with A");
192 mv = Xpetra::IO<real_type, LocalOrdinal, GlobalOrdinal, Node>::ReadMultiVector(color_file,colorMap);
196 mv = Xpetra::IO<real_type, LocalOrdinal, GlobalOrdinal, Node>::ReadMultiVector(color_file,A->getRowMap());
198 TEUCHOS_TEST_FOR_EXCEPTION(mv.is_null(),std::invalid_argument,
"Coloring on disk cannot be read");
199 fc_splitting = LocalOrdinalVectorFactory::Build(A->getRowMap());
200 TEUCHOS_TEST_FOR_EXCEPTION(mv->getLocalLength() != fc_splitting->getLocalLength(),std::invalid_argument,
"Coloring map mismatch");
203 auto boundaryNodes = graph->GetBoundaryNodeMap();
204 ArrayRCP<const real_type> mv_data= mv->getData(0);
205 ArrayRCP<LO> fc_data= fc_splitting->getDataNonConst(0);
206 for(LO i=0; i<(LO)fc_data.size(); i++) {
208 fc_data[i] = DIRICHLET_PT;
210 fc_data[i] = Teuchos::as<LO>(mv_data[i]);
213#ifdef HAVE_MUELU_ZOLTAN2
214 else if(coloringAlgo.find(
"Zoltan2")!=std::string::npos && graph->GetDomainMap()->lib() == Xpetra::UseTpetra) {
216 DoDistributedGraphColoring(graph,myColors,numColors);
219 else if(coloringAlgo ==
"MIS" || graph->GetDomainMap()->lib() == Xpetra::UseTpetra) {
221; TEUCHOS_TEST_FOR_EXCEPTION(A->getRowMap()->getComm()->getSize() != 1, std::invalid_argument,
"MIS on more than 1 MPI rank is not supported");
222 DoMISNaive(*graph,myColors,numColors);
226 TEUCHOS_TEST_FOR_EXCEPTION(A->getRowMap()->getComm()->getSize() != 1, std::invalid_argument,
"KokkosKernels graph coloring on more than 1 MPI rank is not supported");
227 DoGraphColoring(*graph,myColors,numColors);
232 int rank = graph->GetDomainMap()->getComm()->getRank();
234 printf(
"[%d,%d] num colors %d\n",rank,currentLevel.
GetLevelID(),numColors);
236 std::ofstream ofs(std::string(
"m_colors_") + std::to_string(currentLevel.
GetLevelID())+std::string(
"_") + std::to_string(rank) + std::string(
".dat"),std::ofstream::out);
237 RCP<Teuchos::FancyOStream> fancy = Teuchos::fancyOStream(Teuchos::rcpFromRef(ofs));
238 *fancy << myColors();
241 A->getRowMap()->getComm()->barrier();
249 LO num_c_points = 0, num_d_points=0, num_f_points = 0;
250 if(fc_splitting.is_null()) {
252 auto boundaryNodes = graph->GetBoundaryNodeMap();
253 fc_splitting = LocalOrdinalVectorFactory::Build(A->getRowMap());
254 ArrayRCP<LO> myPointType = fc_splitting->getDataNonConst(0);
255 for(LO i=0; i<(LO)myColors.size(); i++) {
256 if(boundaryNodes[i]) {
257 myPointType[i] = DIRICHLET_PT;
260 else if ((LO)myColors[i] == 1) {
261 myPointType[i] = C_PT;
265 myPointType[i] = F_PT;
267 num_f_points = (LO)myColors.size() - num_d_points - num_c_points;
271 ArrayRCP<LO> myPointType = fc_splitting->getDataNonConst(0);
273 for(LO i=0; i<(LO)myPointType.size(); i++) {
274 if(myPointType[i] == DIRICHLET_PT)
276 else if (myPointType[i] == C_PT)
279 num_f_points = (LO)myPointType.size() - num_d_points - num_c_points;
285 GO l_counts[] = {(GO)num_c_points, (GO) num_f_points, (GO) num_d_points};
288 RCP<const Teuchos::Comm<int> > comm = A->getRowMap()->getComm();
289 Teuchos::reduceAll(*comm, Teuchos::REDUCE_SUM, 3, l_counts, g_counts);
290 GetOStream(
Statistics1) <<
"ClassicalMapFactory("<<coloringAlgo<<
"): C/F/D = "<<g_counts[0]<<
"/"<<g_counts[1]<<
"/"<<g_counts[2]<<std::endl;
295 RCP<const Map> coarseMap;
298 GenerateCoarseMap(*A->getRowMap(),num_c_points,coarseMap);
301 Set(currentLevel,
"FC Splitting",fc_splitting);
302 Set(currentLevel,
"CoarseMap", coarseMap);
307 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
309 GenerateCoarseMap(
const Map & fineMap, LO num_c_points, RCP<const Map> & coarseMap)
const {
312 std::vector<size_t> stridingInfo_(1);
314 GO domainGIDOffset = 0;
316 coarseMap = StridedMapFactory::Build(fineMap.lib(),
317 Teuchos::OrdinalTraits<Xpetra::global_size_t>::invalid(),
319 fineMap.getIndexBase(),
326 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
329 const ParameterList& pL = GetParameterList();
330 using graph_t =
typename LWGraph_kokkos::local_graph_type;
331 using KernelHandle = KokkosKernels::Experimental::
332 KokkosKernelsHandle<
typename graph_t::row_map_type::value_type,
333 typename graph_t::entries_type::value_type,
334 typename graph_t::entries_type::value_type,
335 typename graph_t::device_type::execution_space,
336 typename graph_t::device_type::memory_space,
337 typename graph_t::device_type::memory_space>;
341 kh.create_graph_coloring_handle();
344 auto coloringHandle = kh.get_graph_coloring_handle();
347 if(pL.get<
bool>(
"aggregation: deterministic") ==
true) {
348 coloringHandle->set_algorithm( KokkosGraph::COLORING_SERIAL );
350 }
else if(pL.get<std::string>(
"aggregation: coloring algorithm") ==
"serial") {
351 coloringHandle->set_algorithm( KokkosGraph::COLORING_SERIAL );
353 }
else if(pL.get<std::string>(
"aggregation: coloring algorithm") ==
"vertex based") {
354 coloringHandle->set_algorithm( KokkosGraph::COLORING_VB );
356 }
else if(pL.get<std::string>(
"aggregation: coloring algorithm") ==
"vertex based bit array") {
357 coloringHandle->set_algorithm( KokkosGraph::COLORING_VBBIT );
359 }
else if(pL.get<std::string>(
"aggregation: coloring algorithm") ==
"vertex based color set") {
360 coloringHandle->set_algorithm( KokkosGraph::COLORING_VBCS );
362 }
else if(pL.get<std::string>(
"aggregation: coloring algorithm") ==
"vertex based deterministic") {
363 coloringHandle->set_algorithm( KokkosGraph::COLORING_VBD );
364 if(IsPrint(
Statistics1)) GetOStream(
Statistics1) <<
" algorithm: vertex based deterministic" << std::endl;
365 }
else if(pL.get<std::string>(
"aggregation: coloring algorithm") ==
"vertex based deterministic bit array") {
366 coloringHandle->set_algorithm( KokkosGraph::COLORING_VBDBIT );
367 if(IsPrint(
Statistics1)) GetOStream(
Statistics1) <<
" algorithm: vertex based deterministic bit array" << std::endl;
368 }
else if(pL.get<std::string>(
"aggregation: coloring algorithm") ==
"edge based") {
369 coloringHandle->set_algorithm( KokkosGraph::COLORING_EB );
372 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"Unrecognized distance 1 coloring algorithm");
378 auto graphLW =
dynamic_cast<const LWGraph*
>(&graph);
379 auto graphG =
dynamic_cast<const Graph*
>(&graph);
380 TEUCHOS_TEST_FOR_EXCEPTION(!graphLW && !graphLWK && !graphG,std::invalid_argument,
"Graph is not a LWGraph or LWGraph_kokkos object");
385 KokkosGraph::Experimental::graph_color(&kh,
388 graphLWK->getLocalLWGraph().getRowPtrs(),
389 graphLWK->getLocalLWGraph().getEntries(),
393 auto rowptrs = graphLW->getRowPtrs();
394 auto entries = graphLW->getEntries();
396 Teuchos::Array<size_t> rowptrs_s(rowptrs.size());
397 std::copy(rowptrs.begin(),rowptrs.end(),rowptrs_s.begin());
398 Kokkos::View<const size_t*,Kokkos::LayoutLeft,Kokkos::HostSpace> rowptrs_v(rowptrs_s.data(),(
size_t)rowptrs.size());
399 Kokkos::View<const LO*,Kokkos::LayoutLeft,Kokkos::HostSpace> entries_v(entries.getRawPtr(),(
size_t)entries.size());
400 KokkosGraph::Experimental::graph_color(&kh,
409 RCP<const CrsGraph> graphC = graphG->GetGraph();
410 size_t numEntries = graphC->getLocalNumEntries();
411 ArrayView<const LO> indices;
412 graphC->getLocalRowView(0,indices);
413 Kokkos::View<size_t*,Kokkos::LayoutLeft,Kokkos::HostSpace> rowptrs_v(
"rowptrs_v",graphC->getLocalNumRows()+1);
415 for(LO i=0; i<(LO)graphC->getLocalNumRows()+1; i++)
416 rowptrs_v[i+1] = rowptrs_v[i] + graphC->getNumEntriesInLocalRow(i);
417 Kokkos::View<const LO*,Kokkos::LayoutLeft,Kokkos::HostSpace> entries_v(&indices[0],numEntries);
418 KokkosGraph::Experimental::graph_color(&kh,
428 auto myColors_d = coloringHandle->get_vertex_colors();
429 numColors =
static_cast<LO
>(coloringHandle->get_num_colors());
432 auto myColors_h = Kokkos::create_mirror_view(myColors_d);
433 myColors_out.resize(myColors_h.size());
434 Kokkos::View<LO*,Kokkos::LayoutLeft,Kokkos::HostSpace> myColors_v(&myColors_out[0],myColors_h.size());
435 Kokkos::deep_copy(myColors_v,myColors_h);
438 kh.destroy_graph_coloring_handle();
444 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
450 LO LO_INVALID = Teuchos::OrdinalTraits<LO>::invalid();
451 LO MIS = Teuchos::ScalarTraits<LO>::one();
460 for(LO row=0; row < Nrows; row++) {
461 if(boundaryNodes[row])
464 bool has_colored_neighbor=
false;
465 for(LO j=0; !has_colored_neighbor && j<(LO)indices.size(); j++) {
467 if(myColors[indices[j]] == MIS)
468 has_colored_neighbor=
true;
470 if(!has_colored_neighbor)
477 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
480#ifdef HAVE_MUELU_ZOLTAN2
482 Teuchos::ParameterList params;
483 params.set(
"color_choice",
"FirstFit");
484 params.set(
"color_method",
"D1");
493 GraphAdapter z_adapter(graph);
496 Zoltan2::ColoringProblem<GraphAdapter> problem(&z_adapter,¶ms,graph->GetDomainMap()->getComm());
498 Zoltan2::ColoringSolution<GraphAdapter> * soln = problem.getSolution();
499 ArrayRCP<int> colors = soln->getColorsRCP();
500 numColors = (LO)soln->getNumColors();
504 if(std::is_same<LO,int>::value)
505 myColors_out = colors;
507 myColors_out.resize(colors.size());
508 for(LO i=0; i<(LO)myColors_out.size(); i++)
509 myColors_out[i] = (LO) colors[i];
#define SET_VALID_ENTRY(name)
virtual void DoMISNaive(const GraphBase &graph, Teuchos::ArrayRCP< LO > &myColors, LO &numColors) const
RCP< const ParameterList > GetValidParameterList() const override
Return a const parameter list of valid parameters that setParameterList() will accept.
virtual void GenerateCoarseMap(const Map &fineMap, LO num_c_points, Teuchos::RCP< const Map > &coarseMap) const
void DeclareInput(Level ¤tLevel) const override
Specifies the data that this class needs, and the factories that generate that data.
virtual void DoDistributedGraphColoring(RCP< const GraphBase > &graph, Teuchos::ArrayRCP< LO > &myColors, LO &numColors) const
virtual void DoGraphColoring(const GraphBase &graph, Teuchos::ArrayRCP< LO > &myColors, LO &numColors) const
void Build(Level ¤tLevel) const override
Build an object with this factory.
Timer to be used in factories. Similar to Monitor but with additional timers.
MueLu representation of a graph.
virtual const ArrayRCP< const bool > GetBoundaryNodeMap() const =0
virtual Teuchos::ArrayView< const LocalOrdinal > getNeighborVertices(LocalOrdinal v) const =0
Return the list of vertices adjacent to the vertex 'v'.
virtual size_t GetNodeNumVertices() const =0
Return number of vertices owned by the calling node.
MueLu representation of a compressed row storage graph.
Lightweight MueLu representation of a compressed row storage graph.
Lightweight MueLu representation of a compressed row storage graph.
Class that holds all level-specific information.
int GetLevelID() const
Return level number.
Timer to be used in factories. Similar to SubMonitor but adds a timer level by level.
Namespace for MueLu classes and methods.
@ Debug
Print additional debugging information.
@ Statistics1
Print more statistics.