10 #ifndef MUELU_REBALANCETRANSFERFACTORY_DEF_HPP
11 #define MUELU_REBALANCETRANSFERFACTORY_DEF_HPP
14 #include <Teuchos_Tuple.hpp>
17 #include "Xpetra_MultiVectorFactory.hpp"
21 #include <Xpetra_MapFactory.hpp>
29 #include "MueLu_AmalgamationFactory.hpp"
33 #include "MueLu_PerfUtils.hpp"
37 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
40 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
43 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
47 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name))
49 SET_VALID_ENTRY(
"repartition: explicit via new copy rebalance P and R");
54 #undef SET_VALID_ENTRY
58 RCP<validatorType> typeValidator =
rcp(
new validatorType(Teuchos::tuple<std::string>(
"Interpolation",
"Restriction")));
59 validParamList->
set(
"type",
"Interpolation",
"Type of the transfer operator that need to be rebalanced (Interpolation or Restriction)", typeValidator);
62 validParamList->
set<
RCP<const FactoryBase> >(
"P", null,
"Factory of the prolongation operator that need to be rebalanced (only used if type=Interpolation)");
63 validParamList->
set<
RCP<const FactoryBase> >(
"R", null,
"Factory of the restriction operator that need to be rebalanced (only used if type=Restriction)");
64 validParamList->
set<
RCP<const FactoryBase> >(
"Nullspace", null,
"Factory of the nullspace that need to be rebalanced (only used if type=Interpolation)");
65 validParamList->
set<
RCP<const FactoryBase> >(
"Coordinates", null,
"Factory of the coordinates that need to be rebalanced (only used if type=Interpolation)");
66 validParamList->
set<
RCP<const FactoryBase> >(
"Material", null,
"Factory of the material that need to be rebalanced (only used if type=Interpolation)");
67 validParamList->
set<
RCP<const FactoryBase> >(
"BlockNumber", null,
"Factory of the block ids that need to be rebalanced (only used if type=Interpolation)");
68 validParamList->
set<
RCP<const FactoryBase> >(
"Importer", null,
"Factory of the importer object used for the rebalancing");
69 validParamList->
set<
int>(
"write start", -1,
"First level at which coordinates should be written to file");
70 validParamList->
set<
int>(
"write end", -1,
"Last level at which coordinates should be written to file");
76 return validParamList;
79 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
83 if (pL.
get<std::string>(
"type") ==
"Interpolation") {
84 Input(coarseLevel,
"P");
85 if (pL.
get<
bool>(
"repartition: rebalance Nullspace"))
86 Input(coarseLevel,
"Nullspace");
88 Input(coarseLevel,
"Coordinates");
90 Input(coarseLevel,
"Material");
92 Input(coarseLevel,
"BlockNumber");
95 if (pL.
get<
bool>(
"transpose: use implicit") ==
false)
96 Input(coarseLevel,
"R");
99 Input(coarseLevel,
"Importer");
102 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
109 RCP<Matrix> originalP = Get<RCP<Matrix> >(coarseLevel,
"P");
112 if (originalP == Teuchos::null) {
113 Set(coarseLevel,
"P", originalP);
116 int implicit = !pL.
get<
bool>(
"repartition: rebalance P and R");
117 int reallyExplicit = pL.
get<
bool>(
"repartition: explicit via new copy rebalance P and R");
118 int writeStart = pL.
get<
int>(
"write start");
119 int writeEnd = pL.
get<
int>(
"write end");
121 if (writeStart == 0 && fineLevel.
GetLevelID() == 0 && writeStart <= writeEnd && IsAvailable(fineLevel,
"Coordinates")) {
122 std::string fileName =
"coordinates_level_0.m";
124 if (fineCoords != Teuchos::null)
128 if (writeStart == 0 && fineLevel.
GetLevelID() == 0 && writeStart <= writeEnd && IsAvailable(fineLevel,
"BlockNumber")) {
129 std::string fileName =
"BlockNumber_level_0.m";
131 if (fineBlockNumber != Teuchos::null)
143 params->
set(
"printLoadBalancingInfo",
true);
144 params->
set(
"printCommInfo",
true);
147 std::string transferType = pL.
get<std::string>(
"type");
148 if (transferType ==
"Interpolation") {
149 originalP = Get<RCP<Matrix> >(coarseLevel,
"P");
155 if (implicit || importer.
is_null()) {
156 GetOStream(
Runtime0) <<
"Using original prolongator" << std::endl;
157 Set(coarseLevel,
"P", originalP);
184 if (reallyExplicit) {
185 size_t totalMaxPerRow = 0;
186 ArrayRCP<size_t> nnzPerRow(originalP->getRowMap()->getLocalNumElements(), 0);
187 for (
size_t i = 0; i < originalP->getRowMap()->getLocalNumElements(); ++i) {
188 nnzPerRow[i] = originalP->getNumEntriesInLocalRow(i);
189 if (nnzPerRow[i] > totalMaxPerRow) totalMaxPerRow = nnzPerRow[i];
192 rebalancedP = MatrixFactory::Build(originalP->getRowMap(), totalMaxPerRow);
195 RCP<Import> trivialImporter = ImportFactory::Build(originalP->getRowMap(), originalP->getRowMap());
197 std::string monitorLabel{
"Rebalancing prolongator -- import only"};
198 if (
auto sendType = pL.
get<std::string>(
"repartition: send type"); sendType !=
"") {
200 sendTypeList->set(
"Send type", sendType);
201 trivialImporter->setDistributorParameters(sendTypeList);
202 monitorLabel +=
" [" + sendType +
"]";
206 rebalancedP->doImport(*originalP, *trivialImporter,
Xpetra::INSERT);
208 rebalancedP->fillComplete(importer->getTargetMap(), originalP->getRangeMap());
211 rebalancedP = originalP;
216 TEUCHOS_TEST_FOR_EXCEPTION(rebalancedP2 == Teuchos::null, std::runtime_error,
"Xpetra::CrsMatrixWrap doesn't have a CrsMatrix");
219 SubFactoryMonitor subM(*
this,
"Rebalancing prolongator -- fast map replacement", coarseLevel);
224 newImporter = ImportFactory::Build(importer->getTargetMap(), rebalancedP->getColMap());
226 rebalancedP2->replaceDomainMapAndImporter(importer->getTargetMap(), newImporter);
236 std::ostringstream oss;
237 oss <<
"P_" << coarseLevel.GetLevelID();
238 rebalancedP->setObjectLabel(oss.str());
240 Set(coarseLevel,
"P", rebalancedP);
248 if (IsAvailable(coarseLevel,
"Nullspace"))
249 Set(coarseLevel,
"Nullspace", Get<
RCP<MultiVector> >(coarseLevel,
"Nullspace"));
252 if (IsAvailable(coarseLevel,
"Coordinates"))
253 Set(coarseLevel,
"Coordinates", Get<
RCP<xdMV> >(coarseLevel,
"Coordinates"));
256 if (IsAvailable(coarseLevel,
"Material"))
257 Set(coarseLevel,
"Material", Get<
RCP<MultiVector> >(coarseLevel,
"Material"));
260 if (IsAvailable(coarseLevel,
"BlockNumber"))
268 IsAvailable(coarseLevel,
"Coordinates")) {
269 RCP<xdMV> coords = Get<RCP<xdMV> >(coarseLevel,
"Coordinates");
275 LO nodeNumElts = coords->getMap()->getLocalNumElements();
276 LO myBlkSize = 0, blkSize = 0;
278 myBlkSize = importer->getSourceMap()->getLocalNumElements() / nodeNumElts;
279 MueLu_maxAll(coords->getMap()->getComm(), myBlkSize, blkSize);
284 coordImporter = importer;
287 std::vector<size_t> stridingInfo{Teuchos::as<size_t>(blkSize)};
289 importer->getTargetMap()->getLocalElementList(), origMap->getIndexBase(), stridingInfo, origMap->getComm());
293 coordImporter = ImportFactory::Build(origMap, targetVectorMap);
297 permutedCoords->doImport(*coords, *coordImporter,
Xpetra::INSERT);
299 if (pL.
isParameter(
"repartition: use subcommunicators") ==
true && pL.
get<
bool>(
"repartition: use subcommunicators") ==
true)
300 permutedCoords->replaceMap(permutedCoords->getMap()->removeEmptyProcesses());
302 if (permutedCoords->getMap() == Teuchos::null)
303 permutedCoords = Teuchos::null;
305 Set(coarseLevel,
"Coordinates", permutedCoords);
307 std::string fileName =
"rebalanced_coordinates_level_" +
toString(coarseLevel.GetLevelID()) +
".m";
308 if (writeStart <= coarseLevel.GetLevelID() && coarseLevel.GetLevelID() <= writeEnd && permutedCoords->getMap() != Teuchos::null)
312 if (IsAvailable(coarseLevel,
"Material")) {
313 RCP<MultiVector> material = Get<RCP<MultiVector> >(coarseLevel,
"Material");
319 LO nodeNumElts = material->getMap()->getLocalNumElements();
320 LO myBlkSize = 0, blkSize = 0;
322 myBlkSize = importer->getSourceMap()->getLocalNumElements() / nodeNumElts;
323 MueLu_maxAll(material->getMap()->getComm(), myBlkSize, blkSize);
328 materialImporter = importer;
331 std::vector<size_t> stridingInfo{Teuchos::as<size_t>(blkSize)};
333 importer->getTargetMap()->getLocalElementList(), origMap->getIndexBase(), stridingInfo, origMap->getComm());
337 materialImporter = ImportFactory::Build(origMap, targetVectorMap);
340 RCP<MultiVector> permutedMaterial = MultiVectorFactory::Build(materialImporter->getTargetMap(), material->getNumVectors());
341 permutedMaterial->doImport(*material, *materialImporter,
Xpetra::INSERT);
343 if (pL.
get<
bool>(
"repartition: use subcommunicators") ==
true)
344 permutedMaterial->replaceMap(permutedMaterial->getMap()->removeEmptyProcesses());
346 if (permutedMaterial->getMap() == Teuchos::null)
347 permutedMaterial = Teuchos::null;
349 Set(coarseLevel,
"Material", permutedMaterial);
354 IsAvailable(coarseLevel,
"BlockNumber")) {
360 RCP<LocalOrdinalVector> permutedBlockNumber = LocalOrdinalVectorFactory::Build(importer->getTargetMap(),
false);
361 permutedBlockNumber->doImport(*BlockNumber, *importer,
Xpetra::INSERT);
363 if (pL.
isParameter(
"repartition: use subcommunicators") ==
true && pL.
get<
bool>(
"repartition: use subcommunicators") ==
true)
364 permutedBlockNumber->replaceMap(permutedBlockNumber->getMap()->removeEmptyProcesses());
366 if (permutedBlockNumber->getMap() == Teuchos::null)
367 permutedBlockNumber = Teuchos::null;
369 Set(coarseLevel,
"BlockNumber", permutedBlockNumber);
371 std::string fileName =
"rebalanced_BlockNumber_level_" +
toString(coarseLevel.GetLevelID()) +
".m";
372 if (writeStart <= coarseLevel.GetLevelID() && coarseLevel.GetLevelID() <= writeEnd && permutedBlockNumber->getMap() != Teuchos::null)
376 if (IsAvailable(coarseLevel,
"Nullspace")) {
377 RCP<MultiVector> nullspace = Get<RCP<MultiVector> >(coarseLevel,
"Nullspace");
382 RCP<MultiVector> permutedNullspace = MultiVectorFactory::Build(importer->getTargetMap(), nullspace->getNumVectors());
383 permutedNullspace->doImport(*nullspace, *importer,
Xpetra::INSERT);
385 if (pL.
get<
bool>(
"repartition: use subcommunicators") ==
true)
386 permutedNullspace->replaceMap(permutedNullspace->getMap()->removeEmptyProcesses());
388 if (permutedNullspace->getMap() == Teuchos::null)
389 permutedNullspace = Teuchos::null;
391 Set(coarseLevel,
"Nullspace", permutedNullspace);
395 if (pL.
get<
bool>(
"transpose: use implicit") ==
false) {
396 RCP<Matrix> originalR = Get<RCP<Matrix> >(coarseLevel,
"R");
400 if (implicit || importer.
is_null()) {
401 GetOStream(
Runtime0) <<
"Using original restrictor" << std::endl;
402 Set(coarseLevel,
"R", originalR);
407 SubFactoryMonitor subM(*
this,
"Rebalancing restriction -- fusedImport", coarseLevel);
411 listLabel.
set(
"Timer Label",
"MueLu::RebalanceR-" +
Teuchos::toString(coarseLevel.GetLevelID()));
412 rebalancedR = MatrixFactory::Build(originalR, *importer, dummy, importer->getTargetMap(),
Teuchos::rcp(&listLabel,
false));
415 std::ostringstream oss;
416 oss <<
"R_" << coarseLevel.GetLevelID();
417 rebalancedR->setObjectLabel(oss.str());
419 Set(coarseLevel,
"R", rebalancedR);
437 #endif // MUELU_REBALANCETRANSFERFACTORY_DEF_HPP
Exception indicating invalid cast attempted.
T & Get(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Get data without decrementing associated storage counter (i.e., read-only access). Usage: Level->Get< RCP<Matrix> >("A", factory) if factory == NULL => use default factory.
std::string toString(const T &what)
Little helper function to convert non-string types to strings.
#define MueLu_maxAll(rcpComm, in, out)
T & get(const std::string &name, T def_value)
Timer to be used in factories. Similar to Monitor but with additional timers.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
One-liner description of what is happening.
ParameterList & set(std::string const &name, T &&value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
static const NoFactory * get()
virtual ~RebalanceTransferFactory()
Destructor.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
Print even more statistics.
bool isParameter(const std::string &name) const
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Specifies the data that this class needs, and the factories that generate that data.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Class that holds all level-specific information.
Timer to be used in factories. Similar to SubMonitor but adds a timer level by level.
static void AmalgamateMap(const Map &sourceMap, const Matrix &A, RCP< const Map > &amalgamatedMap, Array< LO > &translation)
Method to create merged map for systems of PDEs.
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
#define SET_VALID_ENTRY(name)
static std::string PrintMatrixInfo(const Matrix &A, const std::string &msgTag, RCP< const Teuchos::ParameterList > params=Teuchos::null)
RebalanceTransferFactory()
Constructor.
int GetLevelID() const
Return level number.
std::string toString(const T &t)