47 #include "Teko_BlockInvDiagonalStrategy.hpp"
48 #include "Teko_InverseFactory.hpp"
52 InvFactoryDiagStrategy::InvFactoryDiagStrategy(
const Teuchos::RCP<InverseFactory>& factory) {
54 invDiagFact_.resize(1, factory);
55 defaultInvFact_ = factory;
58 InvFactoryDiagStrategy::InvFactoryDiagStrategy(
59 const std::vector<Teuchos::RCP<InverseFactory> >& factories,
60 const Teuchos::RCP<InverseFactory>& defaultFact) {
61 invDiagFact_ = factories;
63 if (defaultFact == Teuchos::null)
64 defaultInvFact_ = invDiagFact_[0];
66 defaultInvFact_ = defaultFact;
69 InvFactoryDiagStrategy::InvFactoryDiagStrategy(
70 const std::vector<Teuchos::RCP<InverseFactory> >& inverseFactories,
71 const std::vector<Teuchos::RCP<InverseFactory> >& preconditionerFactories,
72 const Teuchos::RCP<InverseFactory>& defaultInverseFact,
73 const Teuchos::RCP<InverseFactory>& defaultPreconditionerFact) {
74 invDiagFact_ = inverseFactories;
75 precDiagFact_ = preconditionerFactories;
77 if (defaultInverseFact == Teuchos::null)
78 defaultInvFact_ = invDiagFact_[0];
80 defaultInvFact_ = defaultInverseFact;
81 defaultPrecFact_ = defaultPreconditionerFact;
88 std::vector<LinearOp>& invDiag)
const {
89 Teko_DEBUG_SCOPE(
"InvFactoryDiagStrategy::getInvD", 10);
92 size_t diagCnt = A->productRange()->numBlocks();
94 Teko_DEBUG_MSG(
"# diags = " << diagCnt <<
", # inverses = " << invDiagFact_.size(), 6);
96 const std::string opPrefix =
"JacobiDiagOp";
97 for (
size_t i = 0; i < diagCnt; i++) {
98 auto precFact = ((i < precDiagFact_.size()) && (!precDiagFact_[i].is_null()))
101 auto invFact = (i < invDiagFact_.size()) ? invDiagFact_[i] : defaultInvFact_;
102 invDiag.push_back(
buildInverse(*invFact, precFact, getBlock(i, i, A), state, opPrefix, i));
107 Teuchos::RCP<InverseFactory>& precFact,
108 const LinearOp& matrix,
110 const std::string& opPrefix,
int i)
const {
111 std::stringstream ss;
112 ss << opPrefix <<
"_" << i;
115 ModifiableLinearOp& precOp = state.
getModifiableOp(
"prec_" + ss.str());
117 if (precFact != Teuchos::null) {
118 if (precOp == Teuchos::null) {
119 precOp = precFact->buildInverse(matrix);
126 if (invOp == Teuchos::null)
127 if (precOp.is_null())
132 if (precOp.is_null())
void rebuildInverse(const InverseFactory &factory, const LinearOp &A, InverseLinearOp &invA)
Abstract class for building an inverse operator.
virtual void getInvD(const BlockedLinearOp &A, BlockPreconditionerState &state, std::vector< LinearOp > &invDiag) const
An implementation of a state object for block preconditioners.
LinearOp buildInverse(const InverseFactory &invFact, Teuchos::RCP< InverseFactory > &precFact, const LinearOp &matrix, BlockPreconditionerState &state, const std::string &opPrefix, int i) const
Conveinence function for building inverse operators.
virtual void addModifiableOp(const std::string &name, const Teko::ModifiableLinearOp &mlo)
Add a named operator to the state object.
InverseLinearOp buildInverse(const InverseFactory &factory, const LinearOp &A)
Build an inverse operator using a factory and a linear operator.
virtual Teko::ModifiableLinearOp & getModifiableOp(const std::string &name)
Add a named operator to the state object.