Thyra Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Thyra_TpetraMultiVector_def.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Thyra: Interfaces and Support for Abstract Numerical Algorithms
4 //
5 // Copyright 2004 NTESS and the Thyra contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef THYRA_TPETRA_MULTIVECTOR_HPP
11 #define THYRA_TPETRA_MULTIVECTOR_HPP
12 
15 #include "Thyra_TpetraVector.hpp"
16 #include "Teuchos_Assert.hpp"
17 #include "Kokkos_Core.hpp"
18 
19 
20 namespace Thyra {
21 
22 
23 // Constructors/initializers/accessors
24 
25 
26 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
28 {}
29 
30 
31 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
33  const RCP<const TpetraVectorSpace<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVectorSpace,
34  const RCP<const ScalarProdVectorSpaceBase<Scalar> > &domainSpace,
35  const RCP<Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraMultiVector
36  )
37 {
38  initializeImpl(tpetraVectorSpace, domainSpace, tpetraMultiVector);
39 }
40 
41 
42 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
44  const RCP<const TpetraVectorSpace<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVectorSpace,
45  const RCP<const ScalarProdVectorSpaceBase<Scalar> > &domainSpace,
46  const RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraMultiVector
47  )
48 {
49  initializeImpl(tpetraVectorSpace, domainSpace, tpetraMultiVector);
50 }
51 
52 
53 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
56 {
57  return tpetraMultiVector_.getNonconstObj();
58 }
59 
60 
61 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
64 {
65  return tpetraMultiVector_;
66 }
67 
68 
69 // Overridden public functions form MultiVectorAdapterBase
70 
71 
72 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
75 {
76  return domainSpace_;
77 }
78 
79 
80 // Overridden protected functions from MultiVectorBase
81 
82 
83 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
84 void
86 {
87  tpetraMultiVector_.getNonconstObj()->putScalar(alpha);
88 }
89 
90 
91 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
93 assignMultiVecImpl(const MultiVectorBase<Scalar>& mv)
94 {
95  auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
96 
97  // If the cast succeeded, call Tpetra directly.
98  // Otherwise, fall back to the RTOp implementation.
99  if (nonnull(tmv)) {
100  tpetraMultiVector_.getNonconstObj()->assign(*tmv);
101  } else {
102  MultiVectorDefaultBase<Scalar>::assignMultiVecImpl(mv);
103  }
104 }
105 
106 
107 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
108 void
110 {
111  tpetraMultiVector_.getNonconstObj()->scale(alpha);
112 }
113 
114 
115 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
117  Scalar alpha,
118  const MultiVectorBase<Scalar>& mv
119  )
120 {
121  auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
122 
123  // If the cast succeeded, call Tpetra directly.
124  // Otherwise, fall back to the RTOp implementation.
125  if (nonnull(tmv)) {
127  tpetraMultiVector_.getNonconstObj()->update(alpha, *tmv, ST::one());
128  } else {
129  MultiVectorDefaultBase<Scalar>::updateImpl(alpha, mv);
130  }
131 }
132 
133 
134 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
136  const ArrayView<const Scalar>& alpha,
137  const ArrayView<const Ptr<const MultiVectorBase<Scalar> > >& mv,
138  const Scalar& beta
139  )
140 {
141 #ifdef TEUCHOS_DEBUG
142  TEUCHOS_ASSERT_EQUALITY(alpha.size(), mv.size());
143 #endif
144 
145  // Try to cast mv to an array of this type
146  typedef Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> TMV;
147  Teuchos::Array<RCP<const TMV> > tmvs(mv.size());
148  RCP<const TMV> tmv;
149  bool allCastsSuccessful = true;
150  {
151  auto mvIter = mv.begin();
152  auto tmvIter = tmvs.begin();
153  for (; mvIter != mv.end(); ++mvIter, ++tmvIter) {
154  tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromPtr(*mvIter));
155  if (nonnull(tmv)) {
156  *tmvIter = tmv;
157  } else {
158  allCastsSuccessful = false;
159  break;
160  }
161  }
162  }
163 
164  // If casts succeeded, or input arrays are size 0, call Tpetra directly.
165  // Otherwise, fall back to the RTOp implementation.
166  auto len = tmvs.size();
167  if (len == 0) {
168  tpetraMultiVector_.getNonconstObj()->scale(beta);
169  } else if (len == 1 && allCastsSuccessful) {
170  tpetraMultiVector_.getNonconstObj()->update(alpha[0], *tmvs[0], beta);
171  } else if (len == 2 && allCastsSuccessful) {
172  tpetraMultiVector_.getNonconstObj()->update(alpha[0], *tmvs[0], alpha[1], *tmvs[1], beta);
173  } else if (allCastsSuccessful) {
175  auto tmvIter = tmvs.begin();
176  auto alphaIter = alpha.begin();
177 
178  // Check if any entry of tmvs aliases this object's wrapped vector.
179  // If so, replace that entry in the array with a copy.
180  tmv = Teuchos::null;
181  for (; tmvIter != tmvs.end(); ++tmvIter) {
182  if (tmvIter->getRawPtr() == tpetraMultiVector_.getConstObj().getRawPtr()) {
183  if (tmv.is_null()) {
184  tmv = Teuchos::rcp(new TMV(*tpetraMultiVector_.getConstObj(), Teuchos::Copy));
185  }
186  *tmvIter = tmv;
187  }
188  }
189  tmvIter = tmvs.begin();
190 
191  // We add two MVs at a time, so only scale if even num MVs,
192  // and additionally do the first addition if odd num MVs.
193  if ((tmvs.size() % 2) == 0) {
194  tpetraMultiVector_.getNonconstObj()->scale(beta);
195  } else {
196  tpetraMultiVector_.getNonconstObj()->update(*alphaIter, *(*tmvIter), beta);
197  ++tmvIter;
198  ++alphaIter;
199  }
200  for (; tmvIter != tmvs.end(); tmvIter+=2, alphaIter+=2) {
201  tpetraMultiVector_.getNonconstObj()->update(
202  *alphaIter, *(*tmvIter), *(alphaIter+1), *(*(tmvIter+1)), ST::one());
203  }
204  } else {
205  MultiVectorDefaultBase<Scalar>::linearCombinationImpl(alpha, mv, beta);
206  }
207 }
208 
209 
210 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
212  const MultiVectorBase<Scalar>& mv,
213  const ArrayView<Scalar>& prods
214  ) const
215 {
216  auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
217 
218  // If the cast succeeded, call Tpetra directly.
219  // Otherwise, fall back to the RTOp implementation.
220  if (nonnull(tmv)) {
221  tpetraMultiVector_.getConstObj()->dot(*tmv, prods);
222  } else {
223  MultiVectorDefaultBase<Scalar>::dotsImpl(mv, prods);
224  }
225 }
226 
227 
228 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
230  const ArrayView<typename ScalarTraits<Scalar>::magnitudeType>& norms
231  ) const
232 {
233  tpetraMultiVector_.getConstObj()->norm1(norms);
234 }
235 
236 
237 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
239  const ArrayView<typename ScalarTraits<Scalar>::magnitudeType>& norms
240  ) const
241 {
242  tpetraMultiVector_.getConstObj()->norm2(norms);
243 }
244 
245 
246 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
248  const ArrayView<typename ScalarTraits<Scalar>::magnitudeType>& norms
249  ) const
250 {
251  tpetraMultiVector_.getConstObj()->normInf(norms);
252 }
253 
254 
255 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
258 {
259 #ifdef TEUCHOS_DEBUG
260  TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE(j, 0, this->domain()->dim());
261 #endif
262  return constTpetraVector<Scalar>(
263  tpetraVectorSpace_,
264  tpetraMultiVector_->getVector(j)
265  );
266 }
267 
268 
269 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
272 {
273 #ifdef TEUCHOS_DEBUG
274  TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE(j, 0, this->domain()->dim());
275 #endif
276  return tpetraVector<Scalar>(
277  tpetraVectorSpace_,
278  tpetraMultiVector_.getNonconstObj()->getVectorNonConst(j)
279  );
280 }
281 
282 
283 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
286  const Range1D& col_rng_in
287  ) const
288 {
289 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
290  std::cerr << "\nTpetraMultiVector::subView(Range1D) const called!\n";
291 #endif
292  const Range1D colRng = this->validateColRange(col_rng_in);
293 
295  this->getConstTpetraMultiVector()->subView(colRng);
296 
297  const RCP<const ScalarProdVectorSpaceBase<Scalar> > viewDomainSpace =
298  tpetraVectorSpace<Scalar>(
299  Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal,Node>(
300  tpetraView->getNumVectors(),
301  tpetraView->getMap()->getComm()
302  )
303  );
304 
305  return constTpetraMultiVector(
306  tpetraVectorSpace_,
307  viewDomainSpace,
308  tpetraView
309  );
310 }
311 
312 
313 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
316  const Range1D& col_rng_in
317  )
318 {
319 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
320  std::cerr << "\nTpetraMultiVector::subView(Range1D) called!\n";
321 #endif
322  const Range1D colRng = this->validateColRange(col_rng_in);
323 
325  this->getTpetraMultiVector()->subViewNonConst(colRng);
326 
327  const RCP<const ScalarProdVectorSpaceBase<Scalar> > viewDomainSpace =
328  tpetraVectorSpace<Scalar>(
329  Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal,Node>(
330  tpetraView->getNumVectors(),
331  tpetraView->getMap()->getComm()
332  )
333  );
334 
335  return tpetraMultiVector(
336  tpetraVectorSpace_,
337  viewDomainSpace,
338  tpetraView
339  );
340 }
341 
342 
343 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
346  const ArrayView<const int>& cols_in
347  ) const
348 {
349 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
350  std::cerr << "\nTpetraMultiVector::subView(ArrayView) const called!\n";
351 #endif
352  // Tpetra wants col indices as size_t
353  Array<std::size_t> cols(cols_in.size());
354  for (Array<std::size_t>::size_type i = 0; i < cols.size(); ++i)
355  cols[i] = static_cast<std::size_t>(cols_in[i]);
356 
358  this->getConstTpetraMultiVector()->subView(cols());
359 
360  const RCP<const ScalarProdVectorSpaceBase<Scalar> > viewDomainSpace =
361  tpetraVectorSpace<Scalar>(
362  Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal,Node>(
363  tpetraView->getNumVectors(),
364  tpetraView->getMap()->getComm()
365  )
366  );
367 
368  return constTpetraMultiVector(
369  tpetraVectorSpace_,
370  viewDomainSpace,
371  tpetraView
372  );
373 }
374 
375 
376 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
379  const ArrayView<const int>& cols_in
380  )
381 {
382 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
383  std::cerr << "\nTpetraMultiVector::subView(ArrayView) called!\n";
384 #endif
385  // Tpetra wants col indices as size_t
386  Array<std::size_t> cols(cols_in.size());
387  for (Array<std::size_t>::size_type i = 0; i < cols.size(); ++i)
388  cols[i] = static_cast<std::size_t>(cols_in[i]);
389 
391  this->getTpetraMultiVector()->subViewNonConst(cols());
392 
393  const RCP<const ScalarProdVectorSpaceBase<Scalar> > viewDomainSpace =
394  tpetraVectorSpace<Scalar>(
395  Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal,Node>(
396  tpetraView->getNumVectors(),
397  tpetraView->getMap()->getComm()
398  )
399  );
400 
401  return tpetraMultiVector(
402  tpetraVectorSpace_,
403  viewDomainSpace,
404  tpetraView
405  );
406 }
407 
408 
409 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
412  const RTOpPack::RTOpT<Scalar> &primary_op,
413  const ArrayView<const Ptr<const MultiVectorBase<Scalar> > > &multi_vecs,
414  const ArrayView<const Ptr<MultiVectorBase<Scalar> > > &targ_multi_vecs,
415  const ArrayView<const Ptr<RTOpPack::ReductTarget> > &reduct_objs,
416  const Ordinal primary_global_offset
417  ) const
418 {
419 
420  MultiVectorAdapterBase<Scalar>::mvMultiReductApplyOpImpl(
421  primary_op, multi_vecs, targ_multi_vecs, reduct_objs, primary_global_offset);
422 }
423 
424 
425 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
428  const Range1D &rowRng,
429  const Range1D &colRng,
431  ) const
432 {
433  SpmdMultiVectorDefaultBase<Scalar>::
434  acquireDetachedMultiVectorViewImpl(rowRng, colRng, sub_mv);
435 }
436 
437 
438 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
441  const Range1D &rowRng,
442  const Range1D &colRng,
444  )
445 {
446  SpmdMultiVectorDefaultBase<Scalar>::
447  acquireNonconstDetachedMultiVectorViewImpl(rowRng, colRng, sub_mv);
448 }
449 
450 
451 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
455  )
456 {
457  SpmdMultiVectorDefaultBase<Scalar>::
458  commitNonconstDetachedMultiVectorViewImpl(sub_mv);
459 
460 }
461 
462 
463 /* ToDo: Implement these?
464 
465 
466 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
467 RCP<const MultiVectorBase<Scalar> >
468 TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::nonContigSubViewImpl(
469  const ArrayView<const int> &cols
470  ) const
471 {
472  THYRA_DEBUG_ASSERT_MV_COLS("nonContigSubViewImpl(cols)", cols);
473  const int numCols = cols.size();
474  const ArrayRCP<Scalar> localValuesView = createContiguousCopy(cols);
475  return defaultSpmdMultiVector<Scalar>(
476  spmdRangeSpace_,
477  createSmallScalarProdVectorSpaceBase<Scalar>(spmdRangeSpace_, numCols),
478  localValuesView
479  );
480 }
481 
482 
483 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
484 RCP<MultiVectorBase<Scalar> >
485 TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::nonconstNonContigSubViewImpl(
486  const ArrayView<const int> &cols )
487 {
488  THYRA_DEBUG_ASSERT_MV_COLS("nonContigSubViewImpl(cols)", cols);
489  const int numCols = cols.size();
490  const ArrayRCP<Scalar> localValuesView = createContiguousCopy(cols);
491  const Ordinal localSubDim = spmdRangeSpace_->localSubDim();
492  RCP<CopyBackSpmdMultiVectorEntries<Scalar> > copyBackView =
493  copyBackSpmdMultiVectorEntries<Scalar>(cols, localValuesView.getConst(),
494  localSubDim, localValues_.create_weak(), leadingDim_);
495  return Teuchos::rcpWithEmbeddedObjPreDestroy(
496  new TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>(
497  spmdRangeSpace_,
498  createSmallScalarProdVectorSpaceBase<Scalar>(spmdRangeSpace_, numCols),
499  localValuesView),
500  copyBackView
501  );
502 }
503 
504 */
505 
506 
507 // Overridden protected members from SpmdMultiVectorBase
508 
509 
510 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
513 {
514  return tpetraVectorSpace_;
515 }
516 
517 
518 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
520  const Ptr<ArrayRCP<Scalar> > &localValues, const Ptr<Ordinal> &leadingDim
521  )
522 {
523  *localValues = tpetraMultiVector_.getNonconstObj()->get1dViewNonConst();
524  *leadingDim = tpetraMultiVector_->getStride();
525 }
526 
527 
528 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
530  const Ptr<ArrayRCP<const Scalar> > &localValues, const Ptr<Ordinal> &leadingDim
531  ) const
532 {
533  *localValues = tpetraMultiVector_->get1dView();
534  *leadingDim = tpetraMultiVector_->getStride();
535 }
536 
537 
538 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
540  const EOpTransp M_trans,
541  const MultiVectorBase<Scalar> &X,
542  const Ptr<MultiVectorBase<Scalar> > &Y,
543  const Scalar alpha,
544  const Scalar beta
545  ) const
546 {
547  // Try to extract Tpetra objects from X and Y
548  typedef Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> TMV;
549  Teuchos::RCP<const TMV> X_tpetra = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(X));
550  Teuchos::RCP<TMV> Y_tpetra = this->getTpetraMultiVector(Teuchos::rcpFromPtr(Y));
551 
552  // If the cast succeeded, call Tpetra directly.
553  // Otherwise, fall back to the default implementation.
554  if (nonnull(X_tpetra) && nonnull(Y_tpetra)) {
556  TEUCHOS_TEST_FOR_EXCEPTION(ST::isComplex && (M_trans == CONJ),
557  std::logic_error,
558  "Error, conjugation without transposition is not allowed for complex scalar types!");
559 
561  switch (M_trans) {
562  case NOTRANS:
563  trans = Teuchos::NO_TRANS;
564  break;
565  case CONJ:
566  trans = Teuchos::NO_TRANS;
567  break;
568  case TRANS:
569  trans = Teuchos::TRANS;
570  break;
571  case CONJTRANS:
572  trans = Teuchos::CONJ_TRANS;
573  break;
574  }
575 
576  Y_tpetra->multiply(trans, Teuchos::NO_TRANS, alpha, *tpetraMultiVector_.getConstObj(), *X_tpetra, beta);
577  Kokkos::fence();
578  } else {
579  SpmdMultiVectorDefaultBase<Scalar>::euclideanApply(M_trans, X, Y, alpha, beta);
580  }
581 
582 }
583 
584 // private
585 
586 
587 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
588 template<class TpetraMultiVector_t>
590  const RCP<const TpetraVectorSpace<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVectorSpace,
591  const RCP<const ScalarProdVectorSpaceBase<Scalar> > &domainSpace,
592  const RCP<TpetraMultiVector_t> &tpetraMultiVector
593  )
594 {
595 #ifdef THYRA_DEBUG
596  TEUCHOS_ASSERT(nonnull(tpetraVectorSpace));
597  TEUCHOS_ASSERT(nonnull(domainSpace));
598  TEUCHOS_ASSERT(nonnull(tpetraMultiVector));
599  // ToDo: Check to make sure that tpetraMultiVector is compatible with
600  // tpetraVectorSpace.
601 #endif
602  tpetraVectorSpace_ = tpetraVectorSpace;
603  domainSpace_ = domainSpace;
604  tpetraMultiVector_.initialize(tpetraMultiVector);
605  this->updateSpmdSpace();
606 }
607 
608 
609 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
612 getTpetraMultiVector(const RCP<MultiVectorBase<Scalar> >& mv) const
613 {
614  using Teuchos::rcp_dynamic_cast;
617 
618  RCP<TMV> tmv = rcp_dynamic_cast<TMV>(mv);
619  if (nonnull(tmv)) {
620  return tmv->getTpetraMultiVector();
621  }
622 
623  RCP<TV> tv = rcp_dynamic_cast<TV>(mv);
624  if (nonnull(tv)) {
625  return tv->getTpetraVector();
626  }
627 
628  return Teuchos::null;
629 }
630 
631 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
634 getConstTpetraMultiVector(const RCP<const MultiVectorBase<Scalar> >& mv) const
635 {
636  using Teuchos::rcp_dynamic_cast;
639 
640  RCP<const TMV> tmv = rcp_dynamic_cast<const TMV>(mv);
641  if (nonnull(tmv)) {
642  return tmv->getConstTpetraMultiVector();
643  }
644 
645  RCP<const TV> tv = rcp_dynamic_cast<const TV>(mv);
646  if (nonnull(tv)) {
647  return tv->getConstTpetraVector();
648  }
649 
650  return Teuchos::null;
651 }
652 
653 
654 } // end namespace Thyra
655 
656 
657 #endif // THYRA_TPETRA_MULTIVECTOR_HPP
virtual void updateImpl(Scalar alpha, const MultiVectorBase< Scalar > &mv)
Concrete implementation of Thyra::MultiVector in terms of Tpetra::MultiVector.
TpetraMultiVector()
Construct to uninitialized.
RCP< const VectorBase< Scalar > > colImpl(Ordinal j) const
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.
virtual void assignImpl(Scalar alpha)
RCP< MultiVectorBase< Scalar > > nonconstContigSubViewImpl(const Range1D &colRng)
iterator begin() const
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual void euclideanApply(const EOpTransp M_trans, const MultiVectorBase< Scalar > &X, const Ptr< MultiVectorBase< Scalar > > &Y, const Scalar alpha, const Scalar beta) const
size_type size() const
void initializeImpl(const RCP< const TpetraVectorSpace< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVectorSpace, const RCP< const ScalarProdVectorSpaceBase< Scalar > > &domainSpace, const RCP< TpetraMultiVector_t > &tpetraMultiVector)
RCP< const MultiVectorBase< Scalar > > nonContigSubViewImpl(const ArrayView< const int > &cols_in) const
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.
#define TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE(index, lower_inclusive, upper_exclusive)
virtual void scaleImpl(Scalar alpha)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
RCP< const MultiVectorBase< Scalar > > contigSubViewImpl(const Range1D &colRng) const
virtual void mvMultiReductApplyOpImpl(const RTOpPack::RTOpT< Scalar > &primary_op, const ArrayView< const Ptr< const MultiVectorBase< Scalar > > > &multi_vecs, const ArrayView< const Ptr< MultiVectorBase< Scalar > > > &targ_multi_vecs, const ArrayView< const Ptr< RTOpPack::ReductTarget > > &reduct_objs, const Ordinal primary_global_offset) const
RCP< MultiVectorBase< Scalar > > nonconstNonContigSubViewImpl(const ArrayView< const int > &cols_in)
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.
virtual void norms2Impl(const ArrayView< typename ScalarTraits< Scalar >::magnitudeType > &norms) const
virtual void linearCombinationImpl(const ArrayView< const Scalar > &alpha, const ArrayView< const Ptr< const MultiVectorBase< Scalar > > > &mv, const Scalar &beta)
Concrete Thyra::SpmdVectorBase using Tpetra::Vector.
virtual void dotsImpl(const MultiVectorBase< Scalar > &mv, const ArrayView< Scalar > &prods) const
RCP< VectorBase< Scalar > > nonconstColImpl(Ordinal j)
bool nonnull(const boost::shared_ptr< T > &p)
void getLocalMultiVectorDataImpl(const Ptr< ArrayRCP< const Scalar > > &localValues, const Ptr< Ordinal > &leadingDim) const
size_type size() const
virtual void assignMultiVecImpl(const MultiVectorBase< Scalar > &mv)
virtual void norms1Impl(const ArrayView< typename ScalarTraits< Scalar >::magnitudeType > &norms) const
RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getConstTpetraMultiVector() const
Extract the underlying const Tpetra::MultiVector object.
RCP< const SpmdVectorSpaceBase< Scalar > > spmdSpaceImpl() const
#define TEUCHOS_ASSERT(assertion_test)
#define TEUCHOS_ASSERT_EQUALITY(val1, val2)
virtual void normsInfImpl(const ArrayView< typename ScalarTraits< Scalar >::magnitudeType > &norms) const
Teuchos_Ordinal Ordinal
iterator begin()
RCP< const ScalarProdVectorSpaceBase< Scalar > > domainScalarProdVecSpc() const
void acquireNonconstDetachedMultiVectorViewImpl(const Range1D &rowRng, const Range1D &colRng, RTOpPack::SubMultiVectorView< Scalar > *sub_mv)
void acquireDetachedMultiVectorViewImpl(const Range1D &rowRng, const Range1D &colRng, RTOpPack::ConstSubMultiVectorView< Scalar > *sub_mv) const
void commitNonconstDetachedMultiVectorViewImpl(RTOpPack::SubMultiVectorView< Scalar > *sub_mv)