42 #ifndef THYRA_TPETRA_MULTIVECTOR_HPP 43 #define THYRA_TPETRA_MULTIVECTOR_HPP 48 #include "Teuchos_Assert.hpp" 57 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
62 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
65 const RCP<
const ScalarProdVectorSpaceBase<Scalar> > &domainSpace,
66 const RCP<Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraMultiVector
69 initializeImpl(tpetraVectorSpace, domainSpace, tpetraMultiVector);
73 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
76 const RCP<
const ScalarProdVectorSpaceBase<Scalar> > &domainSpace,
77 const RCP<
const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraMultiVector
80 initializeImpl(tpetraVectorSpace, domainSpace, tpetraMultiVector);
84 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
85 RCP<Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
88 return tpetraMultiVector_.getNonconstObj();
92 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
93 RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
96 return tpetraMultiVector_;
103 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
104 RCP< const ScalarProdVectorSpaceBase<Scalar> >
114 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
118 tpetraMultiVector_.getNonconstObj()->putScalar(alpha);
122 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
123 RCP<const VectorBase<Scalar> >
127 TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE(j, 0, this->domain()->dim());
129 return constTpetraVector<Scalar>(
131 tpetraMultiVector_->getVector(j)
136 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
137 RCP<VectorBase<Scalar> >
141 TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE(j, 0, this->domain()->dim());
143 return tpetraVector<Scalar>(
145 tpetraMultiVector_.getNonconstObj()->getVectorNonConst(j)
150 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
151 RCP<const MultiVectorBase<Scalar> >
153 const Range1D& col_rng_in
156 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT 157 std::cerr <<
"\nTpetraMultiVector::subView(Range1D) const called!\n";
159 const Range1D colRng = this->validateColRange(col_rng_in);
161 const RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tpetraView =
162 this->getConstTpetraMultiVector()->subView(colRng);
164 const RCP<const ScalarProdVectorSpaceBase<Scalar> > viewDomainSpace =
165 tpetraVectorSpace<Scalar>(
166 Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal>(
167 tpetraView->getNumVectors(),
168 tpetraView->getMap()->getComm(),
169 tpetraView->getMap()->getNode()
173 return constTpetraMultiVector(
181 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
182 RCP<MultiVectorBase<Scalar> >
184 const Range1D& col_rng_in
187 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT 188 std::cerr <<
"\nTpetraMultiVector::subView(Range1D) called!\n";
190 const Range1D colRng = this->validateColRange(col_rng_in);
192 const RCP<Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tpetraView =
193 this->getTpetraMultiVector()->subViewNonConst(colRng);
195 const RCP<const ScalarProdVectorSpaceBase<Scalar> > viewDomainSpace =
196 tpetraVectorSpace<Scalar>(
197 Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal>(
198 tpetraView->getNumVectors(),
199 tpetraView->getMap()->getComm(),
200 tpetraView->getMap()->getNode()
204 return tpetraMultiVector(
259 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
260 RCP<const SpmdVectorSpaceBase<Scalar> >
263 return tpetraVectorSpace_;
267 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
269 const Ptr<ArrayRCP<Scalar> > &localValues,
const Ptr<Ordinal> &leadingDim
272 *localValues = tpetraMultiVector_.getNonconstObj()->get1dViewNonConst();
273 *leadingDim = tpetraMultiVector_->getStride();
277 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
279 const Ptr<ArrayRCP<const Scalar> > &localValues,
const Ptr<Ordinal> &leadingDim
282 *localValues = tpetraMultiVector_->get1dView();
283 *leadingDim = tpetraMultiVector_->getStride();
290 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
291 template<
class TpetraMultiVector_t>
294 const RCP<
const ScalarProdVectorSpaceBase<Scalar> > &domainSpace,
295 const RCP<TpetraMultiVector_t> &tpetraMultiVector
299 TEUCHOS_ASSERT(nonnull(tpetraVectorSpace));
300 TEUCHOS_ASSERT(nonnull(domainSpace));
301 TEUCHOS_ASSERT(nonnull(tpetraMultiVector));
305 tpetraVectorSpace_ = tpetraVectorSpace;
306 domainSpace_ = domainSpace;
307 tpetraMultiVector_.initialize(tpetraMultiVector);
308 this->updateSpmdSpace();
316 #endif // THYRA_TPETRA_MULTIVECTOR_HPP RCP< const ScalarProdVectorSpaceBase< Scalar > > domainScalarProdVecSpc() const
TpetraMultiVector()
Construct to uninitialized.
void getNonconstLocalMultiVectorDataImpl(const Ptr< ArrayRCP< Scalar > > &localValues, const Ptr< Ordinal > &leadingDim)
Concrete implementation of an SPMD vector space for Tpetra.
RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetraMultiVector()
Extract the underlying non-const Tpetra::MultiVector object.
void assignImpl(Scalar alpha)
RCP< const SpmdVectorSpaceBase< Scalar > > spmdSpaceImpl() const
RCP< MultiVectorBase< Scalar > > nonconstContigSubViewImpl(const Range1D &colRng)
void getLocalMultiVectorDataImpl(const Ptr< ArrayRCP< const Scalar > > &localValues, const Ptr< Ordinal > &leadingDim) const
void initializeImpl(const RCP< const TpetraVectorSpace< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVectorSpace, const RCP< const ScalarProdVectorSpaceBase< Scalar > > &domainSpace, const RCP< TpetraMultiVector_t > &tpetraMultiVector)
void initialize(const RCP< const TpetraVectorSpace< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVectorSpace, const RCP< const ScalarProdVectorSpaceBase< Scalar > > &domainSpace, const RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraMultiVector)
Initialize.
RCP< const MultiVectorBase< Scalar > > contigSubViewImpl(const Range1D &colRng) const
RCP< const VectorBase< Scalar > > colImpl(Ordinal j) const
void constInitialize(const RCP< const TpetraVectorSpace< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVectorSpace, const RCP< const ScalarProdVectorSpaceBase< Scalar > > &domainSpace, const RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraMultiVector)
Initialize.
RCP< VectorBase< Scalar > > nonconstColImpl(Ordinal j)
RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getConstTpetraMultiVector() const
Extract the underlying const Tpetra::MultiVector object.