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 //
4 // Thyra: Interfaces and Support for Abstract Numerical Algorithms
5 // Copyright (2004) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Roscoe A. Bartlett (bartlettra@ornl.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
42 #ifndef THYRA_TPETRA_MULTIVECTOR_HPP
43 #define THYRA_TPETRA_MULTIVECTOR_HPP
44 
47 #include "Thyra_TpetraVector.hpp"
48 #include "Teuchos_Assert.hpp"
49 
50 
51 namespace Thyra {
52 
53 
54 // Constructors/initializers/accessors
55 
56 
57 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
59 {}
60 
61 
62 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
64  const RCP<const TpetraVectorSpace<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVectorSpace,
65  const RCP<const ScalarProdVectorSpaceBase<Scalar> > &domainSpace,
66  const RCP<Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraMultiVector
67  )
68 {
69  initializeImpl(tpetraVectorSpace, domainSpace, tpetraMultiVector);
70 }
71 
72 
73 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
75  const RCP<const TpetraVectorSpace<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVectorSpace,
76  const RCP<const ScalarProdVectorSpaceBase<Scalar> > &domainSpace,
77  const RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraMultiVector
78  )
79 {
80  initializeImpl(tpetraVectorSpace, domainSpace, tpetraMultiVector);
81 }
82 
83 
84 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
87 {
88  return tpetraMultiVector_.getNonconstObj();
89 }
90 
91 
92 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
95 {
96  return tpetraMultiVector_;
97 }
98 
99 
100 // Overridden public functions form MultiVectorAdapterBase
101 
102 
103 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
106 {
107  return domainSpace_;
108 }
109 
110 
111 // Overridden protected functions from MultiVectorBase
112 
113 
114 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
115 void
117 {
118  tpetraMultiVector_.getNonconstObj()->putScalar(alpha);
119 }
120 
121 
122 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
124 assignMultiVecImpl(const MultiVectorBase<Scalar>& mv)
125 {
126  auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
127 
128  // If the cast succeeded, call Tpetra directly.
129  // Otherwise, fall back to the RTOp implementation.
130  if (nonnull(tmv)) {
131  tpetraMultiVector_.getNonconstObj()->assign(*tmv);
132  } else {
133  // This version will require/modify the host view of this vector.
134 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE
135  tpetraMultiVector_.getNonconstObj()->template sync<Kokkos::HostSpace>();
136  tpetraMultiVector_.getNonconstObj()->template modify<Kokkos::HostSpace>();
137 #else
138  tpetraMultiVector_.getNonconstObj()->sync_host ();
139  tpetraMultiVector_.getNonconstObj()->modify_host ();
140 #endif
141  MultiVectorDefaultBase<Scalar>::assignMultiVecImpl(mv);
142  }
143 }
144 
145 
146 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
147 void
149 {
150  tpetraMultiVector_.getNonconstObj()->scale(alpha);
151 }
152 
153 
154 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
156  Scalar alpha,
157  const MultiVectorBase<Scalar>& mv
158  )
159 {
160  auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
161 
162  // If the cast succeeded, call Tpetra directly.
163  // Otherwise, fall back to the RTOp implementation.
164  if (nonnull(tmv)) {
166  tpetraMultiVector_.getNonconstObj()->update(alpha, *tmv, ST::one());
167  } else {
168  // This version will require/modify the host view of this vector.
169 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE
170  tpetraMultiVector_.getNonconstObj()->template sync<Kokkos::HostSpace>();
171  tpetraMultiVector_.getNonconstObj()->template modify<Kokkos::HostSpace>();
172 #else
173  tpetraMultiVector_.getNonconstObj()->sync_host ();
174  tpetraMultiVector_.getNonconstObj()->modify_host ();
175 #endif
176  MultiVectorDefaultBase<Scalar>::updateImpl(alpha, mv);
177  }
178 }
179 
180 
181 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
183  const ArrayView<const Scalar>& alpha,
184  const ArrayView<const Ptr<const MultiVectorBase<Scalar> > >& mv,
185  const Scalar& beta
186  )
187 {
188 #ifdef TEUCHOS_DEBUG
189  TEUCHOS_ASSERT_EQUALITY(alpha.size(), mv.size());
190 #endif
191 
192  // Try to cast mv to an array of this type
193  typedef Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> TMV;
194  Teuchos::Array<RCP<const TMV> > tmvs(mv.size());
195  RCP<const TMV> tmv;
196  bool allCastsSuccessful = true;
197  {
198  auto mvIter = mv.begin();
199  auto tmvIter = tmvs.begin();
200  for (; mvIter != mv.end(); ++mvIter, ++tmvIter) {
201  tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromPtr(*mvIter));
202  if (nonnull(tmv)) {
203  *tmvIter = tmv;
204  } else {
205  allCastsSuccessful = false;
206  break;
207  }
208  }
209  }
210 
211  // If casts succeeded, or input arrays are size 0, call Tpetra directly.
212  // Otherwise, fall back to the RTOp implementation.
213  auto len = tmvs.size();
214  if (len == 0) {
215  tpetraMultiVector_.getNonconstObj()->scale(beta);
216  } else if (len == 1 && allCastsSuccessful) {
217  tpetraMultiVector_.getNonconstObj()->update(alpha[0], *tmvs[0], beta);
218  } else if (len == 2 && allCastsSuccessful) {
219  tpetraMultiVector_.getNonconstObj()->update(alpha[0], *tmvs[0], alpha[1], *tmvs[1], beta);
220  } else if (allCastsSuccessful) {
222  auto tmvIter = tmvs.begin();
223  auto alphaIter = alpha.begin();
224 
225  // Check if any entry of tmvs aliases this object's wrapped vector.
226  // If so, replace that entry in the array with a copy.
227  tmv = Teuchos::null;
228  for (; tmvIter != tmvs.end(); ++tmvIter) {
229  if (tmvIter->getRawPtr() == tpetraMultiVector_.getConstObj().getRawPtr()) {
230  if (tmv.is_null()) {
231  tmv = Teuchos::rcp(new TMV(*tpetraMultiVector_.getConstObj(), Teuchos::Copy));
232  }
233  *tmvIter = tmv;
234  }
235  }
236  tmvIter = tmvs.begin();
237 
238  // We add two MVs at a time, so only scale if even num MVs,
239  // and additionally do the first addition if odd num MVs.
240  if ((tmvs.size() % 2) == 0) {
241  tpetraMultiVector_.getNonconstObj()->scale(beta);
242  } else {
243  tpetraMultiVector_.getNonconstObj()->update(*alphaIter, *(*tmvIter), beta);
244  ++tmvIter;
245  ++alphaIter;
246  }
247  for (; tmvIter != tmvs.end(); tmvIter+=2, alphaIter+=2) {
248  tpetraMultiVector_.getNonconstObj()->update(
249  *alphaIter, *(*tmvIter), *(alphaIter+1), *(*(tmvIter+1)), ST::one());
250  }
251  } else {
252  // This version will require/modify the host view of this vector.
253 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE
254  tpetraMultiVector_.getNonconstObj()->template sync<Kokkos::HostSpace>();
255  tpetraMultiVector_.getNonconstObj()->template modify<Kokkos::HostSpace>();
256 #else
257  tpetraMultiVector_.getNonconstObj()->sync_host ();
258  tpetraMultiVector_.getNonconstObj()->modify_host ();
259 #endif
260  MultiVectorDefaultBase<Scalar>::linearCombinationImpl(alpha, mv, beta);
261  }
262 }
263 
264 
265 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
267  const MultiVectorBase<Scalar>& mv,
268  const ArrayView<Scalar>& prods
269  ) const
270 {
271  auto tmv = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(mv));
272 
273  // If the cast succeeded, call Tpetra directly.
274  // Otherwise, fall back to the RTOp implementation.
275  if (nonnull(tmv)) {
276  tpetraMultiVector_.getConstObj()->dot(*tmv, prods);
277  } else {
278  // This version will require/modify the host view of this vector.
279 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE
280  tpetraMultiVector_.getNonconstObj()->template sync<Kokkos::HostSpace>();
281  tpetraMultiVector_.getNonconstObj()->template modify<Kokkos::HostSpace>();
282 #else
283  tpetraMultiVector_.getNonconstObj()->sync_host ();
284  tpetraMultiVector_.getNonconstObj()->modify_host ();
285 #endif
286  MultiVectorDefaultBase<Scalar>::dotsImpl(mv, prods);
287  }
288 }
289 
290 
291 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
293  const ArrayView<typename ScalarTraits<Scalar>::magnitudeType>& norms
294  ) const
295 {
296  tpetraMultiVector_.getConstObj()->norm1(norms);
297 }
298 
299 
300 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
302  const ArrayView<typename ScalarTraits<Scalar>::magnitudeType>& norms
303  ) const
304 {
305  tpetraMultiVector_.getConstObj()->norm2(norms);
306 }
307 
308 
309 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
311  const ArrayView<typename ScalarTraits<Scalar>::magnitudeType>& norms
312  ) const
313 {
314  tpetraMultiVector_.getConstObj()->normInf(norms);
315 }
316 
317 
318 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
321 {
322 #ifdef TEUCHOS_DEBUG
323  TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE(j, 0, this->domain()->dim());
324 #endif
325  return constTpetraVector<Scalar>(
326  tpetraVectorSpace_,
327  tpetraMultiVector_->getVector(j)
328  );
329 }
330 
331 
332 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
335 {
336 #ifdef TEUCHOS_DEBUG
337  TEUCHOS_ASSERT_IN_RANGE_UPPER_EXCLUSIVE(j, 0, this->domain()->dim());
338 #endif
339  return tpetraVector<Scalar>(
340  tpetraVectorSpace_,
341  tpetraMultiVector_.getNonconstObj()->getVectorNonConst(j)
342  );
343 }
344 
345 
346 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
349  const Range1D& col_rng_in
350  ) const
351 {
352 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
353  std::cerr << "\nTpetraMultiVector::subView(Range1D) const called!\n";
354 #endif
355  const Range1D colRng = this->validateColRange(col_rng_in);
356 
358  this->getConstTpetraMultiVector()->subView(colRng);
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 Range1D& col_rng_in
380  )
381 {
382 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
383  std::cerr << "\nTpetraMultiVector::subView(Range1D) called!\n";
384 #endif
385  const Range1D colRng = this->validateColRange(col_rng_in);
386 
388  this->getTpetraMultiVector()->subViewNonConst(colRng);
389 
390  const RCP<const ScalarProdVectorSpaceBase<Scalar> > viewDomainSpace =
391  tpetraVectorSpace<Scalar>(
392  Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal,Node>(
393  tpetraView->getNumVectors(),
394  tpetraView->getMap()->getComm()
395  )
396  );
397 
398  return tpetraMultiVector(
399  tpetraVectorSpace_,
400  viewDomainSpace,
401  tpetraView
402  );
403 }
404 
405 
406 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
409  const ArrayView<const int>& cols_in
410  ) const
411 {
412 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
413  std::cerr << "\nTpetraMultiVector::subView(ArrayView) const called!\n";
414 #endif
415  // Tpetra wants col indices as size_t
416  Array<std::size_t> cols(cols_in.size());
417  for (Array<std::size_t>::size_type i = 0; i < cols.size(); ++i)
418  cols[i] = static_cast<std::size_t>(cols_in[i]);
419 
421  this->getConstTpetraMultiVector()->subView(cols());
422 
423  const RCP<const ScalarProdVectorSpaceBase<Scalar> > viewDomainSpace =
424  tpetraVectorSpace<Scalar>(
425  Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal,Node>(
426  tpetraView->getNumVectors(),
427  tpetraView->getMap()->getComm()
428  )
429  );
430 
431  return constTpetraMultiVector(
432  tpetraVectorSpace_,
433  viewDomainSpace,
434  tpetraView
435  );
436 }
437 
438 
439 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
442  const ArrayView<const int>& cols_in
443  )
444 {
445 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
446  std::cerr << "\nTpetraMultiVector::subView(ArrayView) called!\n";
447 #endif
448  // Tpetra wants col indices as size_t
449  Array<std::size_t> cols(cols_in.size());
450  for (Array<std::size_t>::size_type i = 0; i < cols.size(); ++i)
451  cols[i] = static_cast<std::size_t>(cols_in[i]);
452 
454  this->getTpetraMultiVector()->subViewNonConst(cols());
455 
456  const RCP<const ScalarProdVectorSpaceBase<Scalar> > viewDomainSpace =
457  tpetraVectorSpace<Scalar>(
458  Tpetra::createLocalMapWithNode<LocalOrdinal,GlobalOrdinal,Node>(
459  tpetraView->getNumVectors(),
460  tpetraView->getMap()->getComm()
461  )
462  );
463 
464  return tpetraMultiVector(
465  tpetraVectorSpace_,
466  viewDomainSpace,
467  tpetraView
468  );
469 }
470 
471 
472 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
475  const RTOpPack::RTOpT<Scalar> &primary_op,
476  const ArrayView<const Ptr<const MultiVectorBase<Scalar> > > &multi_vecs,
477  const ArrayView<const Ptr<MultiVectorBase<Scalar> > > &targ_multi_vecs,
478  const ArrayView<const Ptr<RTOpPack::ReductTarget> > &reduct_objs,
479  const Ordinal primary_global_offset
480  ) const
481 {
483 
484  // Sync any non-target Tpetra MVs to host space
485  for (auto itr = multi_vecs.begin(); itr != multi_vecs.end(); ++itr) {
486  Ptr<const TMV> tmv = Teuchos::ptr_dynamic_cast<const TMV>(*itr);
487  if (nonnull(tmv)) {
488 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE
489  Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(
490  tmv->getConstTpetraMultiVector())-> template sync<Kokkos::HostSpace>();
491 #else
492  Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(
493  tmv->getConstTpetraMultiVector())-> sync_host ();
494 #endif
495  }
496  }
497 
498  // Sync any target Tpetra MVs and mark modified
499  for (auto itr = targ_multi_vecs.begin(); itr != targ_multi_vecs.end(); ++itr) {
500  Ptr<TMV> tmv = Teuchos::ptr_dynamic_cast<TMV>(*itr);
501  if (nonnull(tmv)) {
502 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE
503  tmv->getTpetraMultiVector()->template sync<Kokkos::HostSpace>();
504  tmv->getTpetraMultiVector()->template modify<Kokkos::HostSpace>();
505 #else
506  tmv->getTpetraMultiVector()->sync_host ();
507  tmv->getTpetraMultiVector()->modify_host ();
508 #endif
509  }
510  }
511 
512  MultiVectorAdapterBase<Scalar>::mvMultiReductApplyOpImpl(
513  primary_op, multi_vecs, targ_multi_vecs, reduct_objs, primary_global_offset);
514 }
515 
516 
517 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
520  const Range1D &rowRng,
521  const Range1D &colRng,
523  ) const
524 {
525  // Only viewing data, so just sync dual view to host space
526  typedef typename Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> TMV;
527 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE
528  Teuchos::rcp_const_cast<TMV>(
529  tpetraMultiVector_.getConstObj())->template sync<Kokkos::HostSpace>();
530 #else
531  Teuchos::rcp_const_cast<TMV>(
532  tpetraMultiVector_.getConstObj())->sync_host ();
533 #endif
534 
535  SpmdMultiVectorDefaultBase<Scalar>::
536  acquireDetachedMultiVectorViewImpl(rowRng, colRng, sub_mv);
537 }
538 
539 
540 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
543  const Range1D &rowRng,
544  const Range1D &colRng,
546  )
547 {
548  // Sync to host and mark as modified
549 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE
550  tpetraMultiVector_.getNonconstObj()->template sync<Kokkos::HostSpace>();
551  tpetraMultiVector_.getNonconstObj()->template modify<Kokkos::HostSpace>();
552 #else
553  tpetraMultiVector_.getNonconstObj()->sync_host ();
554  tpetraMultiVector_.getNonconstObj()->modify_host ();
555 #endif
556 
557  SpmdMultiVectorDefaultBase<Scalar>::
558  acquireNonconstDetachedMultiVectorViewImpl(rowRng, colRng, sub_mv);
559 }
560 
561 
562 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
566  )
567 {
568  SpmdMultiVectorDefaultBase<Scalar>::
569  commitNonconstDetachedMultiVectorViewImpl(sub_mv);
570 
571  // Sync changes from host view to execution space
572  typedef typename Tpetra::MultiVector<
573  Scalar,LocalOrdinal,GlobalOrdinal,Node>::execution_space execution_space;
574  tpetraMultiVector_.getNonconstObj()->template sync<execution_space>();
575 }
576 
577 
578 /* ToDo: Implement these?
579 
580 
581 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
582 RCP<const MultiVectorBase<Scalar> >
583 TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::nonContigSubViewImpl(
584  const ArrayView<const int> &cols
585  ) const
586 {
587  THYRA_DEBUG_ASSERT_MV_COLS("nonContigSubViewImpl(cols)", cols);
588  const int numCols = cols.size();
589  const ArrayRCP<Scalar> localValuesView = createContiguousCopy(cols);
590  return defaultSpmdMultiVector<Scalar>(
591  spmdRangeSpace_,
592  createSmallScalarProdVectorSpaceBase<Scalar>(spmdRangeSpace_, numCols),
593  localValuesView
594  );
595 }
596 
597 
598 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
599 RCP<MultiVectorBase<Scalar> >
600 TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::nonconstNonContigSubViewImpl(
601  const ArrayView<const int> &cols )
602 {
603  THYRA_DEBUG_ASSERT_MV_COLS("nonContigSubViewImpl(cols)", cols);
604  const int numCols = cols.size();
605  const ArrayRCP<Scalar> localValuesView = createContiguousCopy(cols);
606  const Ordinal localSubDim = spmdRangeSpace_->localSubDim();
607  RCP<CopyBackSpmdMultiVectorEntries<Scalar> > copyBackView =
608  copyBackSpmdMultiVectorEntries<Scalar>(cols, localValuesView.getConst(),
609  localSubDim, localValues_.create_weak(), leadingDim_);
610  return Teuchos::rcpWithEmbeddedObjPreDestroy(
611  new TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>(
612  spmdRangeSpace_,
613  createSmallScalarProdVectorSpaceBase<Scalar>(spmdRangeSpace_, numCols),
614  localValuesView),
615  copyBackView
616  );
617 }
618 
619 */
620 
621 
622 // Overridden protected members from SpmdMultiVectorBase
623 
624 
625 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
628 {
629  return tpetraVectorSpace_;
630 }
631 
632 
633 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
635  const Ptr<ArrayRCP<Scalar> > &localValues, const Ptr<Ordinal> &leadingDim
636  )
637 {
638  *localValues = tpetraMultiVector_.getNonconstObj()->get1dViewNonConst();
639  *leadingDim = tpetraMultiVector_->getStride();
640 }
641 
642 
643 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
645  const Ptr<ArrayRCP<const Scalar> > &localValues, const Ptr<Ordinal> &leadingDim
646  ) const
647 {
648  *localValues = tpetraMultiVector_->get1dView();
649  *leadingDim = tpetraMultiVector_->getStride();
650 }
651 
652 
653 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
655  const EOpTransp M_trans,
656  const MultiVectorBase<Scalar> &X,
657  const Ptr<MultiVectorBase<Scalar> > &Y,
658  const Scalar alpha,
659  const Scalar beta
660  ) const
661 {
662  // Try to extract Tpetra objects from X and Y
663  typedef Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> TMV;
664  Teuchos::RCP<const TMV> X_tpetra = this->getConstTpetraMultiVector(Teuchos::rcpFromRef(X));
665  Teuchos::RCP<TMV> Y_tpetra = this->getTpetraMultiVector(Teuchos::rcpFromPtr(Y));
666 
667  // If the cast succeeded, call Tpetra directly.
668  // Otherwise, fall back to the default implementation.
669  if (nonnull(X_tpetra) && nonnull(Y_tpetra)) {
670  // Sync everything to the execution space
671  typedef typename TMV::execution_space execution_space;
672  Teuchos::rcp_const_cast<TMV>(X_tpetra)->template sync<execution_space>();
673  Y_tpetra->template sync<execution_space>();
674  Teuchos::rcp_const_cast<TMV>(
675  tpetraMultiVector_.getConstObj())->template sync<execution_space>();
676 
678  TEUCHOS_TEST_FOR_EXCEPTION(ST::isComplex && (M_trans == CONJ),
679  std::logic_error,
680  "Error, conjugation without transposition is not allowed for complex scalar types!");
681 
683  switch (M_trans) {
684  case NOTRANS:
685  trans = Teuchos::NO_TRANS;
686  break;
687  case CONJ:
688  trans = Teuchos::NO_TRANS;
689  break;
690  case TRANS:
691  trans = Teuchos::TRANS;
692  break;
693  case CONJTRANS:
694  trans = Teuchos::CONJ_TRANS;
695  break;
696  }
697 
698  Y_tpetra->template modify<execution_space>();
699  Y_tpetra->multiply(trans, Teuchos::NO_TRANS, alpha, *tpetraMultiVector_.getConstObj(), *X_tpetra, beta);
700  } else {
701 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE
702  Teuchos::rcp_const_cast<TMV>(
703  tpetraMultiVector_.getConstObj())->template sync<Kokkos::HostSpace>();
704 #else
705  Teuchos::rcp_const_cast<TMV>(
706  tpetraMultiVector_.getConstObj())->sync_host ();
707 #endif
708  SpmdMultiVectorDefaultBase<Scalar>::euclideanApply(M_trans, X, Y, alpha, beta);
709  }
710 
711 }
712 
713 // private
714 
715 
716 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
717 template<class TpetraMultiVector_t>
719  const RCP<const TpetraVectorSpace<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tpetraVectorSpace,
720  const RCP<const ScalarProdVectorSpaceBase<Scalar> > &domainSpace,
721  const RCP<TpetraMultiVector_t> &tpetraMultiVector
722  )
723 {
724 #ifdef THYRA_DEBUG
725  TEUCHOS_ASSERT(nonnull(tpetraVectorSpace));
726  TEUCHOS_ASSERT(nonnull(domainSpace));
727  TEUCHOS_ASSERT(nonnull(tpetraMultiVector));
728  // ToDo: Check to make sure that tpetraMultiVector is compatible with
729  // tpetraVectorSpace.
730 #endif
731  tpetraVectorSpace_ = tpetraVectorSpace;
732  domainSpace_ = domainSpace;
733  tpetraMultiVector_.initialize(tpetraMultiVector);
734  this->updateSpmdSpace();
735 }
736 
737 
738 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
741 getTpetraMultiVector(const RCP<MultiVectorBase<Scalar> >& mv) const
742 {
743  using Teuchos::rcp_dynamic_cast;
746 
747  RCP<TMV> tmv = rcp_dynamic_cast<TMV>(mv);
748  if (nonnull(tmv)) {
749  return tmv->getTpetraMultiVector();
750  }
751 
752  RCP<TV> tv = rcp_dynamic_cast<TV>(mv);
753  if (nonnull(tv)) {
754  return tv->getTpetraVector();
755  }
756 
757  return Teuchos::null;
758 }
759 
760 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
763 getConstTpetraMultiVector(const RCP<const MultiVectorBase<Scalar> >& mv) const
764 {
765  using Teuchos::rcp_dynamic_cast;
768 
769  RCP<const TMV> tmv = rcp_dynamic_cast<const TMV>(mv);
770  if (nonnull(tmv)) {
771  return tmv->getConstTpetraMultiVector();
772  }
773 
774  RCP<const TV> tv = rcp_dynamic_cast<const TV>(mv);
775  if (nonnull(tv)) {
776  return tv->getConstTpetraVector();
777  }
778 
779  return Teuchos::null;
780 }
781 
782 
783 } // end namespace Thyra
784 
785 
786 #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)