58 #ifndef TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_HPP
59 #define TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_HPP
61 #include "Kokkos_Core.hpp"
62 #include "Kokkos_ArithTraits.hpp"
67 namespace KokkosRefactor {
83 template<
class IntegerType,
84 const bool isSigned = std::numeric_limits<IntegerType>::is_signed>
86 static KOKKOS_INLINE_FUNCTION
bool
87 test (
const IntegerType x,
88 const IntegerType exclusiveUpperBound);
92 template<
class IntegerType>
94 static KOKKOS_INLINE_FUNCTION
bool
95 test (
const IntegerType x,
96 const IntegerType exclusiveUpperBound)
98 return x < static_cast<IntegerType> (0) || x >= exclusiveUpperBound;
103 template<
class IntegerType>
104 struct OutOfBounds<IntegerType, false> {
105 static KOKKOS_INLINE_FUNCTION
bool
106 test (
const IntegerType x,
107 const IntegerType exclusiveUpperBound)
109 return x >= exclusiveUpperBound;
115 template<
class IntegerType>
116 KOKKOS_INLINE_FUNCTION
bool
117 outOfBounds (
const IntegerType x,
const IntegerType exclusiveUpperBound)
127 template <
typename DstView,
typename SrcView,
typename IdxView,
128 typename Enabled =
void>
129 struct PackArraySingleColumn {
130 typedef typename DstView::execution_space execution_space;
131 typedef typename execution_space::size_type size_type;
138 PackArraySingleColumn (
const DstView& dst_,
142 dst(dst_), src(src_), idx(idx_), col(col_) {}
144 KOKKOS_INLINE_FUNCTION
void
145 operator() (
const size_type k)
const {
146 dst(k) = src(idx(k), col);
150 pack (
const DstView& dst,
154 const execution_space &space)
156 typedef Kokkos::RangePolicy<execution_space, size_type> range_type;
158 (
"Tpetra::MultiVector pack one col",
159 range_type (space, 0, idx.size ()),
160 PackArraySingleColumn (dst, src, idx, col));
164 template <
typename DstView,
167 typename SizeType =
typename DstView::execution_space::size_type,
168 typename Enabled =
void>
169 class PackArraySingleColumnWithBoundsCheck {
171 static_assert (Kokkos::is_view<DstView>::value,
172 "DstView must be a Kokkos::View.");
173 static_assert (Kokkos::is_view<SrcView>::value,
174 "SrcView must be a Kokkos::View.");
175 static_assert (Kokkos::is_view<IdxView>::value,
176 "IdxView must be a Kokkos::View.");
177 static_assert (static_cast<int> (DstView::rank) == 1,
178 "DstView must be a rank-1 Kokkos::View.");
179 static_assert (static_cast<int> (SrcView::rank) == 2,
180 "SrcView must be a rank-2 Kokkos::View.");
181 static_assert (static_cast<int> (IdxView::rank) == 1,
182 "IdxView must be a rank-1 Kokkos::View.");
183 static_assert (std::is_integral<SizeType>::value,
184 "SizeType must be a built-in integer type.");
186 using execution_space =
typename DstView::execution_space;
189 typedef SizeType size_type;
190 using value_type = size_t;
197 execution_space space;
200 PackArraySingleColumnWithBoundsCheck (
const DstView& dst_,
203 const size_type col_) :
204 dst (dst_), src (src_), idx (idx_), col (col_) {}
206 KOKKOS_INLINE_FUNCTION
void
207 operator() (
const size_type k, value_type& lclErrCount)
const {
208 using index_type =
typename IdxView::non_const_value_type;
210 const index_type lclRow = idx(k);
211 if (lclRow < static_cast<index_type> (0) ||
212 lclRow >= static_cast<index_type> (src.extent (0))) {
216 dst(k) = src(lclRow, col);
220 KOKKOS_INLINE_FUNCTION
221 void init (value_type& initialErrorCount)
const {
222 initialErrorCount = 0;
225 KOKKOS_INLINE_FUNCTION
void
226 join (value_type& dstErrorCount,
227 const value_type& srcErrorCount)
const
229 dstErrorCount += srcErrorCount;
233 pack (
const DstView& dst,
237 const execution_space &space)
239 typedef Kokkos::RangePolicy<execution_space, size_type> range_type;
240 typedef typename IdxView::non_const_value_type index_type;
242 size_t errorCount = 0;
243 Kokkos::parallel_reduce
244 (
"Tpetra::MultiVector pack one col debug only",
245 range_type (space, 0, idx.size ()),
246 PackArraySingleColumnWithBoundsCheck (dst, src, idx, col),
249 if (errorCount != 0) {
253 auto idx_h = Kokkos::create_mirror_view (idx);
258 std::vector<index_type> badIndices;
259 const size_type numInds = idx_h.extent (0);
260 for (size_type k = 0; k < numInds; ++k) {
261 if (idx_h(k) < static_cast<index_type> (0) ||
262 idx_h(k) >= static_cast<index_type> (src.extent (0))) {
263 badIndices.push_back (idx_h(k));
267 TEUCHOS_TEST_FOR_EXCEPTION
268 (errorCount != badIndices.size (), std::logic_error,
269 "PackArraySingleColumnWithBoundsCheck: errorCount = " << errorCount
270 <<
" != badIndices.size() = " << badIndices.size () <<
". This sho"
271 "uld never happen. Please report this to the Tpetra developers.");
273 std::ostringstream os;
274 os <<
"MultiVector single-column pack kernel had "
275 << badIndices.size () <<
" out-of bounds index/ices. "
277 for (
size_t k = 0; k < badIndices.size (); ++k) {
279 if (k + 1 < badIndices.size ()) {
284 throw std::runtime_error (os.str ());
290 template <
typename DstView,
typename SrcView,
typename IdxView>
292 pack_array_single_column (
const DstView& dst,
297 const typename DstView::execution_space &space)
299 static_assert (Kokkos::is_view<DstView>::value,
300 "DstView must be a Kokkos::View.");
301 static_assert (Kokkos::is_view<SrcView>::value,
302 "SrcView must be a Kokkos::View.");
303 static_assert (Kokkos::is_view<IdxView>::value,
304 "IdxView must be a Kokkos::View.");
305 static_assert (static_cast<int> (DstView::rank) == 1,
306 "DstView must be a rank-1 Kokkos::View.");
307 static_assert (static_cast<int> (SrcView::rank) == 2,
308 "SrcView must be a rank-2 Kokkos::View.");
309 static_assert (static_cast<int> (IdxView::rank) == 1,
310 "IdxView must be a rank-1 Kokkos::View.");
312 using execution_space =
typename DstView::execution_space;
314 static_assert (Kokkos::SpaceAccessibility<execution_space,
315 typename DstView::memory_space>::accessible,
316 "DstView not accessible from execution space");
317 static_assert (Kokkos::SpaceAccessibility<execution_space,
318 typename SrcView::memory_space>::accessible,
319 "SrcView not accessible from execution space");
320 static_assert (Kokkos::SpaceAccessibility<execution_space,
321 typename IdxView::memory_space>::accessible,
322 "IdxView not accessible from execution space");
325 typedef PackArraySingleColumnWithBoundsCheck<DstView,SrcView,IdxView> impl_type;
326 impl_type::pack (dst, src, idx, col, space);
329 typedef PackArraySingleColumn<DstView,SrcView,IdxView> impl_type;
330 impl_type::pack (dst, src, idx, col, space);
336 template <
typename DstView,
typename SrcView,
typename IdxView>
338 pack_array_single_column (
const DstView& dst,
342 const bool debug =
true)
344 pack_array_single_column(dst, src, idx, col, debug,
typename DstView::execution_space());
347 template <
typename DstView,
typename SrcView,
typename IdxView,
348 typename Enabled =
void>
349 struct PackArrayMultiColumn {
350 using execution_space =
typename DstView::execution_space;
351 typedef typename execution_space::size_type size_type;
358 PackArrayMultiColumn (
const DstView& dst_,
361 const size_t numCols_) :
362 dst(dst_), src(src_), idx(idx_), numCols(numCols_) {}
364 KOKKOS_INLINE_FUNCTION
void
365 operator() (
const size_type k)
const {
366 const typename IdxView::value_type localRow = idx(k);
367 const size_t offset = k*numCols;
368 for (
size_t j = 0; j < numCols; ++j) {
369 dst(offset + j) = src(localRow, j);
373 static void pack(
const DstView& dst,
377 const execution_space &space) {
378 typedef Kokkos::RangePolicy<execution_space, size_type> range_type;
380 (
"Tpetra::MultiVector pack multicol const stride",
381 range_type (space, 0, idx.size ()),
382 PackArrayMultiColumn (dst, src, idx, numCols));
386 template <
typename DstView,
389 typename SizeType =
typename DstView::execution_space::size_type,
390 typename Enabled =
void>
391 class PackArrayMultiColumnWithBoundsCheck {
393 using size_type = SizeType;
394 using value_type = size_t;
395 using execution_space =
typename DstView::execution_space;
404 PackArrayMultiColumnWithBoundsCheck (
const DstView& dst_,
407 const size_type numCols_) :
408 dst (dst_), src (src_), idx (idx_), numCols (numCols_) {}
410 KOKKOS_INLINE_FUNCTION
void
411 operator() (
const size_type k, value_type& lclErrorCount)
const {
412 typedef typename IdxView::non_const_value_type index_type;
414 const index_type lclRow = idx(k);
415 if (lclRow < static_cast<index_type> (0) ||
416 lclRow >= static_cast<index_type> (src.extent (0))) {
420 const size_type offset = k*numCols;
421 for (size_type j = 0; j < numCols; ++j) {
422 dst(offset + j) = src(lclRow, j);
427 KOKKOS_INLINE_FUNCTION
428 void init (value_type& initialErrorCount)
const {
429 initialErrorCount = 0;
432 KOKKOS_INLINE_FUNCTION
void
433 join (value_type& dstErrorCount,
434 const value_type& srcErrorCount)
const
436 dstErrorCount += srcErrorCount;
440 pack (
const DstView& dst,
443 const size_type numCols,
444 const execution_space &space)
446 typedef Kokkos::RangePolicy<execution_space, size_type> range_type;
447 typedef typename IdxView::non_const_value_type index_type;
449 size_t errorCount = 0;
450 Kokkos::parallel_reduce
451 (
"Tpetra::MultiVector pack multicol const stride debug only",
452 range_type (space, 0, idx.size ()),
453 PackArrayMultiColumnWithBoundsCheck (dst, src, idx, numCols),
455 if (errorCount != 0) {
459 auto idx_h = Kokkos::create_mirror_view (idx);
464 std::vector<index_type> badIndices;
465 const size_type numInds = idx_h.extent (0);
466 for (size_type k = 0; k < numInds; ++k) {
467 if (idx_h(k) < static_cast<index_type> (0) ||
468 idx_h(k) >= static_cast<index_type> (src.extent (0))) {
469 badIndices.push_back (idx_h(k));
473 TEUCHOS_TEST_FOR_EXCEPTION
474 (errorCount != badIndices.size (), std::logic_error,
475 "PackArrayMultiColumnWithBoundsCheck: errorCount = " << errorCount
476 <<
" != badIndices.size() = " << badIndices.size () <<
". This sho"
477 "uld never happen. Please report this to the Tpetra developers.");
479 std::ostringstream os;
480 os <<
"Tpetra::MultiVector multiple-column pack kernel had "
481 << badIndices.size () <<
" out-of bounds index/ices (errorCount = "
482 << errorCount <<
"): [";
483 for (
size_t k = 0; k < badIndices.size (); ++k) {
485 if (k + 1 < badIndices.size ()) {
490 throw std::runtime_error (os.str ());
496 template <
typename DstView,
500 pack_array_multi_column (
const DstView& dst,
503 const size_t numCols,
505 const typename DstView::execution_space &space)
507 static_assert (Kokkos::is_view<DstView>::value,
508 "DstView must be a Kokkos::View.");
509 static_assert (Kokkos::is_view<SrcView>::value,
510 "SrcView must be a Kokkos::View.");
511 static_assert (Kokkos::is_view<IdxView>::value,
512 "IdxView must be a Kokkos::View.");
513 static_assert (static_cast<int> (DstView::rank) == 1,
514 "DstView must be a rank-1 Kokkos::View.");
515 static_assert (static_cast<int> (SrcView::rank) == 2,
516 "SrcView must be a rank-2 Kokkos::View.");
517 static_assert (static_cast<int> (IdxView::rank) == 1,
518 "IdxView must be a rank-1 Kokkos::View.");
520 using execution_space =
typename DstView::execution_space;
522 static_assert (Kokkos::SpaceAccessibility<execution_space,
523 typename DstView::memory_space>::accessible,
524 "DstView not accessible from execution space");
525 static_assert (Kokkos::SpaceAccessibility<execution_space,
526 typename SrcView::memory_space>::accessible,
527 "SrcView not accessible from execution space");
528 static_assert (Kokkos::SpaceAccessibility<execution_space,
529 typename IdxView::memory_space>::accessible,
530 "IdxView not accessible from execution space");
533 typedef PackArrayMultiColumnWithBoundsCheck<DstView,
534 SrcView, IdxView> impl_type;
535 impl_type::pack (dst, src, idx, numCols, space);
538 typedef PackArrayMultiColumn<DstView, SrcView, IdxView> impl_type;
539 impl_type::pack (dst, src, idx, numCols, space);
543 template <
typename DstView,
547 pack_array_multi_column (
const DstView& dst,
550 const size_t numCols,
551 const bool debug =
true) {
552 pack_array_multi_column(dst, src, idx, numCols, debug,
typename DstView::execution_space());
555 template <
typename DstView,
typename SrcView,
typename IdxView,
556 typename ColView,
typename Enabled =
void>
557 struct PackArrayMultiColumnVariableStride {
558 using execution_space =
typename DstView::execution_space;
559 typedef typename execution_space::size_type size_type;
567 PackArrayMultiColumnVariableStride (
const DstView& dst_,
571 const size_t numCols_) :
572 dst(dst_), src(src_), idx(idx_), col(col_), numCols(numCols_) {}
574 KOKKOS_INLINE_FUNCTION
575 void operator() (
const size_type k)
const {
576 const typename IdxView::value_type localRow = idx(k);
577 const size_t offset = k*numCols;
578 for (
size_t j = 0; j < numCols; ++j) {
579 dst(offset + j) = src(localRow, col(j));
583 static void pack(
const DstView& dst,
588 const execution_space &space) {
589 typedef Kokkos::RangePolicy<execution_space, size_type> range_type;
591 (
"Tpetra::MultiVector pack multicol var stride",
592 range_type (space, 0, idx.size ()),
593 PackArrayMultiColumnVariableStride (dst, src, idx, col, numCols));
597 template <
typename DstView,
601 typename SizeType =
typename DstView::execution_space::size_type,
602 typename Enabled =
void>
603 class PackArrayMultiColumnVariableStrideWithBoundsCheck {
605 using size_type = SizeType;
606 using value_type = size_t;
607 using execution_space =
typename DstView::execution_space;
617 PackArrayMultiColumnVariableStrideWithBoundsCheck (
const DstView& dst_,
621 const size_type numCols_) :
622 dst (dst_), src (src_), idx (idx_), col (col_), numCols (numCols_) {}
624 KOKKOS_INLINE_FUNCTION
void
625 operator() (
const size_type k, value_type& lclErrorCount)
const {
626 typedef typename IdxView::non_const_value_type row_index_type;
627 typedef typename ColView::non_const_value_type col_index_type;
629 const row_index_type lclRow = idx(k);
630 if (lclRow < static_cast<row_index_type> (0) ||
631 lclRow >= static_cast<row_index_type> (src.extent (0))) {
635 const size_type offset = k*numCols;
636 for (size_type j = 0; j < numCols; ++j) {
637 const col_index_type lclCol = col(j);
638 if (Impl::outOfBounds<col_index_type> (lclCol, src.extent (1))) {
642 dst(offset + j) = src(lclRow, lclCol);
648 KOKKOS_INLINE_FUNCTION
void
649 init (value_type& initialErrorCount)
const {
650 initialErrorCount = 0;
653 KOKKOS_INLINE_FUNCTION
void
654 join (value_type& dstErrorCount,
655 const value_type& srcErrorCount)
const
657 dstErrorCount += srcErrorCount;
661 pack (
const DstView& dst,
665 const size_type numCols,
666 const execution_space &space)
668 using range_type = Kokkos::RangePolicy<execution_space, size_type>;
669 using row_index_type =
typename IdxView::non_const_value_type;
670 using col_index_type =
typename ColView::non_const_value_type;
672 size_t errorCount = 0;
673 Kokkos::parallel_reduce
674 (
"Tpetra::MultiVector pack multicol var stride debug only",
675 range_type (space, 0, idx.size ()),
676 PackArrayMultiColumnVariableStrideWithBoundsCheck (dst, src, idx,
679 if (errorCount != 0) {
680 constexpr
size_t maxNumBadIndicesToPrint = 100;
682 std::ostringstream os;
683 os <<
"Tpetra::MultiVector multicolumn variable stride pack kernel "
684 "found " << errorCount
685 <<
" error" << (errorCount != size_t (1) ?
"s" :
"") <<
". ";
690 auto idx_h = Kokkos::create_mirror_view (idx);
695 std::vector<row_index_type> badRows;
696 const size_type numRowInds = idx_h.extent (0);
697 for (size_type k = 0; k < numRowInds; ++k) {
698 if (Impl::outOfBounds<row_index_type> (idx_h(k), src.extent (0))) {
699 badRows.push_back (idx_h(k));
703 if (badRows.size () != 0) {
704 os << badRows.size () <<
" out-of-bounds row ind"
705 << (badRows.size () != size_t (1) ?
"ices" :
"ex");
706 if (badRows.size () <= maxNumBadIndicesToPrint) {
708 for (
size_t k = 0; k < badRows.size (); ++k) {
710 if (k + 1 < badRows.size ()) {
721 os <<
"No out-of-bounds row indices. ";
726 auto col_h = Kokkos::create_mirror_view (col);
731 std::vector<col_index_type> badCols;
732 const size_type numColInds = col_h.extent (0);
733 for (size_type k = 0; k < numColInds; ++k) {
734 if (Impl::outOfBounds<col_index_type> (col_h(k), src.extent (1))) {
735 badCols.push_back (col_h(k));
739 if (badCols.size () != 0) {
740 os << badCols.size () <<
" out-of-bounds column ind"
741 << (badCols.size () != size_t (1) ?
"ices" :
"ex");
742 if (badCols.size () <= maxNumBadIndicesToPrint) {
744 for (
size_t k = 0; k < badCols.size (); ++k) {
746 if (k + 1 < badCols.size ()) {
757 os <<
"No out-of-bounds column indices. ";
760 TEUCHOS_TEST_FOR_EXCEPTION
761 (errorCount != 0 && badRows.size () == 0 && badCols.size () == 0,
762 std::logic_error,
"Tpetra::MultiVector variable stride pack "
763 "kernel reports errorCount=" << errorCount <<
", but we failed "
764 "to find any bad rows or columns. This should never happen. "
765 "Please report this to the Tpetra developers.");
767 throw std::runtime_error (os.str ());
772 template <
typename DstView,
777 pack_array_multi_column_variable_stride (
const DstView& dst,
781 const size_t numCols,
783 const typename DstView::execution_space &space)
785 static_assert (Kokkos::is_view<DstView>::value,
786 "DstView must be a Kokkos::View.");
787 static_assert (Kokkos::is_view<SrcView>::value,
788 "SrcView must be a Kokkos::View.");
789 static_assert (Kokkos::is_view<IdxView>::value,
790 "IdxView must be a Kokkos::View.");
791 static_assert (Kokkos::is_view<ColView>::value,
792 "ColView must be a Kokkos::View.");
793 static_assert (static_cast<int> (DstView::rank) == 1,
794 "DstView must be a rank-1 Kokkos::View.");
795 static_assert (static_cast<int> (SrcView::rank) == 2,
796 "SrcView must be a rank-2 Kokkos::View.");
797 static_assert (static_cast<int> (IdxView::rank) == 1,
798 "IdxView must be a rank-1 Kokkos::View.");
799 static_assert (static_cast<int> (ColView::rank) == 1,
800 "ColView must be a rank-1 Kokkos::View.");
802 using execution_space =
typename DstView::execution_space;
804 static_assert (Kokkos::SpaceAccessibility<execution_space,
805 typename DstView::memory_space>::accessible,
806 "DstView not accessible from execution space");
807 static_assert (Kokkos::SpaceAccessibility<execution_space,
808 typename SrcView::memory_space>::accessible,
809 "SrcView not accessible from execution space");
810 static_assert (Kokkos::SpaceAccessibility<execution_space,
811 typename IdxView::memory_space>::accessible,
812 "IdxView not accessible from execution space");
815 typedef PackArrayMultiColumnVariableStrideWithBoundsCheck<DstView,
816 SrcView, IdxView, ColView> impl_type;
817 impl_type::pack (dst, src, idx, col, numCols, space);
820 typedef PackArrayMultiColumnVariableStride<DstView,
821 SrcView, IdxView, ColView> impl_type;
822 impl_type::pack (dst, src, idx, col, numCols, space);
826 template <
typename DstView,
831 pack_array_multi_column_variable_stride (
const DstView& dst,
835 const size_t numCols,
836 const bool debug =
true) {
837 pack_array_multi_column_variable_stride(dst, src, idx, col, numCols, debug,
838 typename DstView::execution_space());
843 struct atomic_tag {};
844 struct nonatomic_tag {};
848 KOKKOS_INLINE_FUNCTION
849 void operator() (atomic_tag, SC& dest,
const SC& src)
const {
850 Kokkos::atomic_add (&dest, src);
854 KOKKOS_INLINE_FUNCTION
855 void operator() (nonatomic_tag, SC& dest,
const SC& src)
const {
866 KOKKOS_INLINE_FUNCTION
867 void operator() (atomic_tag, SC& dest,
const SC& src)
const {
872 KOKKOS_INLINE_FUNCTION
873 void operator() (nonatomic_tag, SC& dest,
const SC& src)
const {
879 template <
class Scalar>
883 KOKKOS_FUNCTION AbsMaxHelper& operator+=(AbsMaxHelper
const& rhs) {
884 auto lhs_abs_value = Kokkos::ArithTraits<Scalar>::abs(value);
885 auto rhs_abs_value = Kokkos::ArithTraits<Scalar>::abs(rhs.value);
886 value = lhs_abs_value > rhs_abs_value ? lhs_abs_value : rhs_abs_value;
890 KOKKOS_FUNCTION AbsMaxHelper operator+(AbsMaxHelper
const& rhs)
const {
891 AbsMaxHelper ret = *
this;
897 template <
typename SC>
898 KOKKOS_INLINE_FUNCTION
899 void operator() (atomic_tag, SC& dst,
const SC& src)
const {
900 Kokkos::atomic_add(
reinterpret_cast<AbsMaxHelper<SC>*
>(&dst), AbsMaxHelper<SC>{src});
903 template <
typename SC>
904 KOKKOS_INLINE_FUNCTION
905 void operator() (nonatomic_tag, SC& dst,
const SC& src)
const {
906 auto dst_abs_value = Kokkos::ArithTraits<SC>::abs(dst);
907 auto src_abs_value = Kokkos::ArithTraits<SC>::abs(src);
908 dst = dst_abs_value > src_abs_value ? dst_abs_value : src_abs_value;
912 template <
typename ExecutionSpace,
917 typename Enabled =
void>
918 class UnpackArrayMultiColumn {
920 static_assert (Kokkos::is_view<DstView>::value,
921 "DstView must be a Kokkos::View.");
922 static_assert (Kokkos::is_view<SrcView>::value,
923 "SrcView must be a Kokkos::View.");
924 static_assert (Kokkos::is_view<IdxView>::value,
925 "IdxView must be a Kokkos::View.");
926 static_assert (static_cast<int> (DstView::rank) == 2,
927 "DstView must be a rank-2 Kokkos::View.");
928 static_assert (static_cast<int> (SrcView::rank) == 1,
929 "SrcView must be a rank-1 Kokkos::View.");
930 static_assert (static_cast<int> (IdxView::rank) == 1,
931 "IdxView must be a rank-1 Kokkos::View.");
934 typedef typename ExecutionSpace::execution_space execution_space;
935 typedef typename execution_space::size_type size_type;
945 UnpackArrayMultiColumn (
const ExecutionSpace& ,
950 const size_t numCols_) :
958 template<
class TagType>
959 KOKKOS_INLINE_FUNCTION
void
960 operator() (TagType tag,
const size_type k)
const
963 (std::is_same<TagType, atomic_tag>::value ||
964 std::is_same<TagType, nonatomic_tag>::value,
965 "TagType must be atomic_tag or nonatomic_tag.");
967 const typename IdxView::value_type localRow = idx(k);
968 const size_t offset = k*numCols;
969 for (
size_t j = 0; j < numCols; ++j) {
970 op (tag, dst(localRow, j), src(offset+j));
975 unpack (
const ExecutionSpace& execSpace,
980 const size_t numCols,
981 const bool use_atomic_updates)
983 if (use_atomic_updates) {
985 Kokkos::RangePolicy<atomic_tag, execution_space, size_type>;
987 (
"Tpetra::MultiVector unpack const stride atomic",
988 range_type (0, idx.size ()),
989 UnpackArrayMultiColumn (execSpace, dst, src, idx, op, numCols));
993 Kokkos::RangePolicy<nonatomic_tag, execution_space, size_type>;
995 (
"Tpetra::MultiVector unpack const stride nonatomic",
996 range_type (0, idx.size ()),
997 UnpackArrayMultiColumn (execSpace, dst, src, idx, op, numCols));
1002 template <
typename ExecutionSpace,
1007 typename SizeType =
typename ExecutionSpace::execution_space::size_type,
1008 typename Enabled =
void>
1009 class UnpackArrayMultiColumnWithBoundsCheck {
1011 static_assert (Kokkos::is_view<DstView>::value,
1012 "DstView must be a Kokkos::View.");
1013 static_assert (Kokkos::is_view<SrcView>::value,
1014 "SrcView must be a Kokkos::View.");
1015 static_assert (Kokkos::is_view<IdxView>::value,
1016 "IdxView must be a Kokkos::View.");
1017 static_assert (static_cast<int> (DstView::rank) == 2,
1018 "DstView must be a rank-2 Kokkos::View.");
1019 static_assert (static_cast<int> (SrcView::rank) == 1,
1020 "SrcView must be a rank-1 Kokkos::View.");
1021 static_assert (static_cast<int> (IdxView::rank) == 1,
1022 "IdxView must be a rank-1 Kokkos::View.");
1023 static_assert (std::is_integral<SizeType>::value,
1024 "SizeType must be a built-in integer type.");
1027 using execution_space =
typename ExecutionSpace::execution_space;
1028 using size_type = SizeType;
1029 using value_type = size_t;
1039 UnpackArrayMultiColumnWithBoundsCheck (
const ExecutionSpace& ,
1040 const DstView& dst_,
1041 const SrcView& src_,
1042 const IdxView& idx_,
1044 const size_type numCols_) :
1052 template<
class TagType>
1053 KOKKOS_INLINE_FUNCTION
void
1054 operator() (TagType tag,
1056 size_t& lclErrCount)
const
1059 (std::is_same<TagType, atomic_tag>::value ||
1060 std::is_same<TagType, nonatomic_tag>::value,
1061 "TagType must be atomic_tag or nonatomic_tag.");
1062 using index_type =
typename IdxView::non_const_value_type;
1064 const index_type lclRow = idx(k);
1065 if (lclRow < static_cast<index_type> (0) ||
1066 lclRow >= static_cast<index_type> (dst.extent (0))) {
1070 const size_type offset = k*numCols;
1071 for (size_type j = 0; j < numCols; ++j) {
1072 op (tag, dst(lclRow,j), src(offset+j));
1077 template<
class TagType>
1078 KOKKOS_INLINE_FUNCTION
void
1079 init (TagType,
size_t& initialErrorCount)
const {
1080 initialErrorCount = 0;
1083 template<
class TagType>
1084 KOKKOS_INLINE_FUNCTION
void
1086 size_t& dstErrorCount,
1087 const size_t& srcErrorCount)
const
1089 dstErrorCount += srcErrorCount;
1093 unpack (
const ExecutionSpace& execSpace,
1098 const size_type numCols,
1099 const bool use_atomic_updates)
1101 using index_type =
typename IdxView::non_const_value_type;
1103 size_t errorCount = 0;
1104 if (use_atomic_updates) {
1106 Kokkos::RangePolicy<atomic_tag, execution_space, size_type>;
1107 Kokkos::parallel_reduce
1108 (
"Tpetra::MultiVector unpack multicol const stride atomic debug only",
1109 range_type (0, idx.size ()),
1110 UnpackArrayMultiColumnWithBoundsCheck (execSpace, dst, src,
1116 Kokkos::RangePolicy<nonatomic_tag, execution_space, size_type>;
1117 Kokkos::parallel_reduce
1118 (
"Tpetra::MultiVector unpack multicol const stride nonatomic debug only",
1119 range_type (0, idx.size ()),
1120 UnpackArrayMultiColumnWithBoundsCheck (execSpace, dst, src,
1125 if (errorCount != 0) {
1129 auto idx_h = Kokkos::create_mirror_view (idx);
1134 std::vector<index_type> badIndices;
1135 const size_type numInds = idx_h.extent (0);
1136 for (size_type k = 0; k < numInds; ++k) {
1137 if (idx_h(k) < static_cast<index_type> (0) ||
1138 idx_h(k) >= static_cast<index_type> (dst.extent (0))) {
1139 badIndices.push_back (idx_h(k));
1143 if (errorCount != badIndices.size ()) {
1144 std::ostringstream os;
1145 os <<
"MultiVector unpack kernel: errorCount = " << errorCount
1146 <<
" != badIndices.size() = " << badIndices.size ()
1147 <<
". This should never happen. "
1148 "Please report this to the Tpetra developers.";
1149 throw std::logic_error (os.str ());
1152 std::ostringstream os;
1153 os <<
"MultiVector unpack kernel had " << badIndices.size ()
1154 <<
" out-of bounds index/ices. Here they are: [";
1155 for (
size_t k = 0; k < badIndices.size (); ++k) {
1156 os << badIndices[k];
1157 if (k + 1 < badIndices.size ()) {
1162 throw std::runtime_error (os.str ());
1167 template <
typename ExecutionSpace,
1173 unpack_array_multi_column (
const ExecutionSpace& execSpace,
1178 const size_t numCols,
1179 const bool use_atomic_updates,
1182 static_assert (Kokkos::is_view<DstView>::value,
1183 "DstView must be a Kokkos::View.");
1184 static_assert (Kokkos::is_view<SrcView>::value,
1185 "SrcView must be a Kokkos::View.");
1186 static_assert (Kokkos::is_view<IdxView>::value,
1187 "IdxView must be a Kokkos::View.");
1188 static_assert (static_cast<int> (DstView::rank) == 2,
1189 "DstView must be a rank-2 Kokkos::View.");
1190 static_assert (static_cast<int> (SrcView::rank) == 1,
1191 "SrcView must be a rank-1 Kokkos::View.");
1192 static_assert (static_cast<int> (IdxView::rank) == 1,
1193 "IdxView must be a rank-1 Kokkos::View.");
1196 typedef UnpackArrayMultiColumnWithBoundsCheck<ExecutionSpace,
1197 DstView, SrcView, IdxView, Op> impl_type;
1198 impl_type::unpack (execSpace, dst, src, idx, op, numCols,
1199 use_atomic_updates);
1202 typedef UnpackArrayMultiColumn<ExecutionSpace,
1203 DstView, SrcView, IdxView, Op> impl_type;
1204 impl_type::unpack (execSpace, dst, src, idx, op, numCols,
1205 use_atomic_updates);
1209 template <
typename ExecutionSpace,
1215 typename Enabled =
void>
1216 class UnpackArrayMultiColumnVariableStride {
1218 static_assert (Kokkos::is_view<DstView>::value,
1219 "DstView must be a Kokkos::View.");
1220 static_assert (Kokkos::is_view<SrcView>::value,
1221 "SrcView must be a Kokkos::View.");
1222 static_assert (Kokkos::is_view<IdxView>::value,
1223 "IdxView must be a Kokkos::View.");
1224 static_assert (Kokkos::is_view<ColView>::value,
1225 "ColView must be a Kokkos::View.");
1226 static_assert (static_cast<int> (DstView::rank) == 2,
1227 "DstView must be a rank-2 Kokkos::View.");
1228 static_assert (static_cast<int> (SrcView::rank) == 1,
1229 "SrcView must be a rank-1 Kokkos::View.");
1230 static_assert (static_cast<int> (IdxView::rank) == 1,
1231 "IdxView must be a rank-1 Kokkos::View.");
1232 static_assert (static_cast<int> (ColView::rank) == 1,
1233 "ColView must be a rank-1 Kokkos::View.");
1236 using execution_space =
typename ExecutionSpace::execution_space;
1237 using size_type =
typename execution_space::size_type;
1248 UnpackArrayMultiColumnVariableStride (
const ExecutionSpace& ,
1249 const DstView& dst_,
1250 const SrcView& src_,
1251 const IdxView& idx_,
1252 const ColView& col_,
1254 const size_t numCols_) :
1263 template<
class TagType>
1264 KOKKOS_INLINE_FUNCTION
void
1265 operator() (TagType tag,
const size_type k)
const
1268 (std::is_same<TagType, atomic_tag>::value ||
1269 std::is_same<TagType, nonatomic_tag>::value,
1270 "TagType must be atomic_tag or nonatomic_tag.");
1272 const typename IdxView::value_type localRow = idx(k);
1273 const size_t offset = k*numCols;
1274 for (
size_t j = 0; j < numCols; ++j) {
1275 op (tag, dst(localRow, col(j)), src(offset+j));
1280 unpack (
const ExecutionSpace& execSpace,
1286 const size_t numCols,
1287 const bool use_atomic_updates)
1289 if (use_atomic_updates) {
1291 Kokkos::RangePolicy<atomic_tag, execution_space, size_type>;
1292 Kokkos::parallel_for
1293 (
"Tpetra::MultiVector unpack var stride atomic",
1294 range_type (0, idx.size ()),
1295 UnpackArrayMultiColumnVariableStride (execSpace, dst, src,
1296 idx, col, op, numCols));
1300 Kokkos::RangePolicy<nonatomic_tag, execution_space, size_type>;
1301 Kokkos::parallel_for
1302 (
"Tpetra::MultiVector unpack var stride nonatomic",
1303 range_type (0, idx.size ()),
1304 UnpackArrayMultiColumnVariableStride (execSpace, dst, src,
1305 idx, col, op, numCols));
1310 template <
typename ExecutionSpace,
1316 typename SizeType =
typename ExecutionSpace::execution_space::size_type,
1317 typename Enabled =
void>
1318 class UnpackArrayMultiColumnVariableStrideWithBoundsCheck {
1320 static_assert (Kokkos::is_view<DstView>::value,
1321 "DstView must be a Kokkos::View.");
1322 static_assert (Kokkos::is_view<SrcView>::value,
1323 "SrcView must be a Kokkos::View.");
1324 static_assert (Kokkos::is_view<IdxView>::value,
1325 "IdxView must be a Kokkos::View.");
1326 static_assert (Kokkos::is_view<ColView>::value,
1327 "ColView must be a Kokkos::View.");
1328 static_assert (static_cast<int> (DstView::rank) == 2,
1329 "DstView must be a rank-2 Kokkos::View.");
1330 static_assert (static_cast<int> (SrcView::rank) == 1,
1331 "SrcView must be a rank-1 Kokkos::View.");
1332 static_assert (static_cast<int> (IdxView::rank) == 1,
1333 "IdxView must be a rank-1 Kokkos::View.");
1334 static_assert (static_cast<int> (ColView::rank) == 1,
1335 "ColView must be a rank-1 Kokkos::View.");
1336 static_assert (std::is_integral<SizeType>::value,
1337 "SizeType must be a built-in integer type.");
1340 using execution_space =
typename ExecutionSpace::execution_space;
1341 using size_type = SizeType;
1342 using value_type = size_t;
1353 UnpackArrayMultiColumnVariableStrideWithBoundsCheck
1354 (
const ExecutionSpace& ,
1355 const DstView& dst_,
1356 const SrcView& src_,
1357 const IdxView& idx_,
1358 const ColView& col_,
1360 const size_t numCols_) :
1369 template<
class TagType>
1370 KOKKOS_INLINE_FUNCTION
void
1371 operator() (TagType tag,
1373 value_type& lclErrorCount)
const
1376 (std::is_same<TagType, atomic_tag>::value ||
1377 std::is_same<TagType, nonatomic_tag>::value,
1378 "TagType must be atomic_tag or nonatomic_tag.");
1379 using row_index_type =
typename IdxView::non_const_value_type;
1380 using col_index_type =
typename ColView::non_const_value_type;
1382 const row_index_type lclRow = idx(k);
1383 if (lclRow < static_cast<row_index_type> (0) ||
1384 lclRow >= static_cast<row_index_type> (dst.extent (0))) {
1388 const size_type offset = k * numCols;
1389 for (size_type j = 0; j < numCols; ++j) {
1390 const col_index_type lclCol = col(j);
1391 if (Impl::outOfBounds<col_index_type> (lclCol, dst.extent (1))) {
1395 op (tag, dst(lclRow, col(j)), src(offset+j));
1401 KOKKOS_INLINE_FUNCTION
void
1402 init (value_type& initialErrorCount)
const {
1403 initialErrorCount = 0;
1406 KOKKOS_INLINE_FUNCTION
void
1407 join (value_type& dstErrorCount,
1408 const value_type& srcErrorCount)
const
1410 dstErrorCount += srcErrorCount;
1414 unpack (
const ExecutionSpace& execSpace,
1420 const size_type numCols,
1421 const bool use_atomic_updates)
1423 using row_index_type =
typename IdxView::non_const_value_type;
1424 using col_index_type =
typename ColView::non_const_value_type;
1426 size_t errorCount = 0;
1427 if (use_atomic_updates) {
1429 Kokkos::RangePolicy<atomic_tag, execution_space, size_type>;
1430 Kokkos::parallel_reduce
1431 (
"Tpetra::MultiVector unpack var stride atomic debug only",
1432 range_type (0, idx.size ()),
1433 UnpackArrayMultiColumnVariableStrideWithBoundsCheck
1434 (execSpace, dst, src, idx, col, op, numCols),
1439 Kokkos::RangePolicy<nonatomic_tag, execution_space, size_type>;
1440 Kokkos::parallel_reduce
1441 (
"Tpetra::MultiVector unpack var stride nonatomic debug only",
1442 range_type (0, idx.size ()),
1443 UnpackArrayMultiColumnVariableStrideWithBoundsCheck
1444 (execSpace, dst, src, idx, col, op, numCols),
1448 if (errorCount != 0) {
1449 constexpr
size_t maxNumBadIndicesToPrint = 100;
1451 std::ostringstream os;
1452 os <<
"Tpetra::MultiVector multicolumn variable stride unpack kernel "
1453 "found " << errorCount
1454 <<
" error" << (errorCount != size_t (1) ?
"s" :
"") <<
". ";
1460 auto idx_h = Kokkos::create_mirror_view (idx);
1465 std::vector<row_index_type> badRows;
1466 const size_type numRowInds = idx_h.extent (0);
1467 for (size_type k = 0; k < numRowInds; ++k) {
1468 if (idx_h(k) < static_cast<row_index_type> (0) ||
1469 idx_h(k) >= static_cast<row_index_type> (dst.extent (0))) {
1470 badRows.push_back (idx_h(k));
1474 if (badRows.size () != 0) {
1475 os << badRows.size () <<
" out-of-bounds row ind"
1476 << (badRows.size () != size_t (1) ?
"ices" :
"ex");
1477 if (badRows.size () <= maxNumBadIndicesToPrint) {
1479 for (
size_t k = 0; k < badRows.size (); ++k) {
1481 if (k + 1 < badRows.size ()) {
1492 os <<
"No out-of-bounds row indices. ";
1497 auto col_h = Kokkos::create_mirror_view (col);
1502 std::vector<col_index_type> badCols;
1503 const size_type numColInds = col_h.extent (0);
1504 for (size_type k = 0; k < numColInds; ++k) {
1505 if (Impl::outOfBounds<col_index_type> (col_h(k), dst.extent (1))) {
1506 badCols.push_back (col_h(k));
1510 if (badCols.size () != 0) {
1511 os << badCols.size () <<
" out-of-bounds column ind"
1512 << (badCols.size () != size_t (1) ?
"ices" :
"ex");
1513 if (badCols.size () <= maxNumBadIndicesToPrint) {
1514 for (
size_t k = 0; k < badCols.size (); ++k) {
1517 if (k + 1 < badCols.size ()) {
1528 os <<
"No out-of-bounds column indices. ";
1531 TEUCHOS_TEST_FOR_EXCEPTION
1532 (errorCount != 0 && badRows.size () == 0 && badCols.size () == 0,
1533 std::logic_error,
"Tpetra::MultiVector variable stride unpack "
1534 "kernel reports errorCount=" << errorCount <<
", but we failed "
1535 "to find any bad rows or columns. This should never happen. "
1536 "Please report this to the Tpetra developers.");
1538 throw std::runtime_error (os.str ());
1543 template <
typename ExecutionSpace,
1550 unpack_array_multi_column_variable_stride (
const ExecutionSpace& execSpace,
1556 const size_t numCols,
1557 const bool use_atomic_updates,
1560 static_assert (Kokkos::is_view<DstView>::value,
1561 "DstView must be a Kokkos::View.");
1562 static_assert (Kokkos::is_view<SrcView>::value,
1563 "SrcView must be a Kokkos::View.");
1564 static_assert (Kokkos::is_view<IdxView>::value,
1565 "IdxView must be a Kokkos::View.");
1566 static_assert (Kokkos::is_view<ColView>::value,
1567 "ColView must be a Kokkos::View.");
1568 static_assert (static_cast<int> (DstView::rank) == 2,
1569 "DstView must be a rank-2 Kokkos::View.");
1570 static_assert (static_cast<int> (SrcView::rank) == 1,
1571 "SrcView must be a rank-1 Kokkos::View.");
1572 static_assert (static_cast<int> (IdxView::rank) == 1,
1573 "IdxView must be a rank-1 Kokkos::View.");
1574 static_assert (static_cast<int> (ColView::rank) == 1,
1575 "ColView must be a rank-1 Kokkos::View.");
1579 UnpackArrayMultiColumnVariableStrideWithBoundsCheck<ExecutionSpace,
1580 DstView, SrcView, IdxView, ColView, Op>;
1581 impl_type::unpack (execSpace, dst, src, idx, col, op, numCols,
1582 use_atomic_updates);
1585 using impl_type = UnpackArrayMultiColumnVariableStride<ExecutionSpace,
1586 DstView, SrcView, IdxView, ColView, Op>;
1587 impl_type::unpack (execSpace, dst, src, idx, col, op, numCols,
1588 use_atomic_updates);
1592 template <
typename DstView,
typename SrcView,
1593 typename DstIdxView,
typename SrcIdxView,
typename Op,
1594 typename Enabled =
void>
1595 struct PermuteArrayMultiColumn {
1596 using size_type =
typename DstView::size_type;
1605 PermuteArrayMultiColumn (
const DstView& dst_,
1606 const SrcView& src_,
1607 const DstIdxView& dst_idx_,
1608 const SrcIdxView& src_idx_,
1609 const size_t numCols_,
1611 dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
1612 numCols(numCols_), op(op_) {}
1614 KOKKOS_INLINE_FUNCTION
void
1615 operator() (
const size_type k)
const {
1616 const typename DstIdxView::value_type toRow = dst_idx(k);
1617 const typename SrcIdxView::value_type fromRow = src_idx(k);
1619 for (
size_t j = 0; j < numCols; ++j) {
1620 op(tag, dst(toRow, j), src(fromRow, j));
1624 template <
typename ExecutionSpace>
1626 permute (
const ExecutionSpace &space,
1629 const DstIdxView& dst_idx,
1630 const SrcIdxView& src_idx,
1631 const size_t numCols,
1635 Kokkos::RangePolicy<ExecutionSpace, size_type>;
1637 const size_type n = std::min (dst_idx.size (), src_idx.size ());
1638 Kokkos::parallel_for
1639 (
"Tpetra::MultiVector permute multicol const stride",
1640 range_type (space, 0, n),
1641 PermuteArrayMultiColumn (dst, src, dst_idx, src_idx, numCols, op));
1648 template <
typename ExecutionSpace,
typename DstView,
typename SrcView,
1649 typename DstIdxView,
typename SrcIdxView,
typename Op>
1650 void permute_array_multi_column(
const ExecutionSpace &space,
const DstView &dst,
1651 const SrcView &src,
const DstIdxView &dst_idx,
1652 const SrcIdxView &src_idx,
size_t numCols,
1654 PermuteArrayMultiColumn<DstView, SrcView, DstIdxView, SrcIdxView,
1655 Op>::permute(space, dst, src, dst_idx, src_idx,
1662 template <
typename DstView,
typename SrcView,
1663 typename DstIdxView,
typename SrcIdxView,
typename Op>
1664 void permute_array_multi_column(
const DstView& dst,
1666 const DstIdxView& dst_idx,
1667 const SrcIdxView& src_idx,
1670 using execution_space =
typename DstView::execution_space;
1671 PermuteArrayMultiColumn<DstView,SrcView,DstIdxView,SrcIdxView,Op>::permute(
1672 execution_space(), dst, src, dst_idx, src_idx, numCols, op);
1675 template <
typename DstView,
typename SrcView,
1676 typename DstIdxView,
typename SrcIdxView,
1677 typename DstColView,
typename SrcColView,
typename Op,
1678 typename Enabled =
void>
1679 struct PermuteArrayMultiColumnVariableStride {
1680 using size_type =
typename DstView::size_type;
1691 PermuteArrayMultiColumnVariableStride(
const DstView& dst_,
1692 const SrcView& src_,
1693 const DstIdxView& dst_idx_,
1694 const SrcIdxView& src_idx_,
1695 const DstColView& dst_col_,
1696 const SrcColView& src_col_,
1697 const size_t numCols_,
1699 dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
1700 dst_col(dst_col_), src_col(src_col_),
1701 numCols(numCols_), op(op_) {}
1703 KOKKOS_INLINE_FUNCTION
void
1704 operator() (
const size_type k)
const {
1705 const typename DstIdxView::value_type toRow = dst_idx(k);
1706 const typename SrcIdxView::value_type fromRow = src_idx(k);
1707 const nonatomic_tag tag;
1708 for (
size_t j = 0; j < numCols; ++j) {
1709 op(tag, dst(toRow, dst_col(j)), src(fromRow, src_col(j)));
1713 template <
typename ExecutionSpace>
1715 permute (
const ExecutionSpace &space,
1718 const DstIdxView& dst_idx,
1719 const SrcIdxView& src_idx,
1720 const DstColView& dst_col,
1721 const SrcColView& src_col,
1722 const size_t numCols,
1726 static_assert(Kokkos::SpaceAccessibility<
1727 ExecutionSpace,
typename DstView::memory_space>::accessible,
1728 "ExecutionSpace must be able to access DstView");
1730 using range_type = Kokkos::RangePolicy<ExecutionSpace, size_type>;
1731 const size_type n = std::min (dst_idx.size (), src_idx.size ());
1732 Kokkos::parallel_for
1733 (
"Tpetra::MultiVector permute multicol var stride",
1734 range_type (space, 0, n),
1735 PermuteArrayMultiColumnVariableStride (dst, src, dst_idx, src_idx,
1736 dst_col, src_col, numCols, op));
1743 template <
typename ExecutionSpace,
typename DstView,
typename SrcView,
1744 typename DstIdxView,
typename SrcIdxView,
typename DstColView,
1745 typename SrcColView,
typename Op>
1746 void permute_array_multi_column_variable_stride(
1747 const ExecutionSpace &space,
const DstView &dst,
const SrcView &src,
1748 const DstIdxView &dst_idx,
const SrcIdxView &src_idx,
1749 const DstColView &dst_col,
const SrcColView &src_col,
size_t numCols,
1751 PermuteArrayMultiColumnVariableStride<DstView, SrcView, DstIdxView,
1752 SrcIdxView, DstColView, SrcColView,
1753 Op>::permute(space, dst, src, dst_idx,
1754 src_idx, dst_col, src_col,
1761 template <
typename DstView,
typename SrcView,
1762 typename DstIdxView,
typename SrcIdxView,
1763 typename DstColView,
typename SrcColView,
typename Op>
1764 void permute_array_multi_column_variable_stride(
const DstView& dst,
1766 const DstIdxView& dst_idx,
1767 const SrcIdxView& src_idx,
1768 const DstColView& dst_col,
1769 const SrcColView& src_col,
1772 using execution_space =
typename DstView::execution_space;
1773 PermuteArrayMultiColumnVariableStride<DstView,SrcView,
1774 DstIdxView,SrcIdxView,DstColView,SrcColView,Op>::permute(
1775 execution_space(), dst, src, dst_idx, src_idx, dst_col, src_col, numCols, op);
1782 #endif // TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_HPP
KOKKOS_INLINE_FUNCTION bool outOfBounds(const IntegerType x, const IntegerType exclusiveUpperBound)
Is x out of bounds? That is, is x less than zero, or greater than or equal to the given exclusive upp...
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.
Is x out of bounds? That is, is x less than zero, or greater than or equal to the given exclusive upp...