42 #ifndef THYRA_MULTI_VECTOR_LINEAR_OP_WITH_SOLVE_HPP 43 #define THYRA_MULTI_VECTOR_LINEAR_OP_WITH_SOLVE_HPP 45 #include "Thyra_DefaultMultiVectorLinearOpWithSolve_decl.hpp" 46 #include "Thyra_DefaultDiagonalLinearOp.hpp" 47 #include "Thyra_LinearOpWithSolveBase.hpp" 48 #include "Thyra_DefaultMultiVectorProductVectorSpace.hpp" 49 #include "Thyra_DefaultMultiVectorProductVector.hpp" 50 #include "Thyra_AssertOp.hpp" 51 #include "Teuchos_dyn_cast.hpp" 60 template<
class Scalar>
65 template<
class Scalar>
72 validateInitialize(lows,multiVecRange,multiVecDomain);
74 multiVecRange_ = multiVecRange;
75 multiVecDomain_ = multiVecDomain;
79 template<
class Scalar>
86 validateInitialize(lows,multiVecRange,multiVecDomain);
88 multiVecRange_ = multiVecRange;
89 multiVecDomain_ = multiVecDomain;
93 template<
class Scalar>
94 RCP<LinearOpWithSolveBase<Scalar> >
97 return lows_.getNonconstObj();
101 template<
class Scalar>
102 RCP<const LinearOpWithSolveBase<Scalar> >
105 return lows_.getConstObj();
109 template<
class Scalar>
112 lows_.uninitialize();
113 multiVecRange_ = Teuchos::null;
114 multiVecDomain_ = Teuchos::null;
121 template<
class Scalar>
122 RCP< const VectorSpaceBase<Scalar> >
125 return multiVecRange_;
129 template<
class Scalar>
130 RCP< const VectorSpaceBase<Scalar> >
133 return multiVecDomain_;
137 template<
class Scalar>
138 RCP<const LinearOpBase<Scalar> >
141 return Teuchos::null;
151 template<
class Scalar>
156 return Thyra::opSupported(*lows_.getConstObj(),M_trans);
160 template<
class Scalar>
170 using Teuchos::dyn_cast;
175 for (
Ordinal col_j = 0; col_j < numCols; ++col_j) {
177 const RCP<const VectorBase<Scalar> > x = XX.
col(col_j);
178 const RCP<VectorBase<Scalar> > y = YY->col(col_j);
180 RCP<const MultiVectorBase<Scalar> >
181 X = dyn_cast<
const MVPV>(*x).getMultiVector().assert_not_null();
182 RCP<MultiVectorBase<Scalar> >
183 Y = dyn_cast<MVPV>(*y).getNonconstMultiVector().assert_not_null();
185 Thyra::apply( *lows_.getConstObj(), M_trans, *X, Y.ptr(), alpha, beta );
195 template<
class Scalar>
201 return Thyra::solveSupports(*lows_.getConstObj(),M_trans);
205 template<
class Scalar>
211 return Thyra::solveSupportsSolveMeasureType(
212 *lows_.getConstObj(),M_trans,solveMeasureType);
216 template<
class Scalar>
226 using Teuchos::dyn_cast;
227 using Teuchos::outArg;
228 using Teuchos::inOutArg;
234 accumulateSolveStatusInit(outArg(overallSolveStatus));
236 for (
Ordinal col_j = 0; col_j < numCols; ++col_j) {
238 const RCP<const VectorBase<Scalar> > b = BB.
col(col_j);
239 const RCP<VectorBase<Scalar> > x = XX->col(col_j);
241 RCP<const MultiVectorBase<Scalar> >
242 B = dyn_cast<
const MVPV>(*b).getMultiVector().assert_not_null();
243 RCP<MultiVectorBase<Scalar> >
244 X = dyn_cast<MVPV>(*x).getNonconstMultiVector().assert_not_null();
247 Thyra::solve(*lows_.getConstObj(), transp, *B, X.ptr(), solveCriteria);
249 accumulateSolveStatus(
251 solveStatus, inOutArg(overallSolveStatus) );
255 return overallSolveStatus;
263 template<
class Scalar>
271 TEUCHOS_TEST_FOR_EXCEPT(is_null(lows));
272 TEUCHOS_TEST_FOR_EXCEPT(is_null(multiVecRange));
273 TEUCHOS_TEST_FOR_EXCEPT(is_null(multiVecDomain));
274 TEUCHOS_TEST_FOR_EXCEPT( multiVecRange->numBlocks() != multiVecDomain->numBlocks() );
275 if (lows->range() != Teuchos::null)
277 "DefaultMultiVectorLinearOpWithSolve<Scalar>::initialize(lows,multiVecRange,multiVecDomain)",
278 *lows->range(), *multiVecRange->getBlock(0) );
279 if (lows->domain() != Teuchos::null)
281 "DefaultMultiVectorLinearOpWithSolve<Scalar>::initialize(lows,multiVecRange,multiVecDomain)",
282 *lows->domain(), *multiVecDomain->getBlock(0) );
290 #endif // THYRA_MULTI_VECTOR_LINEAR_OP_WITH_SOLVE_HPP #define THYRA_ASSERT_VEC_SPACES(FUNC_NAME, VS1, VS2)
This is a very useful macro that should be used to validate that two vector spaces are compatible...
Implicit concrete LinearOpWithSolveBase subclass that takes a flattended out multi-vector and perform...
Base class for all linear operators that can support a high-level solve operation.
EOpTransp
Enumeration for determining how a linear operator is applied. `*.
bool solveSupportsImpl(EOpTransp M_trans) const
DefaultMultiVectorLinearOpWithSolve()
Construct to uninitialized.
bool opSupportedImpl(EOpTransp M_trans) const
bool solveSupportsSolveMeasureTypeImpl(EOpTransp M_trans, const SolveMeasureType &solveMeasureType) const
Concrete implementation of a product vector which is really composed out of the columns of a multi-ve...
Teuchos::Ordinal Ordinal
Type for the dimension of a vector space. `*.
Interface for a collection of column vectors called a multi-vector.
RCP< const VectorBase< Scalar > > col(Ordinal j) const
Calls colImpl().
void applyImpl(const EOpTransp M_trans, const MultiVectorBase< Scalar > &X, const Ptr< MultiVectorBase< Scalar > > &Y, const Scalar alpha, const Scalar beta) const
void nonconstInitialize(const RCP< LinearOpWithSolveBase< Scalar > > &lows, const RCP< const DefaultMultiVectorProductVectorSpace< Scalar > > &multiVecRange, const RCP< const DefaultMultiVectorProductVectorSpace< Scalar > > &multiVecDomain)
RCP< const LinearOpWithSolveBase< Scalar > > getLinearOpWithSolve() const
Simple struct for the return status from a solve.
Standard concrete implementation of a product vector space that creates product vectors fromed implic...
virtual RCP< const VectorSpaceBase< Scalar > > domain() const =0
Return a smart pointer for the domain space for this operator.
RCP< const LinearOpBase< Scalar > > clone() const
RCP< const VectorSpaceBase< Scalar > > domain() const
RCP< LinearOpWithSolveBase< Scalar > > getNonconstLinearOpWithSolve()
SolveStatus< Scalar > solveImpl(const EOpTransp transp, const MultiVectorBase< Scalar > &B, const Ptr< MultiVectorBase< Scalar > > &X, const Ptr< const SolveCriteria< Scalar > > solveCriteria) const
Simple struct that defines the requested solution criteria for a solve.
void initialize(const RCP< const LinearOpWithSolveBase< Scalar > > &lows, const RCP< const DefaultMultiVectorProductVectorSpace< Scalar > > &multiVecRange, const RCP< const DefaultMultiVectorProductVectorSpace< Scalar > > &multiVecDomain)
RCP< const VectorSpaceBase< Scalar > > range() const