46 #ifndef BELOS_XPETRA_ADAPTER_MULTIVECTOR_HPP 47 #define BELOS_XPETRA_ADAPTER_MULTIVECTOR_HPP 49 #include <Xpetra_ConfigDefs.hpp> 51 #include <Xpetra_MultiVector.hpp> 53 #ifdef HAVE_XPETRA_EPETRA 54 #include <Xpetra_EpetraMultiVector.hpp> 57 #ifdef HAVE_XPETRA_TPETRA 58 #include <Xpetra_TpetraMultiVector.hpp> 61 #include <BelosConfigDefs.hpp> 62 #include <BelosTypes.hpp> 63 #include <BelosMultiVecTraits.hpp> 64 #include <BelosOperatorTraits.hpp> 66 #ifdef HAVE_XPETRA_EPETRA 67 #include <BelosEpetraAdapter.hpp> 70 #ifdef HAVE_XPETRA_TPETRA 71 #include <BelosTpetraAdapter.hpp> 89 template<
class Scalar,
class LO,
class GO,
class Node>
90 class MultiVecTraits<Scalar, Xpetra::MultiVector<Scalar,LO,GO,Node> >
94 #ifdef HAVE_XPETRA_TPETRA 95 typedef Xpetra::TpetraMultiVector<Scalar,LO,GO,Node> TpetraMultiVector;
96 typedef MultiVecTraits <Scalar,Tpetra::MultiVector<Scalar,LO,GO,Node> > MultiVecTraitsTpetra;
101 #ifdef HAVE_BELOS_XPETRA_TIMERS 102 static RCP<Teuchos::Time> mvTimesMatAddMvTimer_, mvTransMvTimer_;
105 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
Clone(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
const int numvecs )
108 #ifdef HAVE_XPETRA_TPETRA 109 if (mv.getMap()->lib() == Xpetra::UseTpetra)
110 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::Clone(toTpetra(mv), numvecs)));
113 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
117 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
CloneCopy(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
120 #ifdef HAVE_XPETRA_TPETRA 121 if (mv.getMap()->lib() == Xpetra::UseTpetra)
122 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv))));
125 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
129 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
CloneCopy(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
const std::vector<int>& index )
132 #ifdef HAVE_XPETRA_TPETRA 133 if (mv.getMap()->lib() == Xpetra::UseTpetra)
134 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv), index)));
137 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
141 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
142 CloneCopy (
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
143 const Teuchos::Range1D& index)
146 #ifdef HAVE_XPETRA_TPETRA 147 if (mv.getMap()->lib() == Xpetra::UseTpetra)
148 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv), index)));
151 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
155 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
CloneViewNonConst( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
const std::vector<int>& index )
158 #ifdef HAVE_XPETRA_TPETRA 159 if (mv.getMap()->lib() == Xpetra::UseTpetra)
160 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::CloneViewNonConst(toTpetra(mv), index)));
163 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
167 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
169 const Teuchos::Range1D& index)
172 #ifdef HAVE_XPETRA_TPETRA 173 if (mv.getMap()->lib() == Xpetra::UseTpetra)
174 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::CloneViewNonConst(toTpetra(mv), index)));
177 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
181 static RCP<const Xpetra::MultiVector<Scalar,LO,GO,Node> >
CloneView(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
const std::vector<int>& index )
184 #ifdef HAVE_XPETRA_TPETRA 185 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
187 RCP<const Tpetra::MultiVector<Scalar,LO,GO,Node> > r = MultiVecTraitsTpetra::CloneView(toTpetra(mv), index);
188 return rcp(
new TpetraMultiVector(Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LO,GO,Node> >(r)));
192 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
196 static RCP<const Xpetra::MultiVector<Scalar,LO,GO,Node> >
197 CloneView (
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
198 const Teuchos::Range1D& index)
201 #ifdef HAVE_XPETRA_TPETRA 202 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
204 RCP<const Tpetra::MultiVector<Scalar,LO,GO,Node> > r = MultiVecTraitsTpetra::CloneView(toTpetra(mv), index);
205 return rcp(
new TpetraMultiVector(Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LO,GO,Node> >(r)));
209 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
216 #ifdef HAVE_XPETRA_TPETRA 217 if (mv.getMap()->lib() == Xpetra::UseTpetra)
218 return MultiVecTraitsTpetra::GetGlobalLength(toTpetra(mv));
221 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
225 static int GetNumberVecs(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
228 #ifdef HAVE_XPETRA_TPETRA 229 if (mv.getMap()->lib() == Xpetra::UseTpetra)
230 return MultiVecTraitsTpetra::GetNumberVecs(toTpetra(mv));
233 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
240 #ifdef HAVE_XPETRA_TPETRA 241 if (mv.getMap()->lib() == Xpetra::UseTpetra)
242 return MultiVecTraitsTpetra::HasConstantStride(toTpetra(mv));
245 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
249 static void MvTimesMatAddMv( Scalar alpha,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
250 const Teuchos::SerialDenseMatrix<int,Scalar>& B,
251 Scalar beta, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
254 #ifdef HAVE_BELOS_XPETRA_TIMERS 255 Teuchos::TimeMonitor lcltimer(*mvTimesMatAddMvTimer_);
258 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
262 #ifdef HAVE_XPETRA_TPETRA 263 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
264 MultiVecTraitsTpetra::MvTimesMatAddMv(alpha, toTpetra(A), B, beta, toTpetra(mv));
269 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
273 static void MvAddMv( Scalar alpha,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, Scalar beta,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
276 #ifdef HAVE_XPETRA_TPETRA 277 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
278 MultiVecTraitsTpetra::MvAddMv(alpha, toTpetra(A), beta, toTpetra(B), toTpetra(mv));
284 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
288 static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, Scalar alpha )
291 #ifdef HAVE_XPETRA_TPETRA 292 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
293 MultiVecTraitsTpetra::MvScale(toTpetra(mv), alpha);
298 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
302 static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
const std::vector<Scalar>& alphas )
305 #ifdef HAVE_XPETRA_TPETRA 306 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
307 MultiVecTraitsTpetra::MvScale(toTpetra(mv), alphas);
312 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
316 static void MvTransMv( Scalar alpha,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, Teuchos::SerialDenseMatrix<int,Scalar>& C)
319 #ifdef HAVE_BELOS_XPETRA_TIMERS 320 Teuchos::TimeMonitor lcltimer(*mvTransMvTimer_);
323 #ifdef HAVE_XPETRA_TPETRA 324 if (A.getMap()->lib() == Xpetra::UseTpetra) {
325 MultiVecTraitsTpetra::MvTransMv(alpha, toTpetra(A), toTpetra(B), C);
330 XPETRA_FACTORY_ERROR_IF_EPETRA(A.getMap()->lib());
334 static void MvDot(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, std::vector<Scalar> &dots)
337 #ifdef HAVE_XPETRA_TPETRA 338 if (A.getMap()->lib() == Xpetra::UseTpetra) {
339 MultiVecTraitsTpetra::MvDot(toTpetra(A), toTpetra(B), dots);
344 XPETRA_FACTORY_ERROR_IF_EPETRA(A.getMap()->lib());
348 static void MvNorm(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, std::vector<
typename Teuchos::ScalarTraits<Scalar>::magnitudeType> &normvec, NormType type=TwoNorm)
351 #ifdef HAVE_XPETRA_TPETRA 352 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
353 MultiVecTraitsTpetra::MvNorm(toTpetra(mv), normvec, type);
358 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
362 static void SetBlock(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
const std::vector<int>& index, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
365 #ifdef HAVE_XPETRA_TPETRA 366 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
367 MultiVecTraitsTpetra::SetBlock(toTpetra(A), index, toTpetra(mv));
372 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
377 SetBlock (
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
378 const Teuchos::Range1D& index,
379 Xpetra::MultiVector<Scalar,LO,GO,Node>& mv)
382 #ifdef HAVE_XPETRA_TPETRA 383 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
384 MultiVecTraitsTpetra::SetBlock(toTpetra(A), index, toTpetra(mv));
389 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
394 Assign (
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
395 Xpetra::MultiVector<Scalar,LO,GO,Node>& mv)
398 #ifdef HAVE_XPETRA_TPETRA 399 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
400 MultiVecTraitsTpetra::Assign(toTpetra(A), toTpetra(mv));
405 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
409 static void MvRandom( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
412 #ifdef HAVE_XPETRA_TPETRA 413 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
414 MultiVecTraitsTpetra::MvRandom(toTpetra(mv));
419 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
423 static void MvInit( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, Scalar alpha = Teuchos::ScalarTraits<Scalar>::zero() )
426 #ifdef HAVE_XPETRA_TPETRA 427 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
428 MultiVecTraitsTpetra::MvInit(toTpetra(mv), alpha);
433 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
437 static void MvPrint(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, std::ostream& os )
440 #ifdef HAVE_XPETRA_TPETRA 441 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
442 MultiVecTraitsTpetra::MvPrint(toTpetra(mv), os);
447 XPETRA_FACTORY_ERROR_IF_EPETRA(mv.getMap()->lib());
455 class MultiVecTraits<double, Xpetra::MultiVector<double,int,int,KokkosClassic::DefaultNode::DefaultNodeType> >
463 typedef KokkosClassic::DefaultNode::DefaultNodeType
Node;
465 #ifdef HAVE_XPETRA_TPETRA 466 typedef Xpetra::TpetraMultiVector<Scalar,LO,GO,Node> TpetraMultiVector;
467 typedef MultiVecTraits <Scalar, Tpetra::MultiVector<Scalar,LO,GO,Node> > MultiVecTraitsTpetra;
470 #ifdef HAVE_XPETRA_EPETRA 471 typedef Xpetra::EpetraMultiVectorT<GO,Node> EpetraMultiVector;
472 typedef MultiVecTraits <Scalar, Epetra_MultiVector> MultiVecTraitsEpetra;
477 #ifdef HAVE_BELOS_XPETRA_TIMERS 478 static RCP<Teuchos::Time> mvTimesMatAddMvTimer_, mvTransMvTimer_;
481 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
Clone(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
const int numvecs )
484 #ifdef HAVE_XPETRA_TPETRA 485 if (mv.getMap()->lib() == Xpetra::UseTpetra)
486 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::Clone(toTpetra(mv), numvecs)));
489 #ifdef HAVE_XPETRA_EPETRA 490 if (mv.getMap()->lib() == Xpetra::UseEpetra)
491 return rcp(
new EpetraMultiVector(MultiVecTraitsEpetra::Clone(toEpetra(mv), numvecs)));
497 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
CloneCopy(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
500 #ifdef HAVE_XPETRA_TPETRA 501 if (mv.getMap()->lib() == Xpetra::UseTpetra)
502 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv))));
505 #ifdef HAVE_XPETRA_EPETRA 506 if (mv.getMap()->lib() == Xpetra::UseEpetra)
507 return rcp(
new EpetraMultiVector(MultiVecTraitsEpetra::CloneCopy(toEpetra(mv))));
513 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
CloneCopy(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
const std::vector<int>& index )
516 #ifdef HAVE_XPETRA_TPETRA 517 if (mv.getMap()->lib() == Xpetra::UseTpetra)
518 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv), index)));
521 #ifdef HAVE_XPETRA_EPETRA 522 if (mv.getMap()->lib() == Xpetra::UseEpetra)
523 return rcp(
new EpetraMultiVector(MultiVecTraitsEpetra::CloneCopy(toEpetra(mv), index)));
529 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
530 CloneCopy (
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
531 const Teuchos::Range1D& index)
534 #ifdef HAVE_XPETRA_TPETRA 535 if (mv.getMap()->lib() == Xpetra::UseTpetra)
536 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv), index)));
539 #ifdef HAVE_XPETRA_EPETRA 540 if (mv.getMap()->lib() == Xpetra::UseEpetra)
541 return rcp(
new EpetraMultiVector(MultiVecTraitsEpetra::CloneCopy(toEpetra(mv), index)));
547 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
CloneViewNonConst( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
const std::vector<int>& index )
550 #ifdef HAVE_XPETRA_TPETRA 551 if (mv.getMap()->lib() == Xpetra::UseTpetra)
552 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::CloneViewNonConst(toTpetra(mv), index)));
555 #ifdef HAVE_XPETRA_EPETRA 556 if (mv.getMap()->lib() == Xpetra::UseEpetra)
557 return rcp(
new EpetraMultiVector(MultiVecTraitsEpetra::CloneViewNonConst(toEpetra(mv), index)));
563 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
565 const Teuchos::Range1D& index)
568 #ifdef HAVE_XPETRA_TPETRA 569 if (mv.getMap()->lib() == Xpetra::UseTpetra)
570 return rcp(
new TpetraMultiVector(MultiVecTraitsTpetra::CloneViewNonConst(toTpetra(mv), index)));
573 #ifdef HAVE_XPETRA_EPETRA 574 if (mv.getMap()->lib() == Xpetra::UseEpetra)
575 return rcp(
new EpetraMultiVector(MultiVecTraitsEpetra::CloneViewNonConst(toEpetra(mv), index)));
581 static RCP<const Xpetra::MultiVector<Scalar,LO,GO,Node> >
CloneView(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
const std::vector<int>& index )
584 #ifdef HAVE_XPETRA_TPETRA 585 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
587 RCP<const Tpetra::MultiVector<Scalar,LO,GO,Node> > r = MultiVecTraitsTpetra::CloneView(toTpetra(mv), index);
588 return rcp(
new TpetraMultiVector(Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LO,GO,Node> >(r)));
592 #ifdef HAVE_XPETRA_EPETRA 593 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
595 RCP<const Epetra_MultiVector > r = MultiVecTraitsEpetra::CloneView(toEpetra(mv), index);
596 return rcp(
new EpetraMultiVector(Teuchos::rcp_const_cast<Epetra_MultiVector>(r)));
603 static RCP<const Xpetra::MultiVector<Scalar,LO,GO,Node> >
604 CloneView (
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
605 const Teuchos::Range1D& index)
608 #ifdef HAVE_XPETRA_TPETRA 609 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
611 RCP<const Tpetra::MultiVector<Scalar,LO,GO,Node> > r = MultiVecTraitsTpetra::CloneView(toTpetra(mv), index);
612 return rcp(
new TpetraMultiVector(Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LO,GO,Node> >(r)));
616 #ifdef HAVE_XPETRA_EPETRA 617 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
619 RCP<const Epetra_MultiVector > r = MultiVecTraitsEpetra::CloneView(toEpetra(mv), index);
620 return rcp(
new EpetraMultiVector(Teuchos::rcp_const_cast<Epetra_MultiVector>(r)));
630 #ifdef HAVE_XPETRA_TPETRA 631 if (mv.getMap()->lib() == Xpetra::UseTpetra)
632 return MultiVecTraitsTpetra::GetGlobalLength(toTpetra(mv));
635 #ifdef HAVE_XPETRA_EPETRA 636 if (mv.getMap()->lib() == Xpetra::UseEpetra)
637 return MultiVecTraitsEpetra::GetGlobalLength(toEpetra(mv));
643 static int GetNumberVecs(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
646 #ifdef HAVE_XPETRA_TPETRA 647 if (mv.getMap()->lib() == Xpetra::UseTpetra)
648 return MultiVecTraitsTpetra::GetNumberVecs(toTpetra(mv));
651 #ifdef HAVE_XPETRA_EPETRA 652 if (mv.getMap()->lib() == Xpetra::UseEpetra)
653 return MultiVecTraitsEpetra::GetNumberVecs(toEpetra(mv));
662 #ifdef HAVE_XPETRA_TPETRA 663 if (mv.getMap()->lib() == Xpetra::UseTpetra)
664 return MultiVecTraitsTpetra::HasConstantStride(toTpetra(mv));
667 #ifdef HAVE_XPETRA_EPETRA 668 if (mv.getMap()->lib() == Xpetra::UseEpetra)
669 return MultiVecTraitsEpetra::HasConstantStride(toEpetra(mv));
675 static void MvTimesMatAddMv( Scalar alpha,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
676 const Teuchos::SerialDenseMatrix<int,Scalar>& B,
677 Scalar beta, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
680 #ifdef HAVE_BELOS_XPETRA_TIMERS 681 Teuchos::TimeMonitor lcltimer(*mvTimesMatAddMvTimer_);
684 #ifdef HAVE_XPETRA_TPETRA 685 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
686 MultiVecTraitsTpetra::MvTimesMatAddMv(alpha, toTpetra(A), B, beta, toTpetra(mv));
691 #ifdef HAVE_XPETRA_EPETRA 692 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
693 MultiVecTraitsEpetra::MvTimesMatAddMv(alpha, toEpetra(A), B, beta, toEpetra(mv));
701 static void MvAddMv( Scalar alpha,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, Scalar beta,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
704 #ifdef HAVE_XPETRA_TPETRA 705 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
706 MultiVecTraitsTpetra::MvAddMv(alpha, toTpetra(A), beta, toTpetra(B), toTpetra(mv));
711 #ifdef HAVE_XPETRA_EPETRA 712 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
713 MultiVecTraitsEpetra::MvAddMv(alpha, toEpetra(A), beta, toEpetra(B), toEpetra(mv));
721 static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, Scalar alpha )
724 #ifdef HAVE_XPETRA_TPETRA 725 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
726 MultiVecTraitsTpetra::MvScale(toTpetra(mv), alpha);
731 #ifdef HAVE_XPETRA_EPETRA 732 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
733 MultiVecTraitsEpetra::MvScale(toEpetra(mv), alpha);
741 static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
const std::vector<Scalar>& alphas )
744 #ifdef HAVE_XPETRA_TPETRA 745 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
746 MultiVecTraitsTpetra::MvScale(toTpetra(mv), alphas);
751 #ifdef HAVE_XPETRA_EPETRA 752 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
753 MultiVecTraitsEpetra::MvScale(toEpetra(mv), alphas);
761 static void MvTransMv( Scalar alpha,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, Teuchos::SerialDenseMatrix<int,Scalar>& C)
764 #ifdef HAVE_BELOS_XPETRA_TIMERS 765 Teuchos::TimeMonitor lcltimer(*mvTransMvTimer_);
768 #ifdef HAVE_XPETRA_TPETRA 769 if (A.getMap()->lib() == Xpetra::UseTpetra) {
770 MultiVecTraitsTpetra::MvTransMv(alpha, toTpetra(A), toTpetra(B), C);
775 #ifdef HAVE_XPETRA_EPETRA 776 if (A.getMap()->lib() == Xpetra::UseEpetra) {
777 MultiVecTraitsEpetra::MvTransMv(alpha, toEpetra(A), toEpetra(B), C);
785 static void MvDot(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, std::vector<Scalar> &dots)
788 #ifdef HAVE_XPETRA_TPETRA 789 if (A.getMap()->lib() == Xpetra::UseTpetra) {
790 MultiVecTraitsTpetra::MvDot(toTpetra(A), toTpetra(B), dots);
795 #ifdef HAVE_XPETRA_EPETRA 796 if (A.getMap()->lib() == Xpetra::UseEpetra) {
797 MultiVecTraitsEpetra::MvDot(toEpetra(A), toEpetra(B), dots);
805 static void MvNorm(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, std::vector<Teuchos::ScalarTraits<Scalar>::magnitudeType> &normvec, NormType type=TwoNorm)
808 #ifdef HAVE_XPETRA_TPETRA 809 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
810 MultiVecTraitsTpetra::MvNorm(toTpetra(mv), normvec, type);
815 #ifdef HAVE_XPETRA_EPETRA 816 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
817 MultiVecTraitsEpetra::MvNorm(toEpetra(mv), normvec, type);
825 static void SetBlock(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
const std::vector<int>& index, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
828 #ifdef HAVE_XPETRA_TPETRA 829 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
830 MultiVecTraitsTpetra::SetBlock(toTpetra(A), index, toTpetra(mv));
835 #ifdef HAVE_XPETRA_EPETRA 836 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
837 MultiVecTraitsEpetra::SetBlock(toEpetra(A), index, toEpetra(mv));
846 SetBlock (
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
847 const Teuchos::Range1D& index,
848 Xpetra::MultiVector<Scalar,LO,GO,Node>& mv)
851 #ifdef HAVE_XPETRA_TPETRA 852 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
853 MultiVecTraitsTpetra::SetBlock(toTpetra(A), index, toTpetra(mv));
858 #ifdef HAVE_XPETRA_EPETRA 859 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
860 MultiVecTraitsEpetra::SetBlock(toEpetra(A), index, toEpetra(mv));
869 Assign (
const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
870 Xpetra::MultiVector<Scalar,LO,GO,Node>& mv)
873 #ifdef HAVE_XPETRA_TPETRA 874 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
875 MultiVecTraitsTpetra::Assign(toTpetra(A), toTpetra(mv));
880 #ifdef HAVE_XPETRA_EPETRA 881 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
882 MultiVecTraitsEpetra::Assign(toEpetra(A), toEpetra(mv));
890 static void MvRandom( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
893 #ifdef HAVE_XPETRA_TPETRA 894 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
895 MultiVecTraitsTpetra::MvRandom(toTpetra(mv));
900 #ifdef HAVE_XPETRA_EPETRA 901 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
902 MultiVecTraitsEpetra::MvRandom(toEpetra(mv));
910 static void MvInit( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, Scalar alpha = Teuchos::ScalarTraits<Scalar>::zero() )
913 #ifdef HAVE_XPETRA_TPETRA 914 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
915 MultiVecTraitsTpetra::MvInit(toTpetra(mv), alpha);
920 #ifdef HAVE_XPETRA_EPETRA 921 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
922 MultiVecTraitsEpetra::MvInit(toEpetra(mv), alpha);
930 static void MvPrint(
const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, std::ostream& os )
933 #ifdef HAVE_XPETRA_TPETRA 934 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
935 MultiVecTraitsTpetra::MvPrint(toTpetra(mv), os);
940 #ifdef HAVE_XPETRA_EPETRA 941 if (mv.getMap()->lib() == Xpetra::UseEpetra) {
942 MultiVecTraitsEpetra::MvPrint(toEpetra(mv), os);
static bool HasConstantStride(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneViewNonConst(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static void MvPrint(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, std::ostream &os)
static ptrdiff_t GetGlobalLength(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvRandom(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvInit(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, Scalar alpha=Teuchos::ScalarTraits< Scalar >::zero())
static void MvTransMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, Teuchos::SerialDenseMatrix< int, Scalar > &C)
static RCP< const Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneView(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static void MvTimesMatAddMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Teuchos::SerialDenseMatrix< int, Scalar > &B, Scalar beta, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > Clone(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const int numvecs)
static void MvPrint(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, std::ostream &os)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static void MvInit(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, Scalar alpha=Teuchos::ScalarTraits< Scalar >::zero())
static int GetNumberVecs(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void SetBlock(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const std::vector< int > &index, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvTransMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, Teuchos::SerialDenseMatrix< int, Scalar > &C)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneViewNonConst(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static void MvRandom(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneViewNonConst(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static void MvScale(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, Scalar alpha)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static void MvDot(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, std::vector< Scalar > &dots)
static RCP< const Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneView(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static void MvAddMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, Scalar beta, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< const Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneView(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static int GetNumberVecs(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvScale(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, Scalar alpha)
static void MvScale(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< Scalar > &alphas)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static void SetBlock(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Teuchos::Range1D &index, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
KokkosClassic::DefaultNode::DefaultNodeType Node
static void SetBlock(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const std::vector< int > &index, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvScale(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< Scalar > &alphas)
static ptrdiff_t GetGlobalLength(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvDot(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, std::vector< Scalar > &dots)
static RCP< const Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneView(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< int > &index)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneViewNonConst(Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > Clone(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const int numvecs)
static void Assign(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void SetBlock(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Teuchos::Range1D &index, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvAddMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, Scalar beta, const Xpetra::MultiVector< Scalar, LO, GO, Node > &B, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvTimesMatAddMv(Scalar alpha, const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, const Teuchos::SerialDenseMatrix< int, Scalar > &B, Scalar beta, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvNorm(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, std::vector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType > &normvec, NormType type=TwoNorm)
static bool HasConstantStride(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static RCP< Xpetra::MultiVector< Scalar, LO, GO, Node > > CloneCopy(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Teuchos::Range1D &index)
static void Assign(const Xpetra::MultiVector< Scalar, LO, GO, Node > &A, Xpetra::MultiVector< Scalar, LO, GO, Node > &mv)
static void MvNorm(const Xpetra::MultiVector< Scalar, LO, GO, Node > &mv, std::vector< Teuchos::ScalarTraits< Scalar >::magnitudeType > &normvec, NormType type=TwoNorm)