42#ifndef THYRA_TPETRA_VECTOR_HPP
43#define THYRA_TPETRA_VECTOR_HPP
55template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
60template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
63 const RCP<Tpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVector
66 initializeImpl(tpetraVectorSpace, tpetraVector);
70template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
73 const RCP<
const Tpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVector
76 initializeImpl(tpetraVectorSpace, tpetraVector);
80template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
84 return tpetraVector_.getNonconstObj();
88template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
97template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
101 if (domainSpace_.is_null()) {
103 Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal,Node>(
105 tpetraVector_.getConstObj()->getMap()->getComm()
116template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
120 return tpetraVectorSpace_;
127template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
131 *
localValues = tpetraVector_.getNonconstObj()->get1dViewNonConst();
135template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
146template <
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);
167template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
172 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
177 tpetraVector_.getNonconstObj()->abs(*
tx);
184template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
189 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
194 tpetraVector_.getNonconstObj()->reciprocal(*
tx);
201template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
206 auto tx = this->getConstTpetraVector(Teuchos::rcpFromRef(x));
212 tpetraVector_.getNonconstObj()->elementWiseMultiply(
213 ST::one(), *
tx, *tpetraVector_.getConstObj(), ST::zero());
220template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
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)));
244template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
257template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
264 SpmdVectorDefaultBase<Scalar>::acquireDetachedVectorViewImpl(
rng,
sub_vec);
268template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
276 SpmdVectorDefaultBase<Scalar>::acquireNonconstDetachedVectorViewImpl(
rng,
sub_vec);
280template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
286 SpmdVectorDefaultBase<Scalar>::commitNonconstDetachedVectorViewImpl(
sub_vec);
293template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
296 tpetraVector_.getNonconstObj()->putScalar(alpha);
300template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
304 auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(
mv));
309 tpetraVector_.getNonconstObj()->assign(*
tmv);
316template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
319 tpetraVector_.getNonconstObj()->scale(alpha);
323template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
329 auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(
mv));
335 tpetraVector_.getNonconstObj()->update(alpha, *
tmv, ST::one());
342template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
361 tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromPtr(*
mvIter));
373 auto len =
mv.size();
375 tpetraVector_.getNonconstObj()->scale(beta);
377 tpetraVector_.getNonconstObj()->update(alpha[0], *
tmvs[0], beta);
379 tpetraVector_.getNonconstObj()->update(alpha[0], *
tmvs[0], alpha[1], *
tmvs[1], beta);
401 if ((
tmvs.size() % 2) == 0) {
402 tpetraVector_.getNonconstObj()->scale(beta);
409 tpetraVector_.getNonconstObj()->update(
418template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
424 auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(
mv));
429 tpetraVector_.getConstObj()->dot(*
tmv,
prods);
436template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
441 tpetraVector_.getConstObj()->norm1(
norms);
445template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
450 tpetraVector_.getConstObj()->norm2(
norms);
454template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
459 tpetraVector_.getConstObj()->normInf(
norms);
463template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
473 typedef Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>
TMV;
483 "Error, conjugation without transposition is not allowed for complex scalar types!");
513template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
514template<
class TpetraVector_t>
524 tpetraVectorSpace_ = tpetraVectorSpace;
525 tpetraVector_.initialize(tpetraVector);
530template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
535 using Teuchos::rcp_dynamic_cast;
541 return tmv->getTpetraMultiVector();
546 return tv->getTpetraVector();
549 return Teuchos::null;
553template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
558 using Teuchos::rcp_dynamic_cast;
564 return tmv->getConstTpetraMultiVector();
569 return tv->getConstTpetraVector();
572 return Teuchos::null;
576template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
584 return tv->getTpetraVector();
586 return Teuchos::null;
591template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
599 return tv->getConstTpetraVector();
601 return Teuchos::null;
bool is_null(const RCP< T > &p)
RCP(T *p, const RCPNodeHandle &node)
bool nonnull(const RCP< T > &p)
Concrete Thyra::SpmdVectorBase using Tpetra::Vector.
Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraMultiVector_t
TpetraVector()
Construct to uninitialized.
#define TEUCHOS_ASSERT(assertion_test)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
#define TEUCHOS_ASSERT_EQUALITY(val1, val2)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)