42#ifndef BELOS_TPETRA_ADAPTER_MP_VECTOR_HPP
43#define BELOS_TPETRA_ADAPTER_MP_VECTOR_HPP
45#include "BelosTpetraAdapter.hpp"
48#include "KokkosBlas.hpp"
72 template<
class Storage,
class LO,
class GO,
class Node>
73 class MultiVecTraits<typename
Storage::value_type,
74 Tpetra::MultiVector< Sacado::MP::Vector<Storage>,
80 typedef Tpetra::MultiVector<Scalar, LO, GO, Node>
MV;
82 typedef typename Tpetra::MultiVector<Scalar,LO,GO,Node>::dot_type
dot_type;
83 typedef typename Tpetra::MultiVector<Scalar,LO,GO,Node>::mag_type
mag_type;
90 static Teuchos::RCP<MV>
Clone (
const MV& X,
const int numVecs) {
91 Teuchos::RCP<MV> Y (
new MV (X.getMap (), numVecs,
false));
92 Y->setCopyOrView (Teuchos::View);
102 Teuchos::RCP<MV> X_copy (
new MV (X, Teuchos::Copy));
109 X_copy->setCopyOrView (Teuchos::View);
124 static Teuchos::RCP<MV>
127#ifdef HAVE_TPETRA_DEBUG
128 const char fnName[] =
"Belos::MultiVecTraits::CloneCopy(mv,index)";
129 const size_t inNumVecs = mv.getNumVectors ();
130 TEUCHOS_TEST_FOR_EXCEPTION(
131 index.size () > 0 && *std::min_element (index.begin (), index.end ()) < 0,
132 std::runtime_error, fnName <<
": All indices must be nonnegative.");
133 TEUCHOS_TEST_FOR_EXCEPTION(
135 static_cast<size_t> (*std::max_element (index.begin (), index.end ())) >= inNumVecs,
137 fnName <<
": All indices must be strictly less than the number of "
138 "columns " << inNumVecs <<
" of the input multivector mv.");
142 Teuchos::Array<size_t> columns (index.size ());
143 for (std::vector<int>::size_type
j = 0;
j < index.size (); ++
j) {
144 columns[
j] = index[
j];
149 Teuchos::RCP<MV> X_copy = mv.subCopy (columns ());
150 X_copy->setCopyOrView (Teuchos::View);
160 static Teuchos::RCP<MV>
163 const bool validRange = index.size() > 0 &&
164 index.lbound() >= 0 &&
165 index.ubound() < GetNumberVecs(mv);
167 std::ostringstream os;
168 os <<
"Belos::MultiVecTraits::CloneCopy(mv,index=["
169 << index.lbound() <<
"," << index.ubound() <<
"]): ";
170 TEUCHOS_TEST_FOR_EXCEPTION(
171 index.size() == 0, std::invalid_argument,
172 os.str() <<
"Empty index range is not allowed.");
173 TEUCHOS_TEST_FOR_EXCEPTION(
174 index.lbound() < 0, std::invalid_argument,
175 os.str() <<
"Index range includes negative index/ices, which is not "
177 TEUCHOS_TEST_FOR_EXCEPTION(
178 index.ubound() >= GetNumberVecs(mv), std::invalid_argument,
179 os.str() <<
"Index range exceeds number of vectors "
180 << mv.getNumVectors() <<
" in the input multivector.");
181 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
182 os.str() <<
"Should never get here!");
184 Teuchos::RCP<MV> X_copy = mv.subCopy (index);
185 X_copy->setCopyOrView (Teuchos::View);
190 static Teuchos::RCP<MV>
193#ifdef HAVE_TPETRA_DEBUG
194 const char fnName[] =
"Belos::MultiVecTraits::CloneViewNonConst(mv,index)";
195 const size_t numVecs = mv.getNumVectors ();
196 TEUCHOS_TEST_FOR_EXCEPTION(
197 index.size () > 0 && *std::min_element (index.begin (), index.end ()) < 0,
198 std::invalid_argument,
199 fnName <<
": All indices must be nonnegative.");
200 TEUCHOS_TEST_FOR_EXCEPTION(
202 static_cast<size_t> (*std::max_element (index.begin (), index.end ())) >= numVecs,
203 std::invalid_argument,
204 fnName <<
": All indices must be strictly less than the number of "
205 "columns " << numVecs <<
" in the input MultiVector mv.");
209 Teuchos::Array<size_t> columns (index.size ());
210 for (std::vector<int>::size_type
j = 0;
j < index.size (); ++
j) {
211 columns[
j] = index[
j];
216 Teuchos::RCP<MV> X_view = mv.subViewNonConst (columns ());
217 X_view->setCopyOrView (Teuchos::View);
221 static Teuchos::RCP<MV>
227 const int numCols =
static_cast<int> (mv.getNumVectors());
228 const bool validRange = index.size() > 0 &&
229 index.lbound() >= 0 && index.ubound() < numCols;
231 std::ostringstream os;
232 os <<
"Belos::MultiVecTraits::CloneViewNonConst(mv,index=["
233 << index.lbound() <<
", " << index.ubound() <<
"]): ";
234 TEUCHOS_TEST_FOR_EXCEPTION(
235 index.size() == 0, std::invalid_argument,
236 os.str() <<
"Empty index range is not allowed.");
237 TEUCHOS_TEST_FOR_EXCEPTION(
238 index.lbound() < 0, std::invalid_argument,
239 os.str() <<
"Index range includes negative inde{x,ices}, which is "
241 TEUCHOS_TEST_FOR_EXCEPTION(
242 index.ubound() >= numCols, std::invalid_argument,
243 os.str() <<
"Index range exceeds number of vectors " << numCols
244 <<
" in the input multivector.");
245 TEUCHOS_TEST_FOR_EXCEPTION(
246 true, std::logic_error,
247 os.str() <<
"Should never get here!");
249 Teuchos::RCP<MV> X_view = mv.subViewNonConst (index);
250 X_view->setCopyOrView (Teuchos::View);
254 static Teuchos::RCP<const MV>
257#ifdef HAVE_TPETRA_DEBUG
258 const char fnName[] =
"Belos::MultiVecTraits<Scalar, "
259 "Tpetra::MultiVector<...> >::CloneView(mv,index)";
260 const size_t numVecs = mv.getNumVectors ();
261 TEUCHOS_TEST_FOR_EXCEPTION(
262 *std::min_element (index.begin (), index.end ()) < 0,
263 std::invalid_argument,
264 fnName <<
": All indices must be nonnegative.");
265 TEUCHOS_TEST_FOR_EXCEPTION(
266 static_cast<size_t> (*std::max_element (index.begin (), index.end ())) >= numVecs,
267 std::invalid_argument,
268 fnName <<
": All indices must be strictly less than the number of "
269 "columns " << numVecs <<
" in the input MultiVector mv.");
273 Teuchos::Array<size_t> columns (index.size ());
274 for (std::vector<int>::size_type
j = 0;
j < index.size (); ++
j) {
275 columns[
j] = index[
j];
280 Teuchos::RCP<const MV> X_view = mv.subView (columns);
281 Teuchos::rcp_const_cast<MV> (X_view)->setCopyOrView (Teuchos::View);
285 static Teuchos::RCP<const MV>
291 const int numCols =
static_cast<int> (mv.getNumVectors());
292 const bool validRange = index.size() > 0 &&
293 index.lbound() >= 0 && index.ubound() < numCols;
295 std::ostringstream os;
296 os <<
"Belos::MultiVecTraits::CloneView(mv, index=["
297 << index.lbound () <<
", " << index.ubound() <<
"]): ";
298 TEUCHOS_TEST_FOR_EXCEPTION(index.size() == 0, std::invalid_argument,
299 os.str() <<
"Empty index range is not allowed.");
300 TEUCHOS_TEST_FOR_EXCEPTION(index.lbound() < 0, std::invalid_argument,
301 os.str() <<
"Index range includes negative index/ices, which is not "
303 TEUCHOS_TEST_FOR_EXCEPTION(
304 index.ubound() >= numCols, std::invalid_argument,
305 os.str() <<
"Index range exceeds number of vectors " << numCols
306 <<
" in the input multivector.");
307 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
308 os.str() <<
"Should never get here!");
310 Teuchos::RCP<const MV> X_view = mv.subView (index);
311 Teuchos::rcp_const_cast<MV> (X_view)->setCopyOrView (Teuchos::View);
316 return static_cast<ptrdiff_t
> (mv.getGlobalLength ());
320 return static_cast<int> (mv.getNumVectors ());
324 return mv.isConstantStride ();
327 template <
class DstType,
class SrcType>
333 if (std::is_same<typename SrcType::memory_space, Kokkos::HostSpace>::value)
336 for (
size_t i = 0; i < src.extent(0); i++)
338 for (
size_t j = 0;
j < src.extent(1);
j++)
340 dst_i(i,
j) = src(i,
j);
349 for (
size_t i = 0; i < src.extent(0); i++)
351 for (
size_t j = 0;
j < src.extent(1);
j++)
353 dst(i,
j) = src_i(i,
j);
361 const Tpetra::MultiVector<Scalar,LO,GO,Node>& A,
362 const Teuchos::SerialDenseMatrix<int,dot_type>& B,
364 Tpetra::MultiVector<Scalar,LO,GO,Node>& C)
370 const int numRowsB = B.numRows ();
371 const int numColsB = B.numCols ();
372 const int strideB = B.stride ();
373 if (numRowsB == 1 && numColsB == 1) {
374 C.update (alpha*B(0,0), A, beta);
381 if (A.isConstantStride() ==
false) Atmp = rcp (
new MV (A, Teuchos::Copy));
382 else Atmp = rcp(&A,
false);
384 if (C.isConstantStride() ==
false) Ctmp = rcp (
new MV (C, Teuchos::Copy));
385 else Ctmp = rcp(&C,
false);
388 typedef Tpetra::MultiVector<dot_type,LO,GO,Node> FMV;
389 typedef typename FMV::dual_view_type::t_dev flat_view_type;
391 typename flat_view_type::const_type flat_A_view = Atmp->getLocalViewDevice(Tpetra::Access::ReadOnly);
392 flat_view_type flat_C_view = Ctmp->getLocalViewDevice(Tpetra::Access::OverwriteAll);
395 typedef Kokkos::View<dot_type**, Kokkos::LayoutLeft, Kokkos::HostSpace> b_host_view_type;
396 b_host_view_type B_view_host_input( B.values(), strideB, numColsB);
397 auto B_view_host = Kokkos::subview(B_view_host_input,
398 Kokkos::pair<int,int>(0,numRowsB),
399 Kokkos::pair<int,int>(0,numColsB));
403 typedef Kokkos::View<dot_type**, Kokkos::LayoutLeft, execution_space> b_view_type;
404 typedef Kokkos::View<dot_type*, Kokkos::LayoutLeft, execution_space> b_1d_view_type;
405 b_1d_view_type B_1d_view_dev(Kokkos::ViewAllocateWithoutInitializing(
"B"), numRowsB*numColsB);
406 b_view_type B_view_dev( B_1d_view_dev.data(), numRowsB, numColsB);
408 if (Kokkos::SpaceAccessibility<Kokkos::HostSpace, typename execution_space::memory_space>::accessible)
414 deep_copy_2d_view_with_intercessory_space(B_view_dev, B_view_host);
419 const char ctransA =
'N', ctransB =
'N';
422 alpha, flat_A_view, B_view_dev, beta, flat_C_view);
425 if (C.isConstantStride() ==
false)
438 const Tpetra::MultiVector<Scalar,LO,GO,Node>& A,
440 const Tpetra::MultiVector<Scalar,LO,GO,Node>& B,
441 Tpetra::MultiVector<Scalar,LO,GO,Node>& mv)
443 mv.update (alpha, A, beta, B, Teuchos::ScalarTraits<Scalar>::zero ());
447 MvScale (Tpetra::MultiVector<Scalar,LO,GO,Node>& mv,
454 MvScale (Tpetra::MultiVector<Scalar,LO,GO,Node>& mv,
455 const std::vector<BaseScalar>& alphas)
457 std::vector<Scalar> alphas_mp(alphas.size());
458 const size_t sz = alphas.size();
459 for (
size_t i=0; i<sz; ++i)
460 alphas_mp[i] = alphas[i];
461 mv.scale (alphas_mp);
465 MvScale (Tpetra::MultiVector<Scalar,LO,GO,Node>& mv,
466 const std::vector<Scalar>& alphas)
473 const Tpetra::MultiVector<Scalar,LO,GO,Node>& A,
474 const Tpetra::MultiVector<Scalar,LO,GO,Node>& B,
475 Teuchos::SerialDenseMatrix<int,dot_type>& C)
480 using Teuchos::REDUCE_SUM;
481 using Teuchos::reduceAll;
484 const int numRowsC = C.numRows ();
485 const int numColsC = C.numCols ();
486 const int strideC = C.stride ();
487 if (numRowsC == 1 && numColsC == 1) {
488 if (alpha == Teuchos::ScalarTraits<Scalar>::zero ()) {
493 A.dot (B, Teuchos::ArrayView<dot_type> (C.values (), 1));
494 if (alpha != Teuchos::ScalarTraits<Scalar>::one ()) {
501 RCP<const MV> Atmp, Btmp;
502 if (A.isConstantStride() ==
false) Atmp = rcp (
new MV (A, Teuchos::Copy));
503 else Atmp = rcp(&A,
false);
505 if (B.isConstantStride() ==
false) Btmp = rcp (
new MV (B, Teuchos::Copy));
506 else Btmp = rcp(&B,
false);
509 typedef Tpetra::MultiVector<dot_type,LO,GO,Node> FMV;
510 typedef typename FMV::dual_view_type::t_dev flat_view_type;
512 typename flat_view_type::const_type flat_A_view = Atmp->getLocalViewDevice(Tpetra::Access::ReadOnly);
513 typename flat_view_type::const_type flat_B_view = Btmp->getLocalViewDevice(Tpetra::Access::ReadOnly);
516 typedef Kokkos::View<dot_type**, Kokkos::LayoutLeft, Kokkos::HostSpace> c_host_view_type;
517 c_host_view_type C_view_host_input( C.values(), strideC, numColsC);
518 auto C_view_host = Kokkos::subview(C_view_host_input,
519 Kokkos::pair<int,int>(0,numRowsC),
520 Kokkos::pair<int,int>(0,numColsC));
524 typedef Kokkos::View<dot_type**, Kokkos::LayoutLeft, execution_space> c_view_type;
525 typedef Kokkos::View<dot_type*, Kokkos::LayoutLeft, execution_space> c_1d_view_type;
526 c_1d_view_type C_1d_view_dev(
"C", numRowsC*numColsC);
527 c_view_type C_view_dev( C_1d_view_dev.data(), numRowsC, numColsC);
531 const char ctransA =
'C', ctransB =
'N';
534 alpha, flat_A_view, flat_B_view,
535 Kokkos::Details::ArithTraits<dot_type>::zero(),
540 RCP<const Comm<int> > pcomm = A.getMap()->getComm ();
541 if (pcomm->getSize () == 1)
543 if (Kokkos::SpaceAccessibility<execution_space, Kokkos::HostSpace>::accessible)
549 deep_copy_2d_view_with_intercessory_space(C_view_host, C_view_dev);
554 typedef Kokkos::View<dot_type*, Kokkos::LayoutLeft, Kokkos::HostSpace> c_1d_host_view_type;
555 c_1d_host_view_type C_1d_view_tmp(Kokkos::ViewAllocateWithoutInitializing(
"C_tmp"), strideC*numColsC);
556 c_host_view_type C_view_tmp_input( C_1d_view_tmp.data(), strideC, numColsC);
557 auto C_view_tmp = Kokkos::subview(C_view_tmp_input,
558 Kokkos::pair<int,int>(0,numRowsC),
559 Kokkos::pair<int,int>(0,numColsC));
561 if (Kokkos::SpaceAccessibility<execution_space, Kokkos::HostSpace>::accessible)
567 deep_copy_2d_view_with_intercessory_space(C_view_tmp, C_view_dev);
570 reduceAll<int> (*pcomm, REDUCE_SUM, strideC*numColsC,
578 MvDot (
const Tpetra::MultiVector<Scalar,LO,GO,Node>& A,
579 const Tpetra::MultiVector<Scalar,LO,GO,Node>& B,
580 std::vector<dot_type>& dots)
582 const size_t numVecs = A.getNumVectors ();
584 TEUCHOS_TEST_FOR_EXCEPTION(
585 numVecs != B.getNumVectors (), std::invalid_argument,
586 "Belos::MultiVecTraits<Scalar,Tpetra::MultiVector>::MvDot(A,B,dots): "
587 "A and B must have the same number of columns. "
588 "A has " << numVecs <<
" column(s), "
589 "but B has " << B.getNumVectors () <<
" column(s).");
590#ifdef HAVE_TPETRA_DEBUG
591 TEUCHOS_TEST_FOR_EXCEPTION(
592 dots.size() < numVecs, std::invalid_argument,
593 "Belos::MultiVecTraits<Scalar,Tpetra::MultiVector>::MvDot(A,B,dots): "
594 "The output array 'dots' must have room for all dot products. "
595 "A and B each have " << numVecs <<
" column(s), "
596 "but 'dots' only has " << dots.size() <<
" entry(/ies).");
599 Teuchos::ArrayView<dot_type> av (dots);
600 A.dot (B, av (0, numVecs));
605 MvNorm (
const Tpetra::MultiVector<Scalar,LO,GO,Node>& mv,
606 std::vector<mag_type> &normvec,
607 NormType type=TwoNorm)
610#ifdef HAVE_TPETRA_DEBUG
611 typedef std::vector<int>::size_type size_type;
612 TEUCHOS_TEST_FOR_EXCEPTION(
613 normvec.size () <
static_cast<size_type
> (mv.getNumVectors ()),
614 std::invalid_argument,
615 "Belos::MultiVecTraits::MvNorm(mv,normvec): The normvec output "
616 "argument must have at least as many entries as the number of vectors "
617 "(columns) in the MultiVector mv. normvec.size() = " << normvec.size ()
618 <<
" < mv.getNumVectors() = " << mv.getNumVectors () <<
".");
621 Teuchos::ArrayView<mag_type> av(normvec);
624 mv.norm1(av(0,mv.getNumVectors()));
627 mv.norm2(av(0,mv.getNumVectors()));
630 mv.normInf(av(0,mv.getNumVectors()));
636 TEUCHOS_TEST_FOR_EXCEPTION(
637 true, std::logic_error,
638 "Belos::MultiVecTraits::MvNorm: Invalid NormType value " << type
639 <<
". Valid values are OneNorm=" << OneNorm <<
", TwoNorm="
640 << TwoNorm <<
", and InfNorm=" << InfNorm <<
". If you are a Belos "
641 "user and have not modified Belos in any way, and you get this "
642 "message, then this is probably a bug in the Belos solver you were "
643 "using. Please report this to the Belos developers.");
650 using Teuchos::Range1D;
652 const size_t inNumVecs = A.getNumVectors ();
653#ifdef HAVE_TPETRA_DEBUG
654 TEUCHOS_TEST_FOR_EXCEPTION(
655 inNumVecs <
static_cast<size_t> (index.size ()), std::invalid_argument,
656 "Belos::MultiVecTraits::SetBlock(A,index,mv): 'index' argument must "
657 "have no more entries as the number of columns in the input MultiVector"
658 " A. A.getNumVectors() = " << inNumVecs <<
" < index.size () = "
659 << index.size () <<
".");
661 RCP<MV> mvsub = CloneViewNonConst (mv, index);
662 if (inNumVecs >
static_cast<size_t> (index.size ())) {
663 RCP<const MV> Asub = A.subView (Range1D (0, index.size () - 1));
664 ::Tpetra::deep_copy (*mvsub, *Asub);
666 ::Tpetra::deep_copy (*mvsub, A);
680 const size_t maxInt =
681 static_cast<size_t> (Teuchos::OrdinalTraits<int>::max ());
682 const bool overflow =
683 maxInt < A.getNumVectors () && maxInt < mv.getNumVectors ();
685 std::ostringstream os;
686 os <<
"Belos::MultiVecTraits::SetBlock(A, index=[" << index.lbound ()
687 <<
", " << index.ubound () <<
"], mv): ";
688 TEUCHOS_TEST_FOR_EXCEPTION(
689 maxInt < A.getNumVectors (), std::range_error, os.str () <<
"Number "
690 "of columns (size_t) in the input MultiVector 'A' overflows int.");
691 TEUCHOS_TEST_FOR_EXCEPTION(
692 maxInt < mv.getNumVectors (), std::range_error, os.str () <<
"Number "
693 "of columns (size_t) in the output MultiVector 'mv' overflows int.");
696 const int numColsA =
static_cast<int> (A.getNumVectors ());
697 const int numColsMv =
static_cast<int> (mv.getNumVectors ());
699 const bool validIndex =
700 index.lbound () >= 0 && index.ubound () < numColsMv;
702 const bool validSource = index.size () <= numColsA;
704 if (! validIndex || ! validSource) {
705 std::ostringstream os;
706 os <<
"Belos::MultiVecTraits::SetBlock(A, index=[" << index.lbound ()
707 <<
", " << index.ubound () <<
"], mv): ";
708 TEUCHOS_TEST_FOR_EXCEPTION(
709 index.lbound() < 0, std::invalid_argument,
710 os.str() <<
"Range lower bound must be nonnegative.");
711 TEUCHOS_TEST_FOR_EXCEPTION(
712 index.ubound() >= numColsMv, std::invalid_argument,
713 os.str() <<
"Range upper bound must be less than the number of "
714 "columns " << numColsA <<
" in the 'mv' output argument.");
715 TEUCHOS_TEST_FOR_EXCEPTION(
716 index.size() > numColsA, std::invalid_argument,
717 os.str() <<
"Range must have no more elements than the number of "
718 "columns " << numColsA <<
" in the 'A' input argument.");
719 TEUCHOS_TEST_FOR_EXCEPTION(
720 true, std::logic_error,
"Should never get here!");
726 Teuchos::RCP<MV> mv_view;
727 if (index.lbound () == 0 && index.ubound () + 1 == numColsMv) {
728 mv_view = Teuchos::rcpFromRef (mv);
730 mv_view = CloneViewNonConst (mv, index);
736 Teuchos::RCP<const MV> A_view;
737 if (index.size () == numColsA) {
738 A_view = Teuchos::rcpFromRef (A);
740 A_view = CloneView (A, Teuchos::Range1D (0, index.size () - 1));
743 ::Tpetra::deep_copy (*mv_view, *A_view);
748 const char errPrefix[] =
"Belos::MultiVecTraits::Assign(A, mv): ";
757 const size_t maxInt =
758 static_cast<size_t> (Teuchos::OrdinalTraits<int>::max ());
759 const bool overflow =
760 maxInt < A.getNumVectors () && maxInt < mv.getNumVectors ();
762 TEUCHOS_TEST_FOR_EXCEPTION(
763 maxInt < A.getNumVectors(), std::range_error,
764 errPrefix <<
"Number of columns in the input multivector 'A' "
765 "(a size_t) overflows int.");
766 TEUCHOS_TEST_FOR_EXCEPTION(
767 maxInt < mv.getNumVectors(), std::range_error,
768 errPrefix <<
"Number of columns in the output multivector 'mv' "
769 "(a size_t) overflows int.");
770 TEUCHOS_TEST_FOR_EXCEPTION(
771 true, std::logic_error,
"Should never get here!");
774 const int numColsA =
static_cast<int> (A.getNumVectors ());
775 const int numColsMv =
static_cast<int> (mv.getNumVectors ());
776 if (numColsA > numColsMv) {
777 TEUCHOS_TEST_FOR_EXCEPTION(
778 numColsA > numColsMv, std::invalid_argument,
779 errPrefix <<
"Input multivector 'A' has " << numColsA <<
" columns, "
780 "but output multivector 'mv' has only " << numColsMv <<
" columns.");
781 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
"Should never get here!");
783 if (numColsA == numColsMv) {
784 ::Tpetra::deep_copy (mv, A);
786 Teuchos::RCP<MV> mv_view =
787 CloneViewNonConst (mv, Teuchos::Range1D (0, numColsA - 1));
788 ::Tpetra::deep_copy (*mv_view, A);
800 mv.putScalar (alpha);
804 Teuchos::FancyOStream fos (Teuchos::rcpFromRef (os));
805 mv.describe (fos, Teuchos::VERB_EXTREME);
808#ifdef HAVE_BELOS_TSQR
812 typedef Tpetra::TsqrAdaptor< Tpetra::MultiVector< Scalar, LO, GO, Node > > tsqr_adaptor_type;
823 template <
class Storage,
class LO,
class GO,
class Node>
824 class OperatorTraits <typename
Storage::value_type,
825 Tpetra::MultiVector<Sacado::MP::Vector<Storage>,
827 Tpetra::Operator<Sacado::MP::Vector<Storage>,
833 Apply (
const Tpetra::Operator<Scalar,LO,GO,Node>& Op,
834 const Tpetra::MultiVector<Scalar,LO,GO,Node>& X,
835 Tpetra::MultiVector<Scalar,LO,GO,Node>& Y,
836 ETrans trans=NOTRANS)
840 Op.apply(X,Y,Teuchos::NO_TRANS);
843 Op.apply(X,Y,Teuchos::TRANS);
846 Op.apply(X,Y,Teuchos::CONJ_TRANS);
849 const std::string scalarName = Teuchos::TypeNameTraits<Scalar>::name();
850 const std::string loName = Teuchos::TypeNameTraits<LO>::name();
851 const std::string goName = Teuchos::TypeNameTraits<GO>::name();
852 const std::string nodeName = Teuchos::TypeNameTraits<Node>::name();
853 const std::string otName =
"Belos::OperatorTraits<" + scalarName
854 +
"," + loName +
"," + goName +
"," + nodeName +
">";
855 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error, otName <<
": Should never "
856 "get here; fell through a switch statement. "
857 "Please report this bug to the Belos developers.");
864 return Op.hasTransposeApply ();
Kokkos::DefaultExecutionSpace execution_space
static int GetNumberVecs(const MV &mv)
static void MvPrint(const MV &mv, std::ostream &os)
static void MvTimesMatAddMv(const dot_type &alpha, const Tpetra::MultiVector< Scalar, LO, GO, Node > &A, const Teuchos::SerialDenseMatrix< int, dot_type > &B, const dot_type &beta, Tpetra::MultiVector< Scalar, LO, GO, Node > &C)
static void MvRandom(MV &mv)
static void MvInit(MV &mv, const Scalar alpha=Teuchos::ScalarTraits< Scalar >::zero())
static void MvNorm(const Tpetra::MultiVector< Scalar, LO, GO, Node > &mv, std::vector< mag_type > &normvec, NormType type=TwoNorm)
For all columns j of mv, set normvec[j] = norm(mv[j]).
static ptrdiff_t GetGlobalLength(const MV &mv)
Storage::ordinal_type s_ordinal
Tpetra::MultiVector< Scalar, LO, GO, Node >::mag_type mag_type
static void MvAddMv(Scalar alpha, const Tpetra::MultiVector< Scalar, LO, GO, Node > &A, Scalar beta, const Tpetra::MultiVector< Scalar, LO, GO, Node > &B, Tpetra::MultiVector< Scalar, LO, GO, Node > &mv)
mv := alpha*A + beta*B
static void SetBlock(const MV &A, const Teuchos::Range1D &index, MV &mv)
static Teuchos::RCP< MV > CloneCopy(const MV &mv, const std::vector< int > &index)
Create and return a deep copy of the given columns of mv.
static Teuchos::RCP< MV > Clone(const MV &X, const int numVecs)
Create a new MultiVector with numVecs columns.
static Teuchos::RCP< const MV > CloneView(const MV &mv, const Teuchos::Range1D &index)
Tpetra::MultiVector< Scalar, LO, GO, Node >::dot_type dot_type
Storage::value_type BaseScalar
Tpetra::MultiVector< Scalar, LO, GO, Node > MV
Sacado::MP::Vector< Storage > Scalar
static Teuchos::RCP< MV > CloneViewNonConst(MV &mv, const std::vector< int > &index)
static Teuchos::RCP< const MV > CloneView(const MV &mv, const std::vector< int > &index)
static void MvDot(const Tpetra::MultiVector< Scalar, LO, GO, Node > &A, const Tpetra::MultiVector< Scalar, LO, GO, Node > &B, std::vector< dot_type > &dots)
For all columns j of A, set dots[j] := A[j]^T * B[j].
static Teuchos::RCP< MV > CloneCopy(const MV &X)
Create and return a deep copy of X.
static void SetBlock(const MV &A, const std::vector< int > &index, MV &mv)
static void MvTransMv(dot_type alpha, const Tpetra::MultiVector< Scalar, LO, GO, Node > &A, const Tpetra::MultiVector< Scalar, LO, GO, Node > &B, Teuchos::SerialDenseMatrix< int, dot_type > &C)
static void Assign(const MV &A, MV &mv)
static void MvScale(Tpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< Scalar > &alphas)
static void MvScale(Tpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< BaseScalar > &alphas)
static Teuchos::RCP< MV > CloneViewNonConst(MV &mv, const Teuchos::Range1D &index)
static void MvScale(Tpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Scalar &alpha)
static void deep_copy_2d_view_with_intercessory_space(const DstType &dst, const SrcType &src)
static bool HasConstantStride(const MV &mv)
static Teuchos::RCP< MV > CloneCopy(const MV &mv, const Teuchos::Range1D &index)
Create and return a deep copy of the given columns of mv.
static bool HasApplyTranspose(const Tpetra::Operator< Scalar, LO, GO, Node > &Op)
static void Apply(const Tpetra::Operator< Scalar, LO, GO, Node > &Op, const Tpetra::MultiVector< Scalar, LO, GO, Node > &X, Tpetra::MultiVector< Scalar, LO, GO, Node > &Y, ETrans trans=NOTRANS)
Sacado::MP::Vector< Storage > Scalar
std::enable_if< Kokkos::is_view_mp_vector< Kokkos::View< DA, PA... > >::value &&Kokkos::is_view_mp_vector< Kokkos::View< DB, PB... > >::value &&Kokkos::is_view_mp_vector< Kokkos::View< DC, PC... > >::value >::type gemm(const char transA[], const char transB[], typename Kokkos::View< DA, PA... >::const_value_type &alpha, const Kokkos::View< DA, PA... > &A, const Kokkos::View< DB, PB... > &B, typename Kokkos::View< DC, PC... >::const_value_type &beta, const Kokkos::View< DC, PC... > &C)
void deep_copy(const Stokhos::CrsMatrix< ValueType, DstDevice, Layout > &dst, const Stokhos::CrsMatrix< ValueType, SrcDevice, Layout > &src)
Stokhos::CrsMatrix< ValueType, Device, Layout >::HostMirror create_mirror(const Stokhos::CrsMatrix< ValueType, Device, Layout > &A)