46 #ifndef MUELU_ZOLTAN2INTERFACE_DEF_HPP 47 #define MUELU_ZOLTAN2INTERFACE_DEF_HPP 53 #if defined(HAVE_MUELU_ZOLTAN2) && defined(HAVE_MPI) 55 #include <Zoltan2_XpetraMultiVectorAdapter.hpp> 56 #include <Zoltan2_PartitioningProblem.hpp> 58 #include <Teuchos_Utils.hpp> 59 #include <Teuchos_DefaultMpiComm.hpp> 60 #include <Teuchos_OpaqueWrapper.hpp> 65 #include "MueLu_Utilities.hpp" 69 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
71 defaultZoltan2Params = rcp(
new ParameterList());
72 defaultZoltan2Params->set(
"algorithm",
"multijagged");
73 defaultZoltan2Params->set(
"partitioning_approach",
"partition");
76 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
78 RCP<ParameterList> validParamList = rcp(
new ParameterList());
80 validParamList->set<
int > (
"rowWeight", 0,
"Default weight to rows (total weight = nnz + rowWeight");
82 validParamList->set< RCP<const FactoryBase> > (
"A", Teuchos::null,
"Factory of the matrix A");
83 validParamList->set< RCP<const FactoryBase> > (
"Coordinates", Teuchos::null,
"Factory of the coordinates");
84 validParamList->set< RCP<const ParameterList> > (
"ParameterList", Teuchos::null,
"Zoltan2 parameters");
86 return validParamList;
90 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
92 Input(currentLevel,
"A");
93 Input(currentLevel,
"Coordinates");
96 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
100 RCP<Matrix> A = Get<RCP<Matrix> >(level,
"A");
101 RCP<const Map> rowMap = A->getRowMap();
103 typedef Xpetra::MultiVector<double, LocalOrdinal, GlobalOrdinal, Node> dMultiVector;
104 RCP<dMultiVector> coords = Get<RCP<dMultiVector> >(level,
"Coordinates");
105 RCP<const Map> map = coords->getMap();
106 GO numElements = map->getNodeNumElements();
108 LO blkSize = A->GetFixedBlockSize();
111 TEUCHOS_TEST_FOR_EXCEPTION(rowMap->getNodeNumElements()/blkSize != coords->getLocalLength(),
Exceptions::Incompatible,
112 "Coordinate vector length (" +
toString(coords->getLocalLength()) <<
" is incompatible with number of block rows in A (" 113 +
toString(rowMap->getNodeNumElements()/blkSize) +
"The vector length should be the same as the number of mesh points.");
114 #ifdef HAVE_MUELU_DEBUG 115 GO indexBase = rowMap->getIndexBase();
116 GetOStream(
Runtime0) <<
"Checking consistence of row and coordinates maps" << std::endl;
118 ArrayView<const GO> rowElements = rowMap->getNodeElementList();
119 ArrayView<const GO> coordsElements = map ->getNodeElementList();
120 for (LO i = 0; i < Teuchos::as<LO>(numElements); i++)
121 TEUCHOS_TEST_FOR_EXCEPTION((coordsElements[i]-indexBase)*blkSize + indexBase != rowElements[i*blkSize],
123 <<
", row GID = " << rowElements[i*blkSize] <<
", blkSize = " << blkSize << std::endl);
126 GO numParts = level.Get<GO>(
"number of partitions");
129 RCP<Xpetra::Vector<GO,LO,GO,NO> > decomposition = Xpetra::VectorFactory<GO, LO, GO, NO>::Build(rowMap,
true);
130 Set(level,
"Partition", decomposition);
134 const ParameterList& pL = GetParameterList();
136 RCP<const ParameterList> providedList = pL.get<RCP<const ParameterList> >(
"ParameterList");
137 ParameterList Zoltan2Params;
138 if (providedList != Teuchos::null)
139 Zoltan2Params = *providedList;
143 for (ParameterList::ConstIterator param = defaultZoltan2Params->begin(); param != defaultZoltan2Params->end(); param++) {
144 const std::string& pName = defaultZoltan2Params->name(param);
145 if (!Zoltan2Params.isParameter(pName))
146 Zoltan2Params.set(pName, defaultZoltan2Params->get<std::string>(pName));
148 Zoltan2Params.set(
"num_global_parts", Teuchos::as<int>(numParts));
150 GetOStream(
Runtime0) <<
"Zoltan2 parameters:\n----------\n" << Zoltan2Params <<
"----------" << std::endl;
152 const std::string& algo = Zoltan2Params.get<std::string>(
"algorithm");
153 TEUCHOS_TEST_FOR_EXCEPTION(algo !=
"multijagged" && algo !=
"rcb", Exceptions::RuntimeError,
154 "Unknown partitioning algorithm: \"" << algo <<
"\"");
156 typedef Zoltan2::XpetraMultiVectorAdapter<dMultiVector> InputAdapterType;
157 typedef Zoltan2::PartitioningProblem<InputAdapterType> ProblemType;
159 int rowWeight = pL.get<
int>(
"rowWeight");
160 GetOStream(
Runtime0) <<
"Using weights formula: nnz + " << rowWeight << std::endl;
162 Array<double> weightsPerRow(numElements);
163 for (LO i = 0; i < numElements; i++) {
164 weightsPerRow[i] = 0.0;
166 for (LO j = 0; j < blkSize; j++) {
167 weightsPerRow[i] += A->getNumEntriesInLocalRow(i*blkSize+j);
173 weightsPerRow[i] += rowWeight;
177 std::vector<int> strides;
178 std::vector<const double*> weights(1, weightsPerRow.getRawPtr());
180 RCP<const Teuchos::MpiComm<int> > dupMpiComm = rcp_dynamic_cast<
const Teuchos::MpiComm<int> >(rowMap->getComm()->duplicate());
181 RCP<const Teuchos::OpaqueWrapper<MPI_Comm> > zoltanComm = dupMpiComm->getRawMpiComm();
183 InputAdapterType adapter(coords, weights, strides);
184 RCP<ProblemType> problem(
new ProblemType(&adapter, &Zoltan2Params, (*zoltanComm)()));
191 RCP<Xpetra::Vector<GO,LO,GO,NO> > decomposition = Xpetra::VectorFactory<GO,LO,GO,NO>::Build(rowMap,
false);
192 ArrayRCP<GO> decompEntries = decomposition->getDataNonConst(0);
194 const typename InputAdapterType::part_t * parts = problem->getSolution().getPartListView();
196 for (GO i = 0; i < numElements; i++) {
197 int partNum = parts[i];
199 for (LO j = 0; j < blkSize; j++)
200 decompEntries[i*blkSize + j] = partNum;
203 Set(level,
"Partition", decomposition);
208 #endif //if defined(HAVE_MUELU_ZOLTAN2) && defined(HAVE_MPI) 210 #endif // MUELU_ZOLTAN2INTERFACE_DEF_HPP std::string toString(const T &what)
Little helper function to convert non-string types to strings.
Timer to be used in factories. Similar to Monitor but with additional timers.
Zoltan2Interface()
Constructor.
One-liner description of what is happening.
Namespace for MueLu classes and methods.
Exception throws to report incompatible objects (like maps).
Class that holds all level-specific information.
Timer to be used in factories. Similar to SubMonitor but adds a timer level by level.
void DeclareInput(Level &level) const
Specifies the data that this class needs, and the factories that generate that data.
void Build(Level &level) const
Build an object with this factory.
Exception throws to report errors in the internal logical of the program.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.