48 #ifndef MUELU_REORDERBLOCKAFACTORY_DEF_HPP_
49 #define MUELU_REORDERBLOCKAFACTORY_DEF_HPP_
56 #include <Xpetra_MapExtractor.hpp>
66 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
72 validParamList->
set<std::string>(
"Reorder Type",
"",
"String describing the reordering of blocks");
75 validParamList->
set<
RCP<const FactoryBase> >(
"Map1", Teuchos::null,
"Generating factory of the fine level map associated with the (1,1) block in your n x n block matrix.");
76 validParamList->
set<
RCP<const FactoryBase> >(
"Map2", Teuchos::null,
"Generating factory of the fine level map associated with the (2,2) block in your n x n block matrix.");
77 validParamList->
set<
RCP<const FactoryBase> >(
"Map3", Teuchos::null,
"Generating factory of the fine level map associated with the (3,3) block in your n x n block matrix.");
78 validParamList->
set<
RCP<const FactoryBase> >(
"Map4", Teuchos::null,
"Generating factory of the fine level map associated with the (4,4) block in your n x n block matrix.");
79 validParamList->
set<
RCP<const FactoryBase> >(
"Map5", Teuchos::null,
"Generating factory of the fine level map associated with the (5,5) block in your n x n block matrix.");
80 validParamList->
set<
RCP<const FactoryBase> >(
"Map6", Teuchos::null,
"Generating factory of the fine level map associated with the (6,6) block in your n x n block matrix.");
81 validParamList->
set<
RCP<const FactoryBase> >(
"Map7", Teuchos::null,
"Generating factory of the fine level map associated with the (7,7) block in your n x n block matrix.");
82 validParamList->
set<
RCP<const FactoryBase> >(
"Map8", Teuchos::null,
"Generating factory of the fine level map associated with the (8,8) block in your n x n block matrix.");
83 validParamList->
set<
RCP<const FactoryBase> >(
"Map9", Teuchos::null,
"Generating factory of the fine level map associated with the (9,9) block in your n x n block matrix.");
85 return validParamList;
88 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
90 Input(currentLevel,
"A");
93 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
98 std::string reorderStr = pL.
get<std::string>(
"Reorder Type");
100 RCP<Matrix> Ain = Get<RCP<Matrix> >(currentLevel,
"A");
106 if (A == Teuchos::null && currentLevel.GetLevelID() == 0) {
107 GetOStream(
Warnings0) <<
"Split input matrix (Warning: this is a rather expensive operation)" << std::endl;
109 std::vector<Teuchos::RCP<const Map> > xmaps;
111 for (
int it = 1; it < 10; it++) {
112 std::stringstream ss;
116 GetOStream(
Runtime1) <<
"Use user-given submap #" << it <<
": length dimension=" << submap->getGlobalNumElements() << std::endl;
117 xmaps.push_back(submap);
121 bool bThyraMode =
false;
142 GetOStream(
Statistics1) <<
"Got a " << A->Rows() <<
"x" << A->Cols() <<
" blocked operator as input" << std::endl;
145 if (reorderStr.empty()) {
146 GetOStream(
Statistics1) <<
"No reordering information provided. Skipping reordering of A." << std::endl;
149 GetOStream(
Debug) <<
"Reordering A using " << brm->
toString() << std::endl;
152 Teuchos::rcp_dynamic_cast<
const ReorderedBlockedCrsMatrix>(
153 Xpetra::buildReorderedBlockedCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(brm, A));
156 "Block reordering of " << A->Rows() <<
"x" << A->Cols()
157 <<
" blocked operator failed.");
159 GetOStream(
Statistics1) <<
"Reordering A using " << brm->toString() <<
" block gives a " << brop->Rows() <<
"x"
160 << brop->Cols() <<
" blocked operator" << std::endl;
161 GetOStream(
Debug) <<
"Reordered operator has " << brop->getRangeMap()->getGlobalNumElements() <<
" rows and "
162 << brop->getDomainMap()->getGlobalNumElements() <<
" columns" << std::endl;
163 GetOStream(
Debug) <<
"Reordered operator: Use of Thyra style gids = "
164 << brop->getRangeMapExtractor()->getThyraMode() << std::endl;
168 Teuchos::rcp_const_cast<ReorderedBlockedCrsMatrix>(brop);
173 currentLevel.Set(
"A", Teuchos::rcp_dynamic_cast<Matrix>(A),
this);
Important warning messages (one line)
Exception indicating invalid cast attempted.
virtual std::string toString() const
Teuchos::RCP< const Xpetra::BlockReorderManager > blockedReorderFromString(std::string reorder)
T & get(const std::string &name, T def_value)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
Timer to be used in factories. Similar to Monitor but with additional timers.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Print additional debugging information.
static Teuchos::RCP< Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > SplitMatrix(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &input, Teuchos::RCP< const Xpetra::MapExtractor< Scalar, LocalOrdinal, GlobalOrdinal, Node > > rangeMapExtractor, Teuchos::RCP< const Xpetra::MapExtractor< Scalar, LocalOrdinal, GlobalOrdinal, Node > > domainMapExtractor, Teuchos::RCP< const Xpetra::MapExtractor< Scalar, LocalOrdinal, GlobalOrdinal, Node > > columnMapExtractor=Teuchos::null, bool bThyraMode=false)
static const NoFactory * get()
RCP< const ParameterList > GetValidParameterList() const
Input.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Class that holds all level-specific information.
void Build(Level ¤tLevel) const
Build an object with this factory.
void DeclareInput(Level ¤tLevel) const
Specifies the data that this class needs, and the factories that generate that data.
Exception throws to report errors in the internal logical of the program.
Description of what is happening (more verbose)
static const RCP< const NoFactory > getRCP()
Static Get() functions.