10 #ifndef THYRA_TPETRA_VECTOR_HPP
11 #define THYRA_TPETRA_VECTOR_HPP
16 #include "Kokkos_Core.hpp"
24 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
29 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
32 const RCP<Tpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVector
35 initializeImpl(tpetraVectorSpace, tpetraVector);
39 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
42 const RCP<
const Tpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVector
45 initializeImpl(tpetraVectorSpace, tpetraVector);
49 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
53 return tpetraVector_.getNonconstObj();
57 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
66 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
70 if (domainSpace_.is_null()) {
71 domainSpace_ = tpetraVectorSpace<Scalar>(
72 Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal,Node>(
74 tpetraVector_.getConstObj()->getMap()->getComm()
85 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
89 return tpetraVectorSpace_;
96 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
98 const Ptr<ArrayRCP<Scalar> > &localValues )
100 *localValues = tpetraVector_.getNonconstObj()->get1dViewNonConst();
104 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
106 const Ptr<ArrayRCP<const Scalar> > &localValues )
const
108 *localValues = tpetraVector_->get1dView();
115 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
123 if (!tpetraVector_.getNonconstObj()->isDistributed()) {
124 auto comm = tpetraVector_.getNonconstObj()->getMap()->getComm();
125 if (tpetraVector_.getConstObj()->getMap()->getComm()->getRank() == 0)
126 tpetraVector_.getNonconstObj()->randomize(l, u);
129 tpetraVector_.getNonconstObj()->reduce();
131 tpetraVector_.getNonconstObj()->randomize(l, u);
136 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
138 const VectorBase<Scalar>& x
141 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
146 tpetraVector_.getNonconstObj()->abs(*tx);
148 VectorDefaultBase<Scalar>::absImpl(x);
153 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
155 const VectorBase<Scalar>& x
158 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
163 tpetraVector_.getNonconstObj()->reciprocal(*tx);
165 VectorDefaultBase<Scalar>::reciprocalImpl(x);
170 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
172 const VectorBase<Scalar>& x
175 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
181 tpetraVector_.getNonconstObj()->elementWiseMultiply(
182 ST::one(), *tx, *tpetraVector_.getConstObj(), ST::zero());
184 VectorDefaultBase<Scalar>::eleWiseScaleImpl(x);
189 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
192 const VectorBase<Scalar>& x
195 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
203 = Tpetra::createVector<Scalar>(tx->getMap());
204 temp->elementWiseMultiply(
205 ST::one(), *tx, *tpetraVector_.getConstObj(), ST::zero());
206 return ST::magnitude(ST::squareroot(tpetraVector_.getConstObj()->dot(*temp)));
208 return VectorDefaultBase<Scalar>::norm2WeightedImpl(x);
213 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
216 const ArrayView<
const Ptr<
const VectorBase<Scalar> > > &vecs,
217 const ArrayView<
const Ptr<VectorBase<Scalar> > > &targ_vecs,
218 const Ptr<RTOpPack::ReductTarget> &reduct_obj,
222 SpmdVectorDefaultBase<Scalar>::applyOpImpl(op, vecs, targ_vecs, reduct_obj, global_offset);
226 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
233 SpmdVectorDefaultBase<Scalar>::acquireDetachedVectorViewImpl(rng, sub_vec);
237 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
245 SpmdVectorDefaultBase<Scalar>::acquireNonconstDetachedVectorViewImpl(rng, sub_vec);
249 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
255 SpmdVectorDefaultBase<Scalar>::commitNonconstDetachedVectorViewImpl(sub_vec);
262 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
265 tpetraVector_.getNonconstObj()->putScalar(alpha);
269 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
273 auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
278 tpetraVector_.getNonconstObj()->assign(*tmv);
280 MultiVectorDefaultBase<Scalar>::assignMultiVecImpl(mv);
285 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
288 tpetraVector_.getNonconstObj()->scale(alpha);
292 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
295 const MultiVectorBase<Scalar>& mv
298 auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
304 tpetraVector_.getNonconstObj()->update(alpha, *tmv, ST::one());
306 MultiVectorDefaultBase<Scalar>::updateImpl(alpha, mv);
311 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
314 const ArrayView<
const Ptr<
const MultiVectorBase<Scalar> > >& mv,
325 bool allCastsSuccessful =
true;
327 auto mvIter = mv.begin();
328 auto tmvIter = tmvs.begin();
329 for (; mvIter != mv.end(); ++mvIter, ++tmvIter) {
330 tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromPtr(*mvIter));
334 allCastsSuccessful =
false;
342 auto len = mv.size();
344 tpetraVector_.getNonconstObj()->scale(beta);
345 }
else if (len == 1 && allCastsSuccessful) {
346 tpetraVector_.getNonconstObj()->update(alpha[0], *tmvs[0], beta);
347 }
else if (len == 2 && allCastsSuccessful) {
348 tpetraVector_.getNonconstObj()->update(alpha[0], *tmvs[0], alpha[1], *tmvs[1], beta);
349 }
else if (allCastsSuccessful) {
351 auto tmvIter = tmvs.begin();
352 auto alphaIter = alpha.
begin();
357 for (; tmvIter != tmvs.end(); ++tmvIter) {
358 if (tmvIter->getRawPtr() == tpetraVector_.getConstObj().getRawPtr()) {
366 tmvIter = tmvs.begin();
370 if ((tmvs.size() % 2) == 0) {
371 tpetraVector_.getNonconstObj()->scale(beta);
373 tpetraVector_.getNonconstObj()->update(*alphaIter, *(*tmvIter), beta);
377 for (; tmvIter != tmvs.end(); tmvIter+=2, alphaIter+=2) {
378 tpetraVector_.getNonconstObj()->update(
379 *alphaIter, *(*tmvIter), *(alphaIter+1), *(*(tmvIter+1)), ST::one());
382 MultiVectorDefaultBase<Scalar>::linearCombinationImpl(alpha, mv, beta);
387 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
389 const MultiVectorBase<Scalar>& mv,
393 auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
398 tpetraVector_.getConstObj()->dot(*tmv, prods);
400 MultiVectorDefaultBase<Scalar>::dotsImpl(mv, prods);
405 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
407 const ArrayView<
typename ScalarTraits<Scalar>::magnitudeType>& norms
410 tpetraVector_.getConstObj()->norm1(norms);
414 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
416 const ArrayView<
typename ScalarTraits<Scalar>::magnitudeType>& norms
419 tpetraVector_.getConstObj()->norm2(norms);
423 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
425 const ArrayView<
typename ScalarTraits<Scalar>::magnitudeType>& norms
428 tpetraVector_.getConstObj()->normInf(norms);
432 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
434 const EOpTransp M_trans,
435 const MultiVectorBase<Scalar> &X,
436 const Ptr<MultiVectorBase<Scalar> > &Y,
442 typedef Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> TMV;
452 "Error, conjugation without transposition is not allowed for complex scalar types!");
470 Y_tpetra->multiply(trans,
Teuchos::NO_TRANS, alpha, *tpetraVector_.getConstObj(), *X_tpetra, beta);
473 VectorDefaultBase<Scalar>::applyImpl(M_trans, X, Y, alpha, beta);
482 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
483 template<
class TpetraVector_t>
493 tpetraVectorSpace_ = tpetraVectorSpace;
494 tpetraVector_.initialize(tpetraVector);
495 this->updateSpmdSpace();
499 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
504 using Teuchos::rcp_dynamic_cast;
508 RCP<TMV> tmv = rcp_dynamic_cast<TMV>(mv);
510 return tmv->getTpetraMultiVector();
513 RCP<TV> tv = rcp_dynamic_cast<TV>(mv);
515 return tv->getTpetraVector();
522 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
527 using Teuchos::rcp_dynamic_cast;
533 return tmv->getConstTpetraMultiVector();
538 return tv->getConstTpetraVector();
545 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
551 RCP<TV> tv = Teuchos::rcp_dynamic_cast<TV>(v);
553 return tv->getTpetraVector();
560 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
568 return tv->getConstTpetraVector();
578 #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)