42#ifndef BELOS_TPETRA_ADAPTER_UQ_PCE_HPP
43#define BELOS_TPETRA_ADAPTER_UQ_PCE_HPP
45#include "BelosTpetraAdapter.hpp"
47#include "KokkosBlas.hpp"
71 template<
class Storage,
class LO,
class GO,
class Node>
72 class MultiVecTraits<typename
Storage::value_type,
73 Tpetra::MultiVector< Sacado::UQ::PCE<Storage>,
79 typedef Tpetra::MultiVector<Scalar, LO, GO, Node>
MV;
81 typedef typename Tpetra::MultiVector<Scalar,LO,GO,Node>::dot_type
dot_type;
82 typedef typename Tpetra::MultiVector<Scalar,LO,GO,Node>::mag_type
mag_type;
89 static Teuchos::RCP<MV>
Clone (
const MV& X,
const int numVecs) {
90 Teuchos::RCP<MV> Y (
new MV (X.getMap (), numVecs,
false));
91 Y->setCopyOrView (Teuchos::View);
101 Teuchos::RCP<MV> X_copy (
new MV (X, Teuchos::Copy));
108 X_copy->setCopyOrView (Teuchos::View);
123 static Teuchos::RCP<MV>
126#ifdef HAVE_TPETRA_DEBUG
127 const char fnName[] =
"Belos::MultiVecTraits::CloneCopy(mv,index)";
128 const size_t inNumVecs = mv.getNumVectors ();
129 TEUCHOS_TEST_FOR_EXCEPTION(
130 index.size () > 0 && *std::min_element (index.begin (), index.end ()) < 0,
131 std::runtime_error, fnName <<
": All indices must be nonnegative.");
132 TEUCHOS_TEST_FOR_EXCEPTION(
134 static_cast<size_t> (*std::max_element (index.begin (), index.end ())) >= inNumVecs,
136 fnName <<
": All indices must be strictly less than the number of "
137 "columns " << inNumVecs <<
" of the input multivector mv.");
141 Teuchos::Array<size_t> columns (index.size ());
142 for (std::vector<int>::size_type
j = 0;
j < index.size (); ++
j) {
143 columns[
j] = index[
j];
148 Teuchos::RCP<MV> X_copy = mv.subCopy (columns ());
149 X_copy->setCopyOrView (Teuchos::View);
159 static Teuchos::RCP<MV>
162 const bool validRange = index.size() > 0 &&
163 index.lbound() >= 0 &&
164 index.ubound() < GetNumberVecs(mv);
166 std::ostringstream os;
167 os <<
"Belos::MultiVecTraits::CloneCopy(mv,index=["
168 << index.lbound() <<
"," << index.ubound() <<
"]): ";
169 TEUCHOS_TEST_FOR_EXCEPTION(
170 index.size() == 0, std::invalid_argument,
171 os.str() <<
"Empty index range is not allowed.");
172 TEUCHOS_TEST_FOR_EXCEPTION(
173 index.lbound() < 0, std::invalid_argument,
174 os.str() <<
"Index range includes negative index/ices, which is not "
176 TEUCHOS_TEST_FOR_EXCEPTION(
177 index.ubound() >= GetNumberVecs(mv), std::invalid_argument,
178 os.str() <<
"Index range exceeds number of vectors "
179 << mv.getNumVectors() <<
" in the input multivector.");
180 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
181 os.str() <<
"Should never get here!");
183 Teuchos::RCP<MV> X_copy = mv.subCopy (index);
184 X_copy->setCopyOrView (Teuchos::View);
189 static Teuchos::RCP<MV>
192#ifdef HAVE_TPETRA_DEBUG
193 const char fnName[] =
"Belos::MultiVecTraits::CloneViewNonConst(mv,index)";
194 const size_t numVecs = mv.getNumVectors ();
195 TEUCHOS_TEST_FOR_EXCEPTION(
196 index.size () > 0 && *std::min_element (index.begin (), index.end ()) < 0,
197 std::invalid_argument,
198 fnName <<
": All indices must be nonnegative.");
199 TEUCHOS_TEST_FOR_EXCEPTION(
201 static_cast<size_t> (*std::max_element (index.begin (), index.end ())) >= numVecs,
202 std::invalid_argument,
203 fnName <<
": All indices must be strictly less than the number of "
204 "columns " << numVecs <<
" in the input MultiVector mv.");
208 Teuchos::Array<size_t> columns (index.size ());
209 for (std::vector<int>::size_type
j = 0;
j < index.size (); ++
j) {
210 columns[
j] = index[
j];
215 Teuchos::RCP<MV> X_view = mv.subViewNonConst (columns ());
216 X_view->setCopyOrView (Teuchos::View);
220 static Teuchos::RCP<MV>
226 const int numCols =
static_cast<int> (mv.getNumVectors());
227 const bool validRange = index.size() > 0 &&
228 index.lbound() >= 0 && index.ubound() < numCols;
230 std::ostringstream os;
231 os <<
"Belos::MultiVecTraits::CloneViewNonConst(mv,index=["
232 << index.lbound() <<
", " << index.ubound() <<
"]): ";
233 TEUCHOS_TEST_FOR_EXCEPTION(
234 index.size() == 0, std::invalid_argument,
235 os.str() <<
"Empty index range is not allowed.");
236 TEUCHOS_TEST_FOR_EXCEPTION(
237 index.lbound() < 0, std::invalid_argument,
238 os.str() <<
"Index range includes negative inde{x,ices}, which is "
240 TEUCHOS_TEST_FOR_EXCEPTION(
241 index.ubound() >= numCols, std::invalid_argument,
242 os.str() <<
"Index range exceeds number of vectors " << numCols
243 <<
" in the input multivector.");
244 TEUCHOS_TEST_FOR_EXCEPTION(
245 true, std::logic_error,
246 os.str() <<
"Should never get here!");
248 Teuchos::RCP<MV> X_view = mv.subViewNonConst (index);
249 X_view->setCopyOrView (Teuchos::View);
253 static Teuchos::RCP<const MV>
256#ifdef HAVE_TPETRA_DEBUG
257 const char fnName[] =
"Belos::MultiVecTraits<Scalar, "
258 "Tpetra::MultiVector<...> >::CloneView(mv,index)";
259 const size_t numVecs = mv.getNumVectors ();
260 TEUCHOS_TEST_FOR_EXCEPTION(
261 *std::min_element (index.begin (), index.end ()) < 0,
262 std::invalid_argument,
263 fnName <<
": All indices must be nonnegative.");
264 TEUCHOS_TEST_FOR_EXCEPTION(
265 static_cast<size_t> (*std::max_element (index.begin (), index.end ())) >= numVecs,
266 std::invalid_argument,
267 fnName <<
": All indices must be strictly less than the number of "
268 "columns " << numVecs <<
" in the input MultiVector mv.");
272 Teuchos::Array<size_t> columns (index.size ());
273 for (std::vector<int>::size_type
j = 0;
j < index.size (); ++
j) {
274 columns[
j] = index[
j];
279 Teuchos::RCP<const MV> X_view = mv.subView (columns);
280 Teuchos::rcp_const_cast<MV> (X_view)->setCopyOrView (Teuchos::View);
284 static Teuchos::RCP<const MV>
290 const int numCols =
static_cast<int> (mv.getNumVectors());
291 const bool validRange = index.size() > 0 &&
292 index.lbound() >= 0 && index.ubound() < numCols;
294 std::ostringstream os;
295 os <<
"Belos::MultiVecTraits::CloneView(mv, index=["
296 << index.lbound () <<
", " << index.ubound() <<
"]): ";
297 TEUCHOS_TEST_FOR_EXCEPTION(index.size() == 0, std::invalid_argument,
298 os.str() <<
"Empty index range is not allowed.");
299 TEUCHOS_TEST_FOR_EXCEPTION(index.lbound() < 0, std::invalid_argument,
300 os.str() <<
"Index range includes negative index/ices, which is not "
302 TEUCHOS_TEST_FOR_EXCEPTION(
303 index.ubound() >= numCols, std::invalid_argument,
304 os.str() <<
"Index range exceeds number of vectors " << numCols
305 <<
" in the input multivector.");
306 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
307 os.str() <<
"Should never get here!");
309 Teuchos::RCP<const MV> X_view = mv.subView (index);
310 Teuchos::rcp_const_cast<MV> (X_view)->setCopyOrView (Teuchos::View);
315 return static_cast<ptrdiff_t
> (mv.getGlobalLength ());
319 return static_cast<int> (mv.getNumVectors ());
323 return mv.isConstantStride ();
328 const Tpetra::MultiVector<Scalar,LO,GO,Node>& A,
329 const Teuchos::SerialDenseMatrix<int,dot_type>& B,
331 Tpetra::MultiVector<Scalar,LO,GO,Node>& C)
337 const int numRowsB = B.numRows ();
338 const int numColsB = B.numCols ();
339 const int strideB = B.stride ();
340 if (numRowsB == 1 && numColsB == 1) {
341 C.update (alpha*B(0,0), A, beta);
348 if (A.isConstantStride() ==
false) Atmp = rcp (
new MV (A, Teuchos::Copy));
349 else Atmp = rcp(&A,
false);
351 if (C.isConstantStride() ==
false) Ctmp = rcp (
new MV (C, Teuchos::Copy));
352 else Ctmp = rcp(&C,
false);
355 typedef Tpetra::MultiVector<dot_type,LO,GO,Node> FMV;
356 typedef Tpetra::MultiVector<const dot_type,LO,GO,Node> CFMV;
357 typedef typename FMV::dual_view_type::t_dev flat_view_type;
358 typedef typename CFMV::dual_view_type::t_dev const_flat_view_type;
360 const_flat_view_type flat_A_view = Atmp->template getLocalView<execution_space>(Tpetra::Access::ReadOnly);
361 flat_view_type flat_C_view = Ctmp->template getLocalView<execution_space>(Tpetra::Access::ReadWrite);
364 typedef Kokkos::View<dot_type**, Kokkos::LayoutLeft, Kokkos::HostSpace> b_host_view_type;
365 b_host_view_type B_view_host_input( B.values(), strideB, numColsB);
366 auto B_view_host = Kokkos::subview( B_view_host_input,
367 Kokkos::pair<int,int>(0, numRowsB),
368 Kokkos::pair<int,int>(0, numColsB));
372 typedef Kokkos::View<dot_type**, Kokkos::LayoutLeft, execution_space> b_view_type;
373 typedef Kokkos::View<dot_type*, Kokkos::LayoutLeft, execution_space> b_1d_view_type;
374 b_1d_view_type B_1d_view_dev(Kokkos::ViewAllocateWithoutInitializing(
"B"), numRowsB*numColsB);
375 b_view_type B_view_dev( B_1d_view_dev.data(), numRowsB, numColsB);
379 for (
int j=0;
j<numColsB; ++
j)
380 Kokkos::deep_copy(Kokkos::subview(B_view_dev,Kokkos::ALL,
j),
381 Kokkos::subview(B_view_host,Kokkos::ALL,
j));
385 const char ctransA =
'N', ctransB =
'N';
388 alpha, flat_A_view, B_view_dev, beta, flat_C_view);
391 if (C.isConstantStride() ==
false)
404 const Tpetra::MultiVector<Scalar,LO,GO,Node>& A,
406 const Tpetra::MultiVector<Scalar,LO,GO,Node>& B,
407 Tpetra::MultiVector<Scalar,LO,GO,Node>& mv)
409 mv.update (alpha, A, beta, B, Teuchos::ScalarTraits<Scalar>::zero ());
413 MvScale (Tpetra::MultiVector<Scalar,LO,GO,Node>& mv,
420 MvScale (Tpetra::MultiVector<Scalar,LO,GO,Node>& mv,
421 const std::vector<BaseScalar>& alphas)
423 std::vector<Scalar> alphas_mp(alphas.size());
424 const size_t sz = alphas.size();
425 for (
size_t i=0; i<sz; ++i)
426 alphas_mp[i] = alphas[i];
427 mv.scale (alphas_mp);
431 MvScale (Tpetra::MultiVector<Scalar,LO,GO,Node>& mv,
432 const std::vector<Scalar>& alphas)
439 const Tpetra::MultiVector<Scalar,LO,GO,Node>& A,
440 const Tpetra::MultiVector<Scalar,LO,GO,Node>& B,
441 Teuchos::SerialDenseMatrix<int,dot_type>& C)
446 using Teuchos::REDUCE_SUM;
447 using Teuchos::reduceAll;
450 const int numRowsC = C.numRows ();
451 const int numColsC = C.numCols ();
452 const int strideC = C.stride ();
453 if (numRowsC == 1 && numColsC == 1) {
454 if (alpha == Teuchos::ScalarTraits<Scalar>::zero ()) {
459 A.dot (B, Teuchos::ArrayView<dot_type> (C.values (), 1));
460 if (alpha != Teuchos::ScalarTraits<Scalar>::one ()) {
467 RCP<const MV> Atmp, Btmp;
468 if (A.isConstantStride() ==
false) Atmp = rcp (
new MV (A, Teuchos::Copy));
469 else Atmp = rcp(&A,
false);
471 if (B.isConstantStride() ==
false) Btmp = rcp (
new MV (B, Teuchos::Copy));
472 else Btmp = rcp(&B,
false);
475 typedef Tpetra::MultiVector<const dot_type,LO,GO,Node> FMV;
476 typedef typename FMV::dual_view_type::t_dev flat_view_type;
478 flat_view_type flat_A_view = Atmp->template getLocalView<execution_space>(Tpetra::Access::ReadOnly);
479 flat_view_type flat_B_view = Btmp->template getLocalView<execution_space>(Tpetra::Access::ReadOnly);
482 typedef Kokkos::View<dot_type**, Kokkos::LayoutLeft, Kokkos::HostSpace> c_host_view_type;
483 c_host_view_type C_view_host_input( C.values(), strideC, numColsC);
484 auto C_view_host = Kokkos::subview(C_view_host_input,
485 Kokkos::pair<int,int>(0, numRowsC),
486 Kokkos::pair<int,int>(0, numColsC));
490 typedef Kokkos::View<dot_type**, Kokkos::LayoutLeft, execution_space> c_view_type;
491 typedef Kokkos::View<dot_type*, Kokkos::LayoutLeft, execution_space> c_1d_view_type;
492 c_1d_view_type C_1d_view_dev(
"C", numRowsC*numColsC);
493 c_view_type C_view_dev( C_1d_view_dev.data(), numRowsC, numColsC);
497 const char ctransA =
'C', ctransB =
'N';
500 alpha, flat_A_view, flat_B_view,
501 Kokkos::Details::ArithTraits<dot_type>::zero(),
505 RCP<const Comm<int> > pcomm = A.getMap()->getComm ();
506 if (pcomm->getSize () == 1) {
510 for (
int j=0;
j<numColsC; ++
j)
511 Kokkos::deep_copy(Kokkos::subview(C_view_host,Kokkos::ALL,
j),
512 Kokkos::subview(C_view_dev,Kokkos::ALL,
j));
515 typedef Kokkos::View<dot_type*, Kokkos::LayoutLeft, Kokkos::HostSpace> c_1d_host_view_type;
516 c_1d_host_view_type C_1d_view_tmp(Kokkos::ViewAllocateWithoutInitializing(
"C_tmp"), strideC*numColsC);
517 c_host_view_type C_view_tmp( C_1d_view_tmp.data(),
519 for (
int j=0;
j<numColsC; ++
j)
520 Kokkos::deep_copy(Kokkos::subview(C_view_tmp,
521 Kokkos::pair<int,int>(0, numRowsC),
523 Kokkos::subview(C_view_dev, Kokkos::ALL,
j));
524 reduceAll<int> (*pcomm, REDUCE_SUM, strideC*numColsC,
532 MvDot (
const Tpetra::MultiVector<Scalar,LO,GO,Node>& A,
533 const Tpetra::MultiVector<Scalar,LO,GO,Node>& B,
534 std::vector<dot_type>& dots)
536 const size_t numVecs = A.getNumVectors ();
538 TEUCHOS_TEST_FOR_EXCEPTION(
539 numVecs != B.getNumVectors (), std::invalid_argument,
540 "Belos::MultiVecTraits<Scalar,Tpetra::MultiVector>::MvDot(A,B,dots): "
541 "A and B must have the same number of columns. "
542 "A has " << numVecs <<
" column(s), "
543 "but B has " << B.getNumVectors () <<
" column(s).");
544#ifdef HAVE_TPETRA_DEBUG
545 TEUCHOS_TEST_FOR_EXCEPTION(
546 dots.size() < numVecs, std::invalid_argument,
547 "Belos::MultiVecTraits<Scalar,Tpetra::MultiVector>::MvDot(A,B,dots): "
548 "The output array 'dots' must have room for all dot products. "
549 "A and B each have " << numVecs <<
" column(s), "
550 "but 'dots' only has " << dots.size() <<
" entry(/ies).");
553 Teuchos::ArrayView<dot_type> av (dots);
554 A.dot (B, av (0, numVecs));
559 MvNorm (
const Tpetra::MultiVector<Scalar,LO,GO,Node>& mv,
560 std::vector<mag_type> &normvec,
561 NormType type=TwoNorm)
564#ifdef HAVE_TPETRA_DEBUG
565 typedef std::vector<int>::size_type size_type;
566 TEUCHOS_TEST_FOR_EXCEPTION(
567 normvec.size () <
static_cast<size_type
> (mv.getNumVectors ()),
568 std::invalid_argument,
569 "Belos::MultiVecTraits::MvNorm(mv,normvec): The normvec output "
570 "argument must have at least as many entries as the number of vectors "
571 "(columns) in the MultiVector mv. normvec.size() = " << normvec.size ()
572 <<
" < mv.getNumVectors() = " << mv.getNumVectors () <<
".");
575 Teuchos::ArrayView<mag_type> av(normvec);
578 mv.norm1(av(0,mv.getNumVectors()));
581 mv.norm2(av(0,mv.getNumVectors()));
584 mv.normInf(av(0,mv.getNumVectors()));
590 TEUCHOS_TEST_FOR_EXCEPTION(
591 true, std::logic_error,
592 "Belos::MultiVecTraits::MvNorm: Invalid NormType value " << type
593 <<
". Valid values are OneNorm=" << OneNorm <<
", TwoNorm="
594 << TwoNorm <<
", and InfNorm=" << InfNorm <<
". If you are a Belos "
595 "user and have not modified Belos in any way, and you get this "
596 "message, then this is probably a bug in the Belos solver you were "
597 "using. Please report this to the Belos developers.");
604 using Teuchos::Range1D;
606 const size_t inNumVecs = A.getNumVectors ();
607#ifdef HAVE_TPETRA_DEBUG
608 TEUCHOS_TEST_FOR_EXCEPTION(
609 inNumVecs <
static_cast<size_t> (index.size ()), std::invalid_argument,
610 "Belos::MultiVecTraits::SetBlock(A,index,mv): 'index' argument must "
611 "have no more entries as the number of columns in the input MultiVector"
612 " A. A.getNumVectors() = " << inNumVecs <<
" < index.size () = "
613 << index.size () <<
".");
615 RCP<MV> mvsub = CloneViewNonConst (mv, index);
616 if (inNumVecs >
static_cast<size_t> (index.size ())) {
617 RCP<const MV> Asub = A.subView (Range1D (0, index.size () - 1));
618 ::Tpetra::deep_copy (*mvsub, *Asub);
620 ::Tpetra::deep_copy (*mvsub, A);
634 const size_t maxInt =
635 static_cast<size_t> (Teuchos::OrdinalTraits<int>::max ());
636 const bool overflow =
637 maxInt < A.getNumVectors () && maxInt < mv.getNumVectors ();
639 std::ostringstream os;
640 os <<
"Belos::MultiVecTraits::SetBlock(A, index=[" << index.lbound ()
641 <<
", " << index.ubound () <<
"], mv): ";
642 TEUCHOS_TEST_FOR_EXCEPTION(
643 maxInt < A.getNumVectors (), std::range_error, os.str () <<
"Number "
644 "of columns (size_t) in the input MultiVector 'A' overflows int.");
645 TEUCHOS_TEST_FOR_EXCEPTION(
646 maxInt < mv.getNumVectors (), std::range_error, os.str () <<
"Number "
647 "of columns (size_t) in the output MultiVector 'mv' overflows int.");
650 const int numColsA =
static_cast<int> (A.getNumVectors ());
651 const int numColsMv =
static_cast<int> (mv.getNumVectors ());
653 const bool validIndex =
654 index.lbound () >= 0 && index.ubound () < numColsMv;
656 const bool validSource = index.size () <= numColsA;
658 if (! validIndex || ! validSource) {
659 std::ostringstream os;
660 os <<
"Belos::MultiVecTraits::SetBlock(A, index=[" << index.lbound ()
661 <<
", " << index.ubound () <<
"], mv): ";
662 TEUCHOS_TEST_FOR_EXCEPTION(
663 index.lbound() < 0, std::invalid_argument,
664 os.str() <<
"Range lower bound must be nonnegative.");
665 TEUCHOS_TEST_FOR_EXCEPTION(
666 index.ubound() >= numColsMv, std::invalid_argument,
667 os.str() <<
"Range upper bound must be less than the number of "
668 "columns " << numColsA <<
" in the 'mv' output argument.");
669 TEUCHOS_TEST_FOR_EXCEPTION(
670 index.size() > numColsA, std::invalid_argument,
671 os.str() <<
"Range must have no more elements than the number of "
672 "columns " << numColsA <<
" in the 'A' input argument.");
673 TEUCHOS_TEST_FOR_EXCEPTION(
674 true, std::logic_error,
"Should never get here!");
680 Teuchos::RCP<MV> mv_view;
681 if (index.lbound () == 0 && index.ubound () + 1 == numColsMv) {
682 mv_view = Teuchos::rcpFromRef (mv);
684 mv_view = CloneViewNonConst (mv, index);
690 Teuchos::RCP<const MV> A_view;
691 if (index.size () == numColsA) {
692 A_view = Teuchos::rcpFromRef (A);
694 A_view = CloneView (A, Teuchos::Range1D (0, index.size () - 1));
697 ::Tpetra::deep_copy (*mv_view, *A_view);
702 const char errPrefix[] =
"Belos::MultiVecTraits::Assign(A, mv): ";
711 const size_t maxInt =
712 static_cast<size_t> (Teuchos::OrdinalTraits<int>::max ());
713 const bool overflow =
714 maxInt < A.getNumVectors () && maxInt < mv.getNumVectors ();
716 TEUCHOS_TEST_FOR_EXCEPTION(
717 maxInt < A.getNumVectors(), std::range_error,
718 errPrefix <<
"Number of columns in the input multivector 'A' "
719 "(a size_t) overflows int.");
720 TEUCHOS_TEST_FOR_EXCEPTION(
721 maxInt < mv.getNumVectors(), std::range_error,
722 errPrefix <<
"Number of columns in the output multivector 'mv' "
723 "(a size_t) overflows int.");
724 TEUCHOS_TEST_FOR_EXCEPTION(
725 true, std::logic_error,
"Should never get here!");
728 const int numColsA =
static_cast<int> (A.getNumVectors ());
729 const int numColsMv =
static_cast<int> (mv.getNumVectors ());
730 if (numColsA > numColsMv) {
731 TEUCHOS_TEST_FOR_EXCEPTION(
732 numColsA > numColsMv, std::invalid_argument,
733 errPrefix <<
"Input multivector 'A' has " << numColsA <<
" columns, "
734 "but output multivector 'mv' has only " << numColsMv <<
" columns.");
735 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
"Should never get here!");
737 if (numColsA == numColsMv) {
738 ::Tpetra::deep_copy (mv, A);
740 Teuchos::RCP<MV> mv_view =
741 CloneViewNonConst (mv, Teuchos::Range1D (0, numColsA - 1));
742 ::Tpetra::deep_copy (*mv_view, A);
754 mv.putScalar (alpha);
758 Teuchos::FancyOStream fos (Teuchos::rcpFromRef (os));
759 mv.describe (fos, Teuchos::VERB_EXTREME);
762#ifdef HAVE_BELOS_TSQR
766 typedef Tpetra::TsqrAdaptor< Tpetra::MultiVector< Scalar, LO, GO, Node > > tsqr_adaptor_type;
777 template <
class Storage,
class LO,
class GO,
class Node>
778 class OperatorTraits <typename
Storage::value_type,
779 Tpetra::MultiVector<Sacado::UQ::PCE<Storage>,
781 Tpetra::Operator<Sacado::UQ::PCE<Storage>,
787 Apply (
const Tpetra::Operator<Scalar,LO,GO,Node>& Op,
788 const Tpetra::MultiVector<Scalar,LO,GO,Node>& X,
789 Tpetra::MultiVector<Scalar,LO,GO,Node>& Y,
790 ETrans trans=NOTRANS)
794 Op.apply(X,Y,Teuchos::NO_TRANS);
797 Op.apply(X,Y,Teuchos::TRANS);
800 Op.apply(X,Y,Teuchos::CONJ_TRANS);
803 const std::string scalarName = Teuchos::TypeNameTraits<Scalar>::name();
804 const std::string loName = Teuchos::TypeNameTraits<LO>::name();
805 const std::string goName = Teuchos::TypeNameTraits<GO>::name();
806 const std::string nodeName = Teuchos::TypeNameTraits<Node>::name();
807 const std::string otName =
"Belos::OperatorTraits<" + scalarName
808 +
"," + loName +
"," + goName +
"," + nodeName +
">";
809 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error, otName <<
": Should never "
810 "get here; fell through a switch statement. "
811 "Please report this bug to the Belos developers.");
818 return Op.hasTransposeApply ();
Kokkos::DefaultExecutionSpace execution_space
static void SetBlock(const MV &A, const Teuchos::Range1D &index, 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.
Tpetra::MultiVector< Scalar, LO, GO, Node >::dot_type dot_type
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< const MV > CloneView(const MV &mv, const std::vector< int > &index)
static Teuchos::RCP< MV > CloneCopy(const MV &X)
Create and return a deep copy of X.
static Teuchos::RCP< MV > CloneViewNonConst(MV &mv, const std::vector< int > &index)
static Teuchos::RCP< MV > CloneViewNonConst(MV &mv, const Teuchos::Range1D &index)
static int GetNumberVecs(const MV &mv)
static void Assign(const MV &A, MV &mv)
static Teuchos::RCP< MV > Clone(const MV &X, const int numVecs)
Create a new MultiVector with numVecs columns.
static ptrdiff_t GetGlobalLength(const MV &mv)
static Teuchos::RCP< const MV > CloneView(const MV &mv, const Teuchos::Range1D &index)
Tpetra::MultiVector< Scalar, LO, GO, Node > MV
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 void SetBlock(const MV &A, const std::vector< int > &index, MV &mv)
static void MvScale(Tpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< Scalar > &alphas)
Storage::ordinal_type s_ordinal
static void MvPrint(const MV &mv, std::ostream &os)
Storage::value_type BaseScalar
static void MvRandom(MV &mv)
static void MvScale(Tpetra::MultiVector< Scalar, LO, GO, Node > &mv, const std::vector< BaseScalar > &alphas)
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 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)
Sacado::UQ::PCE< Storage > Scalar
static void MvInit(MV &mv, const Scalar alpha=Teuchos::ScalarTraits< Scalar >::zero())
static bool HasConstantStride(const MV &mv)
static void MvScale(Tpetra::MultiVector< Scalar, LO, GO, Node > &mv, const Scalar &alpha)
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].
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 bool HasApplyTranspose(const Tpetra::Operator< Scalar, LO, GO, Node > &Op)
Sacado::UQ::PCE< Storage > Scalar
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)
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)