42 #ifndef THYRA_DEFAULT_INVERSE_LINEAR_OP_DEF_HPP 43 #define THYRA_DEFAULT_INVERSE_LINEAR_OP_DEF_HPP 45 #include "Thyra_DefaultInverseLinearOp_decl.hpp" 46 #include "Thyra_MultiVectorStdOps.hpp" 47 #include "Thyra_AssertOp.hpp" 48 #include "Teuchos_Utils.hpp" 49 #include "Teuchos_TypeNameTraits.hpp" 58 template<
class Scalar>
63 template<
class Scalar>
73 lows,fwdSolveCriteria,throwOnFwdSolveFailure
74 ,adjSolveCriteria,throwOnAdjSolveFailure
79 template<
class Scalar>
81 const Teuchos::RCP<
const LinearOpWithSolveBase<Scalar> > &lows,
82 const SolveCriteria<Scalar> *fwdSolveCriteria,
84 const SolveCriteria<Scalar> *adjSolveCriteria,
89 lows,fwdSolveCriteria,throwOnFwdSolveFailure
90 ,adjSolveCriteria,throwOnAdjSolveFailure
95 template<
class Scalar>
97 const Teuchos::RCP<LinearOpWithSolveBase<Scalar> > &lows,
98 const SolveCriteria<Scalar> *fwdSolveCriteria,
100 const SolveCriteria<Scalar> *adjSolveCriteria,
105 lows,fwdSolveCriteria,throwOnFwdSolveFailure
106 ,adjSolveCriteria,throwOnAdjSolveFailure
111 template<
class Scalar>
113 const Teuchos::RCP<
const LinearOpWithSolveBase<Scalar> > &lows,
114 const SolveCriteria<Scalar> *fwdSolveCriteria,
116 const SolveCriteria<Scalar> *adjSolveCriteria,
121 lows,fwdSolveCriteria,throwOnFwdSolveFailure
122 ,adjSolveCriteria,throwOnAdjSolveFailure
127 template<
class Scalar>
130 lows_.uninitialize();
131 fwdSolveCriteria_ = Teuchos::null;
132 adjSolveCriteria_ = Teuchos::null;
139 template<
class Scalar>
142 return lows_.isConst();
146 template<
class Scalar>
147 Teuchos::RCP<LinearOpWithSolveBase<Scalar> >
150 return lows_.getNonconstObj();
154 template<
class Scalar>
155 Teuchos::RCP<const LinearOpWithSolveBase<Scalar> >
158 return lows_.getConstObj();
165 template<
class Scalar>
166 Teuchos::RCP< const VectorSpaceBase<Scalar> >
170 return lows_.getConstObj()->domain();
174 template<
class Scalar>
175 Teuchos::RCP< const VectorSpaceBase<Scalar> >
179 return lows_.getConstObj()->range();
183 template<
class Scalar>
184 Teuchos::RCP<const LinearOpBase<Scalar> >
187 return Teuchos::null;
194 template<
class Scalar>
198 std::ostringstream oss;
200 << Teuchos::Describable::description() <<
"{" 201 <<
"lows="<<lows_.getConstObj()->description()
202 <<
",fwdSolveCriteria="<<(fwdSolveCriteria_.get()?
"...":
"DEFAULT")
203 <<
",adjSolveCriteria="<<(adjSolveCriteria_.get()?
"...":
"DEFAULT")
209 template<
class Scalar>
211 Teuchos::FancyOStream &out,
212 const Teuchos::EVerbosityLevel verbLevel
216 using Teuchos::OSTab;
220 case Teuchos::VERB_DEFAULT:
221 case Teuchos::VERB_LOW:
222 out << this->description() << std::endl;
224 case Teuchos::VERB_MEDIUM:
225 case Teuchos::VERB_HIGH:
226 case Teuchos::VERB_EXTREME:
229 << Teuchos::Describable::description() <<
"{" 230 <<
"rangeDim=" << this->range()->dim()
231 <<
",domainDim=" << this->domain()->dim() <<
"}:\n";
234 if(!lows_.getConstObj().get()) {
238 out << Teuchos::describe(*lows_.getConstObj(),verbLevel);
243 TEUCHOS_TEST_FOR_EXCEPT(
true);
254 template<
class Scalar>
257 if (nonnull(lows_)) {
258 return solveSupports(*lows_.getConstObj(),M_trans);
264 template<
class Scalar>
273 typedef Teuchos::ScalarTraits<Scalar> ST;
284 Teuchos::RCP<MultiVectorBase<Scalar> > T;
285 if(beta==ST::zero()) {
286 T = Teuchos::rcpFromPtr(Y);
289 T = createMembers(Y->range(),Y->domain()->dim());
293 const Ptr<const SolveCriteria<Scalar> > solveCriteria =
296 ? fwdSolveCriteria_.ptr()
297 : adjSolveCriteria_.ptr()
299 assign(T.ptr(), ST::zero());
301 Thyra::solve<Scalar>(*lows_.getConstObj(), M_trans, X, T.ptr(), solveCriteria);
303 TEUCHOS_TEST_FOR_EXCEPTION(
309 ,
"Error, the LOWS object " << lows_.getConstObj()->description() <<
" returned an unconverged" 310 "status of " <<
toString(solveStatus.solveStatus) <<
" with the message " 311 << solveStatus.message <<
"." 314 if(beta==ST::zero()) {
318 update( alpha, *T, Y );
326 template<
class Scalar>
329 const Teuchos::RCP<LOWS> &lows,
340 fwdSolveCriteria_ = Teuchos::null;
344 adjSolveCriteria_ = Teuchos::null;
345 throwOnFwdSolveFailure_ = throwOnFwdSolveFailure;
346 throwOnAdjSolveFailure_ = throwOnAdjSolveFailure;
347 const std::string lowsLabel = lows_.getConstObj()->getObjectLabel();
348 if(lowsLabel.length())
349 this->setObjectLabel(
"inv("+lowsLabel+
")" );
359 template<
class Scalar>
360 Teuchos::RCP<Thyra::LinearOpBase<Scalar> >
361 Thyra::nonconstInverse(
362 const RCP<LinearOpWithSolveBase<Scalar> > &A,
363 const Ptr<
const SolveCriteria<Scalar> > &fwdSolveCriteria,
365 const Ptr<
const SolveCriteria<Scalar> > &adjSolveCriteria,
370 new DefaultInverseLinearOp<Scalar>(
371 A, fwdSolveCriteria.get(), throwOnFwdSolveFailure,
372 adjSolveCriteria.get(), throwOnAdjSolveFailure
377 template<
class Scalar>
378 Teuchos::RCP<Thyra::LinearOpBase<Scalar> >
380 const RCP<
const LinearOpWithSolveBase<Scalar> > &A,
381 const Ptr<
const SolveCriteria<Scalar> > &fwdSolveCriteria,
383 const Ptr<
const SolveCriteria<Scalar> > &adjSolveCriteria,
388 new DefaultInverseLinearOp<Scalar>(
389 A, fwdSolveCriteria.get(), throwOnFwdSolveFailure,
390 adjSolveCriteria.get(), throwOnAdjSolveFailure
403 #define THYRA_DEFAULT_INVERSE_LINEAR_OP_INSTANT(SCALAR) \ 405 template class DefaultInverseLinearOp<SCALAR >; \ 407 template RCP<LinearOpBase<SCALAR > > \ 409 const RCP<LinearOpWithSolveBase<SCALAR > > &A, \ 410 const Ptr<const SolveCriteria<SCALAR > > &fwdSolveCriteria, \ 411 const EThrowOnSolveFailure throwOnFwdSolveFailure, \ 412 const Ptr<const SolveCriteria<SCALAR > > &adjSolveCriteria, \ 413 const EThrowOnSolveFailure throwOnAdjSolveFailure \ 416 template RCP<LinearOpBase<SCALAR > > \ 418 const RCP<const LinearOpWithSolveBase<SCALAR > > &A, \ 419 const Ptr<const SolveCriteria<SCALAR > > &fwdSolveCriteria, \ 420 const EThrowOnSolveFailure throwOnFwdSolveFailure, \ 421 const Ptr<const SolveCriteria<SCALAR > > &adjSolveCriteria, \ 422 const EThrowOnSolveFailure throwOnAdjSolveFailure \ 426 #endif // THYRA_DEFAULT_INVERSE_LINEAR_OP_DEF_HPP Base class for all linear operators that can support a high-level solve operation.
RCP< const VectorSpaceBase< Scalar > > domain() const
Returns this->getLows()->range() if <t>this->getLows().get()!=NULL and returns Teuchos::null otherwis...
EOpTransp
Enumeration for determining how a linear operator is applied. `*.
void uninitialize()
Set to uninitialized.
DefaultInverseLinearOp()
Constructs to uninitialized (see postconditions for uninitialize()).
Use the non-transposed operator.
EOpTransp real_trans(EOpTransp transp)
Return NOTRANS or TRANS for real scalar valued operators and this also is used for determining struct...
bool opSupportedImpl(EOpTransp M_trans) const
Returns true only if all constituent operators support M_trans.
RCP< const LinearOpBase< Scalar > > clone() const
std::string description() const
RCP< const LinearOpWithSolveBase< Scalar > > getLows() const
RCP< const VectorSpaceBase< Scalar > > range() const
Returns this->getLows()->domain() if <t>this->getLows().get()!=NULL and returns Teuchos::null otherwi...
Throw an exception if a solve fails to converge.
Interface for a collection of column vectors called a multi-vector.
void describe(FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
void applyImpl(const EOpTransp M_trans, const MultiVectorBase< Scalar > &X, const Ptr< MultiVectorBase< Scalar > > &Y, const Scalar alpha, const Scalar beta) const
Simple struct for the return status from a solve.
void initialize(const RCP< LinearOpWithSolveBase< Scalar > > &lows, const SolveCriteria< Scalar > *fwdSolveCriteria=NULL, const EThrowOnSolveFailure throwOnFwdSolveFailure=THROW_ON_SOLVE_FAILURE, const SolveCriteria< Scalar > *adjSolveCriteria=NULL, const EThrowOnSolveFailure throwOnAdjSolveFailure=THROW_ON_SOLVE_FAILURE)
Initialize given a non-const LinearOpWithSolveBase object and an optional .
EThrowOnSolveFailure
Determines what to do if inverse solve fails.
The requested solution criteria has likely been achieved.
Exception type thrown on an catastrophic solve failure.
const char * toString(EConj conj)
Return a string name for a EOpTransp value. `*.
Simple struct that defines the requested solution criteria for a solve.
Concrete LinearOpBase subclass that creates an implicit LinearOpBase object using the inverse action ...
RCP< LinearOpWithSolveBase< Scalar > > getNonconstLows()