46 #ifndef PACKAGES_MUELU_SRC_REBALANCING_MUELU_CLONEREPARTITIONINTERFACE_DEF_HPP_ 47 #define PACKAGES_MUELU_SRC_REBALANCING_MUELU_CLONEREPARTITIONINTERFACE_DEF_HPP_ 49 #include <Xpetra_MultiVectorFactory.hpp> 50 #include <Xpetra_VectorFactory.hpp> 57 #include "MueLu_Utilities.hpp" 62 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
64 Teuchos::RCP<ParameterList> validParamList = rcp(
new ParameterList());
65 validParamList->set< Teuchos::RCP<const FactoryBase> >(
"A", Teuchos::null,
"Factory of the matrix A");
66 validParamList->set< Teuchos::RCP<const FactoryBase> >(
"number of partitions", Teuchos::null,
"Instance of RepartitionHeuristicFactory.");
67 validParamList->set< Teuchos::RCP<const FactoryBase> >(
"Partition", Teuchos::null,
"Factory generating the Partition array (e.g. ZoltanInterface)");
69 return validParamList;
73 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
75 Input(currentLevel,
"A");
76 Input(currentLevel,
"Partition");
79 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
84 Teuchos::RCP<Matrix> A = Get< Teuchos::RCP<Matrix> > (currentLevel,
"A");
85 Teuchos::RCP<const Teuchos::Comm< int > > comm = A->getRowMap()->getComm();
88 if (currentLevel.IsAvailable(
"number of partitions")) {
89 GetOStream(
Warnings0) <<
"Using user-provided \"number of partitions\", the performance is unknown" << std::endl;
96 RCP<GOVector> decomposition = Teuchos::null;
99 decomposition = Get<RCP<GOVector> >(currentLevel,
"Partition");
100 ArrayRCP<const GO> decompEntries = decomposition->getData(0);
102 if (decomposition.is_null()) {
103 GetOStream(
Warnings0) <<
"No repartitioning necessary: partitions were left unchanged by the repartitioner" << std::endl;
104 Set<RCP<const Import> >(currentLevel,
"Importer", Teuchos::null);
109 Teuchos::RCP<GOVector> ret = Xpetra::VectorFactory<GO, LO, GO, NO>::Build(A->getRowMap(),
false);
110 ArrayRCP<GO> retDecompEntries = ret->getDataNonConst(0);
113 LocalOrdinal blkSize = 1;
116 if(A->IsView(
"stridedMaps") &&
117 Teuchos::rcp_dynamic_cast<
const StridedMap>(A->getRowMap(
"stridedMaps")) != Teuchos::null) {
118 Xpetra::viewLabel_t oldView = A->SwitchToView(
"stridedMaps");
119 RCP<const StridedMap> strMap = Teuchos::rcp_dynamic_cast<
const StridedMap>(A->getRowMap());
120 TEUCHOS_TEST_FOR_EXCEPTION(strMap == Teuchos::null,
Exceptions::BadCast,
"MueLu::CloneRepartitionInterface::Build: cast to strided row map failed.");
121 LocalOrdinal stridedBlock = strMap->getStridedBlockId();
122 if (stridedBlock == -1)
123 blkSize = strMap->getFixedBlockSize();
125 std::vector<size_t> strInfo = strMap->getStridingData();
126 blkSize = strInfo[stridedBlock];
128 oldView = A->SwitchToView(oldView);
129 GetOStream(
Statistics1) <<
"CloneRepartitionInterface::Build():" <<
" found blockdim=" << blkSize <<
" from strided maps."<< std::endl;
131 GetOStream(
Statistics1) <<
"CloneRepartitionInterface::Build(): no striding information available. Use blockdim=" << blkSize <<
" (DofsPerNode)." << std::endl;
132 blkSize = A->GetFixedBlockSize();
136 size_t inLocalLength = decomposition->getLocalLength();
137 size_t outLocalLength = A->getRowMap()->getNodeNumElements();
140 size_t numLocalNodes = outLocalLength / blkSize;
141 TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::as<size_t>(outLocalLength % blkSize) != 0,
MueLu::Exceptions::RuntimeError,
"CloneRepartitionInterface: inconsistent number of local DOFs (" << outLocalLength <<
") and degrees of freedoms (" << blkSize <<
")");
143 if (numLocalNodes > 0) {
144 TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::as<size_t>(inLocalLength % numLocalNodes) != 0,
MueLu::Exceptions::RuntimeError,
"CloneRepartitionInterface: inconsistent number of local DOFs (" << inLocalLength <<
") and number of local nodes (" << numLocalNodes <<
")");
145 LocalOrdinal inBlkSize = Teuchos::as<LocalOrdinal>(inLocalLength / numLocalNodes);
148 for(LO i = 0; i<Teuchos::as<LO>(numLocalNodes); i++) {
149 for(LO j = 0; j < blkSize; j++) {
150 retDecompEntries[i*blkSize + j] = Teuchos::as<GO>(decompEntries[i*inBlkSize]);
154 Set(currentLevel,
"Partition", ret);
Important warning messages (one line)
Exception indicating invalid cast attempted.
Timer to be used in factories. Similar to Monitor but with additional timers.
Namespace for MueLu classes and methods.
Print statistics that do not involve significant additional computation.
Class that holds all level-specific information.
void print(std::ostream &out, const VerbLevel verbLevel=Default) const
Printing method.
Exception throws to report errors in the internal logical of the program.
void Build(Level &level) const
Build an object with this factory.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
void DeclareInput(Level &level) const
Specifies the data that this class needs, and the factories that generate that data.