42 #ifndef TPETRA_BLOCKMULTIVECTOR_DEF_HPP
43 #define TPETRA_BLOCKMULTIVECTOR_DEF_HPP
47 #include "Teuchos_OrdinalTraits.hpp"
62 template<
class MultiVectorType>
63 struct RawHostPtrFromMultiVector {
64 typedef typename MultiVectorType::impl_scalar_type impl_scalar_type;
66 static impl_scalar_type* getRawPtr (MultiVectorType& X) {
72 auto X_view_host = X.template getLocalView<typename MultiVectorType::dual_view_type::t_host::device_type> ();
73 impl_scalar_type* X_raw = X_view_host.data ();
90 template<
class S,
class LO,
class GO,
class N>
94 return RawHostPtrFromMultiVector<MV>::getRawPtr (X);
101 template<
class Scalar,
class LO,
class GO,
class Node>
102 typename BlockMultiVector<Scalar, LO, GO, Node>::mv_type
109 template<
class Scalar,
class LO,
class GO,
class Node>
110 Teuchos::RCP<const BlockMultiVector<Scalar, LO, GO, Node> >
115 const BMV* src_bmv =
dynamic_cast<const BMV*
> (&src);
116 TEUCHOS_TEST_FOR_EXCEPTION(
117 src_bmv ==
nullptr, std::invalid_argument,
"Tpetra::"
118 "BlockMultiVector: The source object of an Import or Export to a "
119 "BlockMultiVector, must also be a BlockMultiVector.");
120 return Teuchos::rcp (src_bmv,
false);
123 template<
class Scalar,
class LO,
class GO,
class Node>
126 const Teuchos::DataAccess copyOrView) :
128 meshMap_ (in.meshMap_),
129 pointMap_ (in.pointMap_),
130 mv_ (in.mv_, copyOrView),
131 mvData_ (getRawHostPtrFromMultiVector (mv_)),
132 blockSize_ (in.blockSize_)
135 template<
class Scalar,
class LO,
class GO,
class Node>
142 pointMap_ (makePointMap (meshMap, blockSize)),
143 mv_ (Teuchos::rcpFromRef (pointMap_), numVecs),
144 mvData_ (getRawHostPtrFromMultiVector (mv_)),
145 blockSize_ (blockSize)
148 template<
class Scalar,
class LO,
class GO,
class Node>
156 pointMap_ (pointMap),
157 mv_ (Teuchos::rcpFromRef (pointMap_), numVecs),
158 mvData_ (getRawHostPtrFromMultiVector (mv_)),
159 blockSize_ (blockSize)
162 template<
class Scalar,
class LO,
class GO,
class Node>
166 const LO blockSize) :
170 blockSize_ (blockSize)
186 RCP<const mv_type> X_view_const;
189 Teuchos::Array<size_t> cols (0);
190 X_view_const = X_mv.
subView (cols ());
192 X_view_const = X_mv.
subView (Teuchos::Range1D (0, numCols-1));
194 TEUCHOS_TEST_FOR_EXCEPTION(
195 X_view_const.is_null (), std::logic_error,
"Tpetra::"
196 "BlockMultiVector constructor: X_mv.subView(...) returned null. This "
197 "should never happen. Please report this bug to the Tpetra developers.");
202 RCP<mv_type> X_view = Teuchos::rcp_const_cast<
mv_type> (X_view_const);
203 TEUCHOS_TEST_FOR_EXCEPTION(
204 X_view->getCopyOrView () != Teuchos::View, std::logic_error,
"Tpetra::"
205 "BlockMultiVector constructor: We just set a MultiVector "
206 "to have view semantics, but it claims that it doesn't have view "
207 "semantics. This should never happen. "
208 "Please report this bug to the Tpetra developers.");
213 Teuchos::RCP<const map_type> pointMap =
mv_.
getMap ();
214 if (! pointMap.is_null ()) {
215 pointMap_ = *pointMap;
217 mvData_ = getRawHostPtrFromMultiVector (
mv_);
220 template<
class Scalar,
class LO,
class GO,
class Node>
225 const size_t offset) :
227 meshMap_ (newMeshMap),
228 pointMap_ (newPointMap),
229 mv_ (X.mv_, newPointMap, offset * X.getBlockSize ()),
230 mvData_ (getRawHostPtrFromMultiVector (mv_)),
231 blockSize_ (X.getBlockSize ())
234 template<
class Scalar,
class LO,
class GO,
class Node>
238 const size_t offset) :
240 meshMap_ (newMeshMap),
241 pointMap_ (makePointMap (newMeshMap, X.getBlockSize ())),
242 mv_ (X.mv_, pointMap_, offset * X.getBlockSize ()),
243 mvData_ (getRawHostPtrFromMultiVector (mv_)),
244 blockSize_ (X.getBlockSize ())
247 template<
class Scalar,
class LO,
class GO,
class Node>
255 template<
class Scalar,
class LO,
class GO,
class Node>
261 typedef typename Teuchos::ArrayView<const GO>::size_type size_type;
263 const GST gblNumMeshMapInds =
265 const size_t lclNumMeshMapIndices =
267 const GST gblNumPointMapInds =
268 gblNumMeshMapInds *
static_cast<GST
> (blockSize);
269 const size_t lclNumPointMapInds =
270 lclNumMeshMapIndices *
static_cast<size_t> (blockSize);
274 return map_type (gblNumPointMapInds, lclNumPointMapInds, indexBase,
282 const size_type lclNumMeshGblInds = lclMeshGblInds.size ();
283 Teuchos::Array<GO> lclPointGblInds (lclNumPointMapInds);
284 for (size_type g = 0; g < lclNumMeshGblInds; ++g) {
285 const GO meshGid = lclMeshGblInds[g];
286 const GO pointGidStart = indexBase +
287 (meshGid - indexBase) * static_cast<GO> (blockSize);
288 const size_type offset = g *
static_cast<size_type
> (blockSize);
289 for (LO k = 0; k < blockSize; ++k) {
290 const GO pointGid = pointGidStart +
static_cast<GO
> (k);
291 lclPointGblInds[offset +
static_cast<size_type
> (k)] = pointGid;
294 return map_type (gblNumPointMapInds, lclPointGblInds (), indexBase,
300 template<
class Scalar,
class LO,
class GO,
class Node>
305 const Scalar vals[])
const
307 auto X_dst = getLocalBlock (localRowIndex, colIndex);
308 typename const_little_vec_type::HostMirror::const_type X_src (reinterpret_cast<const impl_scalar_type*> (vals),
314 template<
class Scalar,
class LO,
class GO,
class Node>
319 const Scalar vals[])
const
321 if (! meshMap_.isNodeLocalElement (localRowIndex)) {
324 replaceLocalValuesImpl (localRowIndex, colIndex, vals);
329 template<
class Scalar,
class LO,
class GO,
class Node>
334 const Scalar vals[])
const
336 const LO localRowIndex = meshMap_.getLocalElement (globalRowIndex);
337 if (localRowIndex == Teuchos::OrdinalTraits<LO>::invalid ()) {
340 replaceLocalValuesImpl (localRowIndex, colIndex, vals);
345 template<
class Scalar,
class LO,
class GO,
class Node>
350 const Scalar vals[])
const
352 auto X_dst = getLocalBlock (localRowIndex, colIndex);
353 typename const_little_vec_type::HostMirror::const_type X_src (reinterpret_cast<const impl_scalar_type*> (vals),
355 AXPY (static_cast<impl_scalar_type> (STS::one ()), X_src, X_dst);
358 template<
class Scalar,
class LO,
class GO,
class Node>
363 const Scalar vals[])
const
365 if (! meshMap_.isNodeLocalElement (localRowIndex)) {
368 sumIntoLocalValuesImpl (localRowIndex, colIndex, vals);
373 template<
class Scalar,
class LO,
class GO,
class Node>
378 const Scalar vals[])
const
380 const LO localRowIndex = meshMap_.getLocalElement (globalRowIndex);
381 if (localRowIndex == Teuchos::OrdinalTraits<LO>::invalid ()) {
384 sumIntoLocalValuesImpl (localRowIndex, colIndex, vals);
389 template<
class Scalar,
class LO,
class GO,
class Node>
394 if (! meshMap_.isNodeLocalElement (localRowIndex)) {
397 auto X_ij = getLocalBlock (localRowIndex, colIndex);
398 vals =
reinterpret_cast<Scalar*
> (X_ij.data ());
403 template<
class Scalar,
class LO,
class GO,
class Node>
408 const LO localRowIndex = meshMap_.getLocalElement (globalRowIndex);
409 if (localRowIndex == Teuchos::OrdinalTraits<LO>::invalid ()) {
412 auto X_ij = getLocalBlock (localRowIndex, colIndex);
413 vals =
reinterpret_cast<Scalar*
> (X_ij.data ());
418 template<
class Scalar,
class LO,
class GO,
class Node>
422 const LO colIndex)
const
436 if (! isValidLocalMeshIndex (localRowIndex)) {
437 return typename little_vec_type::HostMirror ();
439 const size_t blockSize = getBlockSize ();
440 const size_t offset = colIndex * this->getStrideY () +
441 localRowIndex * blockSize;
443 return typename little_vec_type::HostMirror (blockRaw, blockSize);
447 template<
class Scalar,
class LO,
class GO,
class Node>
448 Teuchos::RCP<const typename BlockMultiVector<Scalar, LO, GO, Node>::mv_type>
453 using Teuchos::rcpFromRef;
460 const this_type* srcBlkVec =
dynamic_cast<const this_type*
> (&src);
461 if (srcBlkVec ==
nullptr) {
462 const mv_type* srcMultiVec =
dynamic_cast<const mv_type*
> (&src);
463 if (srcMultiVec ==
nullptr) {
467 return rcp (
new mv_type ());
469 return rcp (srcMultiVec,
false);
472 return rcpFromRef (srcBlkVec->mv_);
476 template<
class Scalar,
class LO,
class GO,
class Node>
480 return ! getMultiVectorFromSrcDistObject (src).is_null ();
483 template<
class Scalar,
class LO,
class GO,
class Node>
485 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
487 #else // TPETRA_ENABLE_DEPRECATED_CODE
489 #endif // TPETRA_ENABLE_DEPRECATED_CODE
491 const size_t numSameIDs,
492 const Kokkos::DualView<
const local_ordinal_type*,
493 buffer_device_type>& permuteToLIDs,
494 const Kokkos::DualView<
const local_ordinal_type*,
495 buffer_device_type>& permuteFromLIDs)
497 TEUCHOS_TEST_FOR_EXCEPTION
498 (
true, std::logic_error,
499 "Tpetra::BlockMultiVector::copyAndPermute: Do NOT use this "
500 "instead, create a point importer using makePointMap function.");
503 template<
class Scalar,
class LO,
class GO,
class Node>
504 void BlockMultiVector<Scalar, LO, GO, Node>::
505 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
507 #else // TPETRA_ENABLE_DEPRECATED_CODE
509 #endif // TPETRA_ENABLE_DEPRECATED_CODE
510 (
const SrcDistObject& src,
511 const Kokkos::DualView<
const local_ordinal_type*,
512 buffer_device_type>& exportLIDs,
513 Kokkos::DualView<packet_type*,
514 buffer_device_type>& exports,
515 Kokkos::DualView<
size_t*,
516 buffer_device_type> numPacketsPerLID,
517 size_t& constantNumPackets,
520 TEUCHOS_TEST_FOR_EXCEPTION
521 (
true, std::logic_error,
522 "Tpetra::BlockMultiVector::copyAndPermute: Do NOT use this; "
523 "instead, create a point importer using makePointMap function.");
526 template<
class Scalar,
class LO,
class GO,
class Node>
527 void BlockMultiVector<Scalar, LO, GO, Node>::
528 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
530 #else // TPETRA_ENABLE_DEPRECATED_CODE
532 #endif // TPETRA_ENABLE_DEPRECATED_CODE
534 buffer_device_type>& importLIDs,
535 Kokkos::DualView<packet_type*,
536 buffer_device_type> imports,
537 Kokkos::DualView<
size_t*,
538 buffer_device_type> numPacketsPerLID,
539 const size_t constantNumPackets,
543 TEUCHOS_TEST_FOR_EXCEPTION
544 (
true, std::logic_error,
545 "Tpetra::BlockMultiVector::copyAndPermute: Do NOT use this; "
546 "instead, create a point importer using makePointMap function.");
549 template<
class Scalar,
class LO,
class GO,
class Node>
553 return meshLocalIndex != Teuchos::OrdinalTraits<LO>::invalid () &&
554 meshMap_.isNodeLocalElement (meshLocalIndex);
557 template<
class Scalar,
class LO,
class GO,
class Node>
564 template<
class Scalar,
class LO,
class GO,
class Node>
571 template<
class Scalar,
class LO,
class GO,
class Node>
577 mv_.update (alpha, X.
mv_, beta);
582 template <
typename Scalar,
typename ViewY,
typename ViewD,
typename ViewX>
583 struct BlockWiseMultiply {
584 typedef typename ViewD::size_type Size;
587 typedef typename ViewD::device_type Device;
588 typedef typename ViewD::non_const_value_type ImplScalar;
589 typedef Kokkos::MemoryTraits<Kokkos::Unmanaged> Unmanaged;
591 template <
typename View>
592 using UnmanagedView = Kokkos::View<
typename View::data_type,
typename View::array_layout,
593 typename View::device_type, Unmanaged>;
594 typedef UnmanagedView<ViewY> UnMViewY;
595 typedef UnmanagedView<ViewD> UnMViewD;
596 typedef UnmanagedView<ViewX> UnMViewX;
598 const Size block_size_;
605 BlockWiseMultiply (
const Size block_size,
const Scalar& alpha,
606 const ViewY& Y,
const ViewD& D,
const ViewX& X)
607 : block_size_(block_size), alpha_(alpha), Y_(Y), D_(D), X_(X)
610 KOKKOS_INLINE_FUNCTION
611 void operator() (
const Size k)
const {
612 const auto zero = Kokkos::Details::ArithTraits<Scalar>::zero();
613 auto D_curBlk = Kokkos::subview(D_, k, Kokkos::ALL (), Kokkos::ALL ());
614 const auto num_vecs = X_.extent(1);
615 for (Size i = 0; i < num_vecs; ++i) {
616 Kokkos::pair<Size, Size> kslice(k*block_size_, (k+1)*block_size_);
617 auto X_curBlk = Kokkos::subview(X_, kslice, i);
618 auto Y_curBlk = Kokkos::subview(Y_, kslice, i);
627 template <
typename Scalar,
typename ViewY,
typename ViewD,
typename ViewX>
628 inline BlockWiseMultiply<Scalar, ViewY, ViewD, ViewX>
629 createBlockWiseMultiply (
const int block_size,
const Scalar& alpha,
630 const ViewY& Y,
const ViewD& D,
const ViewX& X) {
631 return BlockWiseMultiply<Scalar, ViewY, ViewD, ViewX>(block_size, alpha, Y, D, X);
634 template <
typename ViewY,
638 typename LO =
typename ViewY::size_type>
639 class BlockJacobiUpdate {
641 typedef typename ViewD::device_type Device;
642 typedef typename ViewD::non_const_value_type ImplScalar;
643 typedef Kokkos::MemoryTraits<Kokkos::Unmanaged> Unmanaged;
645 template <
typename ViewType>
646 using UnmanagedView = Kokkos::View<
typename ViewType::data_type,
647 typename ViewType::array_layout,
648 typename ViewType::device_type,
650 typedef UnmanagedView<ViewY> UnMViewY;
651 typedef UnmanagedView<ViewD> UnMViewD;
652 typedef UnmanagedView<ViewZ> UnMViewZ;
662 BlockJacobiUpdate (
const ViewY& Y,
666 const Scalar& beta) :
667 blockSize_ (D.extent (1)),
675 static_assert (static_cast<int> (ViewY::rank) == 1,
676 "Y must have rank 1.");
677 static_assert (static_cast<int> (ViewD::rank) == 3,
"D must have rank 3.");
678 static_assert (static_cast<int> (ViewZ::rank) == 1,
679 "Z must have rank 1.");
685 KOKKOS_INLINE_FUNCTION
void
686 operator() (
const LO& k)
const
689 using Kokkos::subview;
690 typedef Kokkos::pair<LO, LO> range_type;
691 typedef Kokkos::Details::ArithTraits<Scalar> KAT;
695 auto D_curBlk = subview (D_, k, ALL (), ALL ());
696 const range_type kslice (k*blockSize_, (k+1)*blockSize_);
700 auto Z_curBlk = subview (Z_, kslice);
701 auto Y_curBlk = subview (Y_, kslice);
703 if (beta_ == KAT::zero ()) {
706 else if (beta_ != KAT::one ()) {
713 template<
class ViewY,
717 class LO =
typename ViewD::size_type>
719 blockJacobiUpdate (
const ViewY& Y,
725 static_assert (Kokkos::Impl::is_view<ViewY>::value,
"Y must be a Kokkos::View.");
726 static_assert (Kokkos::Impl::is_view<ViewD>::value,
"D must be a Kokkos::View.");
727 static_assert (Kokkos::Impl::is_view<ViewZ>::value,
"Z must be a Kokkos::View.");
728 static_assert (static_cast<int> (ViewY::rank) == static_cast<int> (ViewZ::rank),
729 "Y and Z must have the same rank.");
730 static_assert (static_cast<int> (ViewD::rank) == 3,
"D must have rank 3.");
732 const auto lclNumMeshRows = D.extent (0);
734 #ifdef HAVE_TPETRA_DEBUG
738 const auto blkSize = D.extent (1);
739 const auto lclNumPtRows = lclNumMeshRows * blkSize;
740 TEUCHOS_TEST_FOR_EXCEPTION
741 (Y.extent (0) != lclNumPtRows, std::invalid_argument,
742 "blockJacobiUpdate: Y.extent(0) = " << Y.extent (0) <<
" != "
743 "D.extent(0)*D.extent(1) = " << lclNumMeshRows <<
" * " << blkSize
744 <<
" = " << lclNumPtRows <<
".");
745 TEUCHOS_TEST_FOR_EXCEPTION
746 (Y.extent (0) != Z.extent (0), std::invalid_argument,
747 "blockJacobiUpdate: Y.extent(0) = " << Y.extent (0) <<
" != "
748 "Z.extent(0) = " << Z.extent (0) <<
".");
749 TEUCHOS_TEST_FOR_EXCEPTION
750 (Y.extent (1) != Z.extent (1), std::invalid_argument,
751 "blockJacobiUpdate: Y.extent(1) = " << Y.extent (1) <<
" != "
752 "Z.extent(1) = " << Z.extent (1) <<
".");
753 #endif // HAVE_TPETRA_DEBUG
755 BlockJacobiUpdate<ViewY, Scalar, ViewD, ViewZ, LO> functor (Y, alpha, D, Z, beta);
756 typedef Kokkos::RangePolicy<typename ViewY::execution_space, LO> range_type;
758 range_type range (0, static_cast<LO> (lclNumMeshRows));
759 Kokkos::parallel_for (range, functor);
764 template<
class Scalar,
class LO,
class GO,
class Node>
773 typedef typename device_type::memory_space memory_space;
774 const LO lclNumMeshRows = meshMap_.getNodeNumElements ();
776 if (alpha == STS::zero ()) {
777 this->putScalar (STS::zero ());
780 const LO blockSize = this->getBlockSize ();
782 auto X_lcl = X.
mv_.template getLocalView<memory_space> ();
783 auto Y_lcl = this->mv_.template getLocalView<memory_space> ();
784 auto bwm = Impl::createBlockWiseMultiply (blockSize, alphaImpl, Y_lcl, D, X_lcl);
791 Kokkos::RangePolicy<execution_space, LO> range (0, lclNumMeshRows);
792 Kokkos::parallel_for (range, bwm);
796 template<
class Scalar,
class LO,
class GO,
class Node>
806 using Kokkos::subview;
807 typedef typename device_type::memory_space memory_space;
810 const IST alphaImpl =
static_cast<IST
> (alpha);
811 const IST betaImpl =
static_cast<IST
> (beta);
812 const LO numVecs = mv_.getNumVectors ();
814 auto X_lcl = X.
mv_.template getLocalView<memory_space> ();
815 auto Y_lcl = this->mv_.template getLocalView<memory_space> ();
816 auto Z_lcl = Z.
mv_.template getLocalView<memory_space> ();
818 if (alpha == STS::zero ()) {
822 Z.
update (STS::one (), X, -STS::one ());
823 for (LO j = 0; j < numVecs; ++j) {
824 auto X_lcl_j = subview (X_lcl, ALL (), j);
825 auto Y_lcl_j = subview (Y_lcl, ALL (), j);
826 auto Z_lcl_j = subview (Z_lcl, ALL (), j);
827 Impl::blockJacobiUpdate (Y_lcl_j, alphaImpl, D, Z_lcl_j, betaImpl);
839 #define TPETRA_BLOCKMULTIVECTOR_INSTANT(S,LO,GO,NODE) \
840 template class BlockMultiVector< S, LO, GO, NODE >;
842 #endif // TPETRA_BLOCKMULTIVECTOR_DEF_HPP
Teuchos::RCP< const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > subView(const Teuchos::Range1D &colRng) const
Return a const MultiVector with const views of selected columns.
void update(const Scalar &alpha, const BlockMultiVector< Scalar, LO, GO, Node > &X, const Scalar &beta)
Update: this = beta*this + alpha*X.
little_vec_type::HostMirror getLocalBlock(const LO localRowIndex, const LO colIndex) const
Get a host view of the degrees of freedom at the given mesh point.
bool replaceLocalValues(const LO localRowIndex, const LO colIndex, const Scalar vals[]) const
Replace all values at the given mesh point, using local row and column indices.
bool getLocalRowView(const LO localRowIndex, const LO colIndex, Scalar *&vals) const
Get a writeable view of the entries at the given mesh point, using a local index. ...
bool isContiguous() const
True if this Map is distributed contiguously, else false.
void putScalar(const Scalar &val)
Fill all entries with the given value val.
size_t getNumVectors() const
Number of columns in the multivector.
void blockJacobiUpdate(const Scalar &alpha, const Kokkos::View< const impl_scalar_type ***, device_type, Kokkos::MemoryUnmanaged > &D, const BlockMultiVector< Scalar, LO, GO, Node > &X, BlockMultiVector< Scalar, LO, GO, Node > &Z, const Scalar &beta)
Block Jacobi update .
void blockWiseMultiply(const Scalar &alpha, const Kokkos::View< const impl_scalar_type ***, device_type, Kokkos::MemoryUnmanaged > &D, const BlockMultiVector< Scalar, LO, GO, Node > &X)
*this := alpha * D * X, where D is a block diagonal matrix.
One or more distributed dense vectors.
virtual bool checkSizes(const Tpetra::SrcDistObject &source)
Compare the source and target (this) objects for compatibility.
Linear algebra kernels for small dense matrices and vectors.
GlobalOrdinal getIndexBase() const
The index base for this Map.
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
Accessors for the Teuchos::Comm and Kokkos Node objects.
KOKKOS_INLINE_FUNCTION void FILL(const ViewType &x, const InputType &val)
Set every entry of x to val.
static map_type makePointMap(const map_type &meshMap, const LO blockSize)
Create and return the point Map corresponding to the given mesh Map and block size.
int local_ordinal_type
Default value of Scalar template parameter.
MultiVector for multiple degrees of freedom per mesh point.
KOKKOS_INLINE_FUNCTION void SCAL(const CoefficientType &alpha, const ViewType &x)
x := alpha*x, where x is either rank 1 (a vector) or rank 2 (a matrix).
size_t global_size_t
Global size_t object.
void deep_copy(MultiVector< DS, DL, DG, DN > &dst, const MultiVector< SS, SL, SG, SN > &src)
Copy the contents of the MultiVector src into dst.
BlockMultiVector()
Default constructor.
typename device_type::execution_space execution_space
The Kokkos execution space.
bool isValidLocalMeshIndex(const LO meshLocalIndex) const
True if and only if meshLocalIndex is a valid local index in the mesh Map.
bool sumIntoGlobalValues(const GO globalRowIndex, const LO colIndex, const Scalar vals[]) const
Sum into all values at the given mesh point, using a global index.
typename mv_type::device_type device_type
The Kokkos Device type.
CombineMode
Rule for combining data in an Import or Export.
mv_type mv_
The Tpetra::MultiVector used to represent the data.
Abstract base class for objects that can be the source of an Import or Export operation.
bool getGlobalRowView(const GO globalRowIndex, const LO colIndex, Scalar *&vals) const
Get a writeable view of the entries at the given mesh point, using a global index.
KOKKOS_INLINE_FUNCTION void GEMV(const CoeffType &alpha, const BlkType &A, const VecType1 &x, const VecType2 &y)
y := y + alpha * A * x (dense matrix-vector multiply)
bool replaceGlobalValues(const GO globalRowIndex, const LO colIndex, const Scalar vals[]) const
Replace all values at the given mesh point, using a global index.
Teuchos::ArrayView< const GlobalOrdinal > getNodeElementList() const
Return a NONOWNING view of the global indices owned by this process.
typename Kokkos::Details::ArithTraits< Scalar >::val_type impl_scalar_type
The type used internally in place of Scalar.
size_t getNodeNumElements() const
The number of elements belonging to the calling process.
mv_type getMultiVectorView() const
Get a Tpetra::MultiVector that views this BlockMultiVector's data.
virtual Teuchos::RCP< const map_type > getMap() const
The Map describing the parallel distribution of this object.
void scale(const Scalar &val)
Multiply all entries in place by the given value val.
Teuchos::DataAccess getCopyOrView() const
Get whether this has copy (copyOrView = Teuchos::Copy) or view (copyOrView = Teuchos::View) semantics...
typename mv_type::impl_scalar_type impl_scalar_type
The implementation type of entries in the object.
bool sumIntoLocalValues(const LO localRowIndex, const LO colIndex, const Scalar vals[]) const
Sum into all values at the given mesh point, using a local index.
KOKKOS_INLINE_FUNCTION void AXPY(const CoefficientType &alpha, const ViewType1 &x, const ViewType2 &y)
y := y + alpha * x (dense vector or matrix update)
global_size_t getGlobalNumElements() const
The number of elements in this Map.
Declaration of Tpetra::Details::Behavior, a class that describes Tpetra's behavior.