47 #ifndef MUELU_MAPTRANSFERFACTORY_DEF_HPP_ 48 #define MUELU_MAPTRANSFERFACTORY_DEF_HPP_ 62 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
68 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
70 RCP<ParameterList> validParamList = rcp(
new ParameterList());
78 validParamList->setEntry(
"map: name", Teuchos::ParameterEntry(std::string(
"")));
79 validParamList->setEntry(
"map: factory", Teuchos::ParameterEntry(std::string(
"null")));
80 validParamList->set< RCP<const FactoryBase> >(
"P", Teuchos::null,
"Tentative prolongator factory");
82 return validParamList;
85 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
87 const ParameterList & pL = GetParameterList();
88 std::string mapFactName = pL.get<std::string>(
"map: factory");
89 std::string mapName = pL.get<std::string>(
"map: name");
92 if (mapFactName ==
"" || mapFactName ==
"NoFactory") {
95 else if (mapFactName !=
"null") {
104 Teuchos::RCP<const FactoryBase> tentPFact = GetFactory(
"P");
105 if (tentPFact == Teuchos::null)
110 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
112 Monitor m(*
this,
"Contact Map transfer factory");
114 const ParameterList & pL = GetParameterList();
115 std::string mapName = pL.get<std::string>(
"map: name");
117 if (fineLevel.
IsAvailable(mapName, mapFact_.get())==
false)
118 GetOStream(
Runtime0) <<
"MapTransferFactory::Build: User provided map " << mapName <<
" not found in Level class." << std::endl;
121 RCP<Map> transferMap = fineLevel.
Get<RCP<Map> >(mapName,mapFact_.get());
126 RCP<const FactoryBase> tentPFact = GetFactory(
"P");
127 if (tentPFact == Teuchos::null)
130 "MueLu::MapTransferFactory::Build(): P (generated by TentativePFactory) not available.");
131 RCP<Matrix> Ptent = coarseLevel.
Get<RCP<Matrix> >(
"P", tentPFact.get());
133 std::vector<GO > coarseMapGids;
136 for (
size_t row=0; row<Ptent->getNodeNumRows(); row++) {
137 GO grid = Ptent->getRowMap()->getGlobalElement(row);
139 if (transferMap->isNodeGlobalElement(grid)) {
140 Teuchos::ArrayView<const LO> indices;
141 Teuchos::ArrayView<const SC> vals;
142 Ptent->getLocalRowView(row, indices, vals);
144 for (
size_t i = 0; i < as<size_t>(indices.size()); i++) {
146 GO gcid = Ptent->getColMap()->getGlobalElement(indices[i]);
147 coarseMapGids.push_back(gcid);
153 const GO INVALID = Teuchos::OrdinalTraits<Xpetra::global_size_t>::invalid();
154 std::sort(coarseMapGids.begin(), coarseMapGids.end());
155 coarseMapGids.erase(std::unique(coarseMapGids.begin(), coarseMapGids.end()), coarseMapGids.end());
156 Teuchos::ArrayView<GO> coarseMapGidsView(&coarseMapGids[0], coarseMapGids.size());
157 Teuchos::RCP<Map> coarseTransferMap = MapFactory::Build(Ptent->getColMap()->lib(), INVALID, coarseMapGidsView, Ptent->getColMap()->getIndexBase(), Ptent->getColMap()->getComm());
160 coarseLevel.
Set(mapName, coarseTransferMap, mapFact_.get());
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.
One-liner description of what is happening.
Namespace for MueLu classes and methods.
Class that holds all level-specific information.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
void Set(const std::string &ename, const T &entry, const FactoryBase *factory=NoFactory::get())
Timer to be used in non-factories.
const RCP< const FactoryManagerBase > GetFactoryManager()
returns the current factory manager
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.
Exception throws to report errors in the internal logical of the program.
MapTransferFactory()
Constructor.
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
static const RCP< const NoFactory > getRCP()
Static Get() functions.