46 #ifndef MUELU_REBALANCEMAPFACTORY_DEF_HPP_ 47 #define MUELU_REBALANCEMAPFACTORY_DEF_HPP_ 51 #include <Teuchos_Utils.hpp> 58 #include "MueLu_Utilities.hpp" 62 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
64 RCP<ParameterList> validParamList = rcp(
new ParameterList());
66 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name)) 68 #undef SET_VALID_ENTRY 71 validParamList->set< std::string > (
"Map name" ,
"",
"Name of map to rebalanced.");
72 validParamList->set< RCP<const FactoryBase> >(
"Map factory",
MueLu::NoFactory::getRCP(),
"Generating factory of map to be rebalanced.");
75 validParamList->set< RCP<const FactoryBase> >(
"Importer", Teuchos::null,
"Factory of the importer object used for the rebalancing");
77 return validParamList;
81 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
83 const Teuchos::ParameterList & pL = GetParameterList();
84 std::string mapName = pL.get<std::string> (
"Map name");
85 Teuchos::RCP<const FactoryBase> mapFactory = GetFactory (
"Map factory");
86 currentLevel.
DeclareInput(mapName,mapFactory.get(),
this);
88 Input(currentLevel,
"Importer");
91 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
98 const Teuchos::ParameterList & pL = GetParameterList();
99 std::string mapName = pL.get<std::string> (
"Map name");
100 Teuchos::RCP<const FactoryBase> mapFactory = GetFactory (
"Map factory");
102 RCP<const Import> rebalanceImporter = Get<RCP<const Import> >(level,
"Importer");
104 if(rebalanceImporter != Teuchos::null) {
106 RCP<const Map> map = level.Get< RCP<const Map> >(mapName,mapFactory.get());
110 RCP<Vector> v = VectorFactory::Build(map);
115 RCP<const Import> blowUpImporter = ImportFactory::Build(map, rebalanceImporter->getSourceMap());
116 RCP<Vector> pv = VectorFactory::Build(rebalanceImporter->getSourceMap());
117 pv->doImport(*v,*blowUpImporter,Xpetra::INSERT);
120 RCP<Vector> ptv = VectorFactory::Build(rebalanceImporter->getTargetMap());
121 ptv->doImport(*pv,*rebalanceImporter,Xpetra::INSERT);
123 if (pL.get<
bool>(
"repartition: use subcommunicators") ==
true)
124 ptv->replaceMap(ptv->getMap()->removeEmptyProcesses());
127 Teuchos::ArrayRCP< const Scalar > ptvData = ptv->getData(0);
128 std::vector<GlobalOrdinal> localGIDs;
130 for (
size_t k = 0; k < ptv->getLocalLength(); k++) {
131 if(ptvData[k] == 1.0) {
132 localGIDs.push_back(ptv->getMap()->getGlobalElement(k));
136 const Teuchos::ArrayView<const GlobalOrdinal> localGIDs_view(&localGIDs[0],localGIDs.size());
138 Teuchos::RCP<const Map> localGIDsMap = MapFactory::Build(
140 Teuchos::OrdinalTraits<int>::invalid(),
142 0, ptv->getMap()->getComm());
146 level.Set(mapName, localGIDsMap, mapFactory.get());
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
#define SET_VALID_ENTRY(name)
Timer to be used in factories. Similar to Monitor but with additional timers.
Namespace for MueLu classes and methods.
Class that holds all level-specific information.
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.
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()
static const RCP< const NoFactory > getRCP()
Static Get() functions.