42 #ifndef THYRA_TPETRA_VECTOR_HPP
43 #define THYRA_TPETRA_VECTOR_HPP
55 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
60 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
63 const RCP<Tpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVector
66 initializeImpl(tpetraVectorSpace, tpetraVector);
70 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
73 const RCP<
const Tpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVector
76 initializeImpl(tpetraVectorSpace, tpetraVector);
80 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
84 return tpetraVector_.getNonconstObj();
88 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
97 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
101 if (domainSpace_.is_null()) {
102 domainSpace_ = tpetraVectorSpace<Scalar>(
103 Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal,Node>(
105 tpetraVector_.getConstObj()->getMap()->getComm()
116 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
120 return tpetraVectorSpace_;
127 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
129 const Ptr<ArrayRCP<Scalar> > &localValues )
131 *localValues = tpetraVector_.getNonconstObj()->get1dViewNonConst();
135 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
137 const Ptr<ArrayRCP<const Scalar> > &localValues )
const
139 *localValues = tpetraVector_->get1dView();
146 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
154 if (!tpetraVector_.getNonconstObj()->isDistributed()) {
155 auto comm = tpetraVector_.getNonconstObj()->getMap()->getComm();
156 if (tpetraVector_.getConstObj()->getMap()->getComm()->getRank() == 0)
157 tpetraVector_.getNonconstObj()->randomize(l, u);
160 tpetraVector_.getNonconstObj()->reduce();
162 tpetraVector_.getNonconstObj()->randomize(l, u);
167 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
169 const VectorBase<Scalar>& x
172 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
177 tpetraVector_.getNonconstObj()->abs(*tx);
179 VectorDefaultBase<Scalar>::absImpl(x);
184 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
186 const VectorBase<Scalar>& x
189 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
194 tpetraVector_.getNonconstObj()->reciprocal(*tx);
196 VectorDefaultBase<Scalar>::reciprocalImpl(x);
201 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
203 const VectorBase<Scalar>& x
206 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
212 tpetraVector_.getNonconstObj()->elementWiseMultiply(
213 ST::one(), *tx, *tpetraVector_.getConstObj(), ST::zero());
215 VectorDefaultBase<Scalar>::eleWiseScaleImpl(x);
220 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
223 const VectorBase<Scalar>& x
226 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
234 = Tpetra::createVector<Scalar>(tx->getMap());
235 temp->elementWiseMultiply(
236 ST::one(), *tx, *tpetraVector_.getConstObj(), ST::zero());
237 return ST::magnitude(ST::squareroot(tpetraVector_.getConstObj()->dot(*temp)));
239 return VectorDefaultBase<Scalar>::norm2WeightedImpl(x);
244 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
247 const ArrayView<
const Ptr<
const VectorBase<Scalar> > > &vecs,
248 const ArrayView<
const Ptr<VectorBase<Scalar> > > &targ_vecs,
249 const Ptr<RTOpPack::ReductTarget> &reduct_obj,
253 SpmdVectorDefaultBase<Scalar>::applyOpImpl(op, vecs, targ_vecs, reduct_obj, global_offset);
257 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
264 SpmdVectorDefaultBase<Scalar>::acquireDetachedVectorViewImpl(rng, sub_vec);
268 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
276 SpmdVectorDefaultBase<Scalar>::acquireNonconstDetachedVectorViewImpl(rng, sub_vec);
280 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
286 SpmdVectorDefaultBase<Scalar>::commitNonconstDetachedVectorViewImpl(sub_vec);
293 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
296 tpetraVector_.getNonconstObj()->putScalar(alpha);
300 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
304 auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
309 tpetraVector_.getNonconstObj()->assign(*tmv);
311 MultiVectorDefaultBase<Scalar>::assignMultiVecImpl(mv);
316 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
319 tpetraVector_.getNonconstObj()->scale(alpha);
323 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
326 const MultiVectorBase<Scalar>& mv
329 auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
335 tpetraVector_.getNonconstObj()->update(alpha, *tmv, ST::one());
337 MultiVectorDefaultBase<Scalar>::updateImpl(alpha, mv);
342 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
345 const ArrayView<
const Ptr<
const MultiVectorBase<Scalar> > >& mv,
356 bool allCastsSuccessful =
true;
358 auto mvIter = mv.begin();
359 auto tmvIter = tmvs.begin();
360 for (; mvIter != mv.end(); ++mvIter, ++tmvIter) {
361 tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromPtr(*mvIter));
365 allCastsSuccessful =
false;
373 auto len = mv.size();
375 tpetraVector_.getNonconstObj()->scale(beta);
376 }
else if (len == 1 && allCastsSuccessful) {
377 tpetraVector_.getNonconstObj()->update(alpha[0], *tmvs[0], beta);
378 }
else if (len == 2 && allCastsSuccessful) {
379 tpetraVector_.getNonconstObj()->update(alpha[0], *tmvs[0], alpha[1], *tmvs[1], beta);
380 }
else if (allCastsSuccessful) {
382 auto tmvIter = tmvs.begin();
383 auto alphaIter = alpha.
begin();
388 for (; tmvIter != tmvs.end(); ++tmvIter) {
389 if (tmvIter->getRawPtr() == tpetraVector_.getConstObj().getRawPtr()) {
397 tmvIter = tmvs.begin();
401 if ((tmvs.size() % 2) == 0) {
402 tpetraVector_.getNonconstObj()->scale(beta);
404 tpetraVector_.getNonconstObj()->update(*alphaIter, *(*tmvIter), beta);
408 for (; tmvIter != tmvs.end(); tmvIter+=2, alphaIter+=2) {
409 tpetraVector_.getNonconstObj()->update(
410 *alphaIter, *(*tmvIter), *(alphaIter+1), *(*(tmvIter+1)), ST::one());
413 MultiVectorDefaultBase<Scalar>::linearCombinationImpl(alpha, mv, beta);
418 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
420 const MultiVectorBase<Scalar>& mv,
424 auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
429 tpetraVector_.getConstObj()->dot(*tmv, prods);
431 MultiVectorDefaultBase<Scalar>::dotsImpl(mv, prods);
436 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
438 const ArrayView<
typename ScalarTraits<Scalar>::magnitudeType>& norms
441 tpetraVector_.getConstObj()->norm1(norms);
445 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
447 const ArrayView<
typename ScalarTraits<Scalar>::magnitudeType>& norms
450 tpetraVector_.getConstObj()->norm2(norms);
454 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
456 const ArrayView<
typename ScalarTraits<Scalar>::magnitudeType>& norms
459 tpetraVector_.getConstObj()->normInf(norms);
463 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
465 const EOpTransp M_trans,
466 const MultiVectorBase<Scalar> &X,
467 const Ptr<MultiVectorBase<Scalar> > &Y,
473 typedef Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> TMV;
483 "Error, conjugation without transposition is not allowed for complex scalar types!");
501 Y_tpetra->multiply(trans,
Teuchos::NO_TRANS, alpha, *tpetraVector_.getConstObj(), *X_tpetra, beta);
504 VectorDefaultBase<Scalar>::applyImpl(M_trans, X, Y, alpha, beta);
513 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
514 template<
class TpetraVector_t>
524 tpetraVectorSpace_ = tpetraVectorSpace;
525 tpetraVector_.initialize(tpetraVector);
526 this->updateSpmdSpace();
530 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
535 using Teuchos::rcp_dynamic_cast;
539 RCP<TMV> tmv = rcp_dynamic_cast<TMV>(mv);
541 return tmv->getTpetraMultiVector();
544 RCP<TV> tv = rcp_dynamic_cast<TV>(mv);
546 return tv->getTpetraVector();
553 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
558 using Teuchos::rcp_dynamic_cast;
564 return tmv->getConstTpetraMultiVector();
569 return tv->getConstTpetraVector();
576 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
582 RCP<TV> tv = Teuchos::rcp_dynamic_cast<TV>(v);
584 return tv->getTpetraVector();
591 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
599 return tv->getConstTpetraVector();
609 #endif // THYRA_TPETRA_VECTOR_HPP
void constInitialize(const RCP< const TpetraVectorSpace< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVectorSpace, const RCP< const Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVector)
Initialize.
Concrete implementation of Thyra::MultiVector in terms of Tpetra::MultiVector.
RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetraMultiVector(const RCP< MultiVectorBase< Scalar > > &mv) const
void applyImpl(const EOpTransp M_trans, const MultiVectorBase< Scalar > &X, const Ptr< MultiVectorBase< Scalar > > &Y, const Scalar alpha, const Scalar beta) const
Concrete implementation of an SPMD vector space for Tpetra.
void commitNonconstDetachedVectorViewImpl(RTOpPack::SubVectorView< Scalar > *sub_vec)
void initialize(const RCP< const TpetraVectorSpace< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVectorSpace, const RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVector)
Initialize.
RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getConstTpetraMultiVector(const RCP< const MultiVectorBase< Scalar > > &mv) const
TpetraVector()
Construct to uninitialized.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual void norms1Impl(const ArrayView< typename ScalarTraits< Scalar >::magnitudeType > &norms) const
virtual void applyOpImpl(const RTOpPack::RTOpT< Scalar > &op, const ArrayView< const Ptr< const VectorBase< Scalar > > > &vecs, const ArrayView< const Ptr< VectorBase< Scalar > > > &targ_vecs, const Ptr< RTOpPack::ReductTarget > &reduct_obj, const Ordinal global_offset) const
virtual void norms2Impl(const ArrayView< typename ScalarTraits< Scalar >::magnitudeType > &norms) const
void getLocalVectorDataImpl(const Ptr< ArrayRCP< const Scalar > > &localValues) const
virtual void assignMultiVecImpl(const MultiVectorBase< Scalar > &mv)
void acquireNonconstDetachedVectorViewImpl(const Range1D &rng, RTOpPack::SubVectorView< Scalar > *sub_vec)
RCP< const Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getConstTpetraVector() const
Get the embedded non-const Tpetra::Vector.
virtual void reciprocalImpl(const VectorBase< Scalar > &x)
void acquireDetachedVectorViewImpl(const Range1D &rng, RTOpPack::ConstSubVectorView< Scalar > *sub_vec) const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
virtual Teuchos::ScalarTraits< Scalar >::magnitudeType norm2WeightedImpl(const VectorBase< Scalar > &x) const
virtual void linearCombinationImpl(const ArrayView< const Scalar > &alpha, const ArrayView< const Ptr< const MultiVectorBase< Scalar > > > &mv, const Scalar &beta)
virtual void absImpl(const VectorBase< Scalar > &x)
void getNonconstLocalVectorDataImpl(const Ptr< ArrayRCP< Scalar > > &localValues)
Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraMultiVector_t
RCP< const SpmdVectorSpaceBase< Scalar > > spmdSpaceImpl() const
void initializeImpl(const RCP< const TpetraVectorSpace< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVectorSpace, const RCP< TpetraVector_t > &tpetraVector)
Concrete Thyra::SpmdVectorBase using Tpetra::Vector.
virtual void dotsImpl(const MultiVectorBase< Scalar > &mv, const ArrayView< Scalar > &prods) const
RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetraVector()
Get the embedded non-const Tpetra::Vector.
bool nonnull(const boost::shared_ptr< T > &p)
virtual void assignImpl(Scalar alpha)
virtual void eleWiseScaleImpl(const VectorBase< Scalar > &x)
virtual void normsInfImpl(const ArrayView< typename ScalarTraits< Scalar >::magnitudeType > &norms) const
RCP< const VectorSpaceBase< Scalar > > domain() const
virtual void randomizeImpl(Scalar l, Scalar u)
virtual void updateImpl(Scalar alpha, const MultiVectorBase< Scalar > &mv)
virtual void scaleImpl(Scalar alpha)
#define TEUCHOS_ASSERT(assertion_test)
#define TEUCHOS_ASSERT_EQUALITY(val1, val2)