44 #ifndef TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_UQ_PCE_HPP 45 #define TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_UQ_PCE_HPP 51 #include "Tpetra_KokkosRefactor_Details_MultiVectorDistObjectKernels.hpp" 59 namespace KokkosRefactor {
65 template <
typename DS,
typename ... DP,
66 typename SS,
typename ... SP,
68 struct PackArraySingleColumn<
70 Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
73 typedef Kokkos::View<Sacado::UQ::PCE<DS>*,DP...>
DstView;
74 typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>
SrcView;
76 typedef typename execution_space::size_type
size_type;
78 static const unsigned BlockSize = 32;
89 dst(dst_), src(src_), idx(idx_), col(col_) {}
91 KOKKOS_INLINE_FUNCTION
93 dst(k) = src(idx(k), col);
96 KOKKOS_INLINE_FUNCTION
99 dst(k).fastAccessCoeff(i) = src(idx(k), col).fastAccessCoeff(i);
107 Kokkos::parallel_for(
109 PackArraySingleColumn(dst,src,idx,col) );
111 Kokkos::parallel_for( idx.size(),
112 PackArraySingleColumn(dst,src,idx,col) );
116 template <
typename DS,
typename ... DP,
117 typename SS,
typename ... SP,
119 struct PackArrayMultiColumn<
121 Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
124 typedef Kokkos::View<Sacado::UQ::PCE<DS>*,DP...>
DstView;
125 typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>
SrcView;
129 static const unsigned BlockSize = 32;
140 dst(dst_), src(src_), idx(idx_), numCols(numCols_) {}
142 KOKKOS_INLINE_FUNCTION
145 const size_t offset = k*numCols;
146 for (
size_t j = 0;
j < numCols; ++
j)
147 dst(offset +
j) = src(localRow,
j);
150 KOKKOS_INLINE_FUNCTION
153 const size_t offset = k*numCols;
154 for (
size_t j = 0;
j < numCols; ++
j)
156 dst(offset +
j).fastAccessCoeff(i) =
157 src(localRow,
j).fastAccessCoeff(i);
165 Kokkos::parallel_for(
167 PackArrayMultiColumn(dst,src,idx,numCols) );
169 Kokkos::parallel_for( idx.size(),
170 PackArrayMultiColumn(dst,src,idx,numCols) );
174 template <
typename DS,
typename ... DP,
175 typename SS,
typename ... SP,
178 struct PackArrayMultiColumnVariableStride<
180 Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
184 typedef Kokkos::View<Sacado::UQ::PCE<DS>*,DP...>
DstView;
185 typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>
SrcView;
189 static const unsigned BlockSize = 32;
202 dst(dst_), src(src_), idx(idx_), col(col_), numCols(numCols_) {}
204 KOKKOS_INLINE_FUNCTION
207 const size_t offset = k*numCols;
208 for (
size_t j = 0;
j < numCols; ++
j)
209 dst(offset +
j) = src(localRow, col(
j));
212 KOKKOS_INLINE_FUNCTION
215 const size_t offset = k*numCols;
216 for (
size_t j = 0;
j < numCols; ++
j)
218 dst(offset +
j).fastAccessCoeff(i) =
219 src(localRow, col(
j)).fastAccessCoeff(i);
228 Kokkos::parallel_for(
230 PackArrayMultiColumnVariableStride(dst,src,idx,col,numCols) );
232 Kokkos::parallel_for( idx.size(),
233 PackArrayMultiColumnVariableStride(
234 dst,src,idx,col,numCols) );
238 template <
typename ExecutionSpace,
245 struct UnpackArrayMultiColumn<
248 Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>,
252 typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>
DstView;
253 typedef Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>
SrcView;
257 static const unsigned BlockSize = 32;
270 const size_t numCols_) :
278 KOKKOS_INLINE_FUNCTION
void 282 const size_t offset = k*numCols;
283 for (
size_t j = 0;
j < numCols; ++
j) {
284 op (dst(localRow,
j), src(offset+
j));
288 KOKKOS_INLINE_FUNCTION
void 292 const size_t offset = k*numCols;
293 for (
size_t j = 0;
j < numCols; ++
j) {
307 const size_t numCols)
311 (
"Tpetra::MultiVector (Sacado::UQ::PCE) unpack (constant stride)",
313 UnpackArrayMultiColumn (execSpace, dst, src, idx, op, numCols));
317 (
"Tpetra::MultiVector (Sacado::UQ::PCE) unpack (constant stride)",
318 Kokkos::RangePolicy<execution_space, size_type> (0, idx.size ()),
319 UnpackArrayMultiColumn (execSpace, dst, src, idx, op, numCols));
324 template <
typename ExecutionSpace,
332 struct UnpackArrayMultiColumnVariableStride<
335 Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>,
340 typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>
DstView;
341 typedef Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>
SrcView;
345 static const unsigned BlockSize = 32;
361 dst(dst_), src(src_), idx(idx_), col(col_), op(op_), numCols(numCols_) {}
363 KOKKOS_INLINE_FUNCTION
366 const size_t offset = k*numCols;
367 for (
size_t j = 0;
j < numCols; ++
j)
368 op( dst(localRow,col(
j)), src(offset+
j) );
371 KOKKOS_INLINE_FUNCTION
374 const size_t offset = k*numCols;
375 for (
size_t j = 0;
j < numCols; ++
j)
377 op( dst(localRow,col(
j)).fastAccessCoeff(i),
378 src(offset+
j).fastAccessCoeff(i) );
388 const size_t numCols)
392 (
"Tpetra::MultiVector unpack (Sacado::UQ::PCE) (nonconstant stride)",
394 UnpackArrayMultiColumnVariableStride (execSpace, dst, src,
395 idx, col, op, numCols));
399 (
"Tpetra::MultiVector unpack (Sacado::UQ::PCE) (nonconstant stride)",
400 Kokkos::RangePolicy<execution_space, size_type> (0, idx.size ()),
401 UnpackArrayMultiColumnVariableStride (execSpace, dst, src,
402 idx, col, op, numCols));
407 template <
typename DS,
typename ... DP,
408 typename SS,
typename ... SP,
409 typename DstIdxView,
typename SrcIdxView>
410 struct PermuteArrayMultiColumn<
412 Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
416 typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>
DstView;
417 typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>
SrcView;
421 static const unsigned BlockSize = 32;
431 const DstIdxView& dst_idx_,
432 const SrcIdxView& src_idx_,
434 dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
437 KOKKOS_INLINE_FUNCTION
441 for (
size_t j = 0;
j < numCols; ++
j)
442 dst(toRow,
j) = src(fromRow,
j);
445 KOKKOS_INLINE_FUNCTION
449 for (
size_t j = 0;
j < numCols; ++
j)
451 dst(toRow,
j).fastAccessCoeff(i) =
452 src(fromRow,
j).fastAccessCoeff(i);
457 const DstIdxView& dst_idx,
458 const SrcIdxView& src_idx,
462 Kokkos::parallel_for(
464 PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
466 Kokkos::parallel_for(
467 n, PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
471 template <
typename DS,
typename ... DP,
472 typename SS,
typename ... SP,
473 typename DstIdxView,
typename SrcIdxView,
474 typename DstColView,
typename SrcColView>
475 struct PermuteArrayMultiColumnVariableStride<
477 Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
478 DstIdxView, SrcIdxView,
479 DstColView, SrcColView >
481 typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>
DstView;
482 typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>
SrcView;
486 static const unsigned BlockSize = 32;
498 const DstIdxView& dst_idx_,
499 const SrcIdxView& src_idx_,
500 const DstColView& dst_col_,
501 const SrcColView& src_col_,
503 dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
504 dst_col(dst_col_), src_col(src_col_),
507 KOKKOS_INLINE_FUNCTION
511 for (
size_t j = 0;
j < numCols; ++
j)
512 dst(toRow, dst_col(
j)) = src(fromRow, src_col(
j));
515 KOKKOS_INLINE_FUNCTION
519 for (
size_t j = 0;
j < numCols; ++
j)
522 src(fromRow, src_col(
j)).fastAccessCoeff(i);
527 const DstIdxView& dst_idx,
528 const SrcIdxView& src_idx,
529 const DstColView& dst_col,
530 const SrcColView& src_col,
534 Kokkos::parallel_for(
536 PermuteArrayMultiColumnVariableStride(
537 dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
539 Kokkos::parallel_for(
540 n, PermuteArrayMultiColumnVariableStride(
541 dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
549 #endif // TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_UQ_PCE_HPP
Kokkos::View< Sacado::UQ::PCE< DS > **, DP... > DstView
static void pack(const DstView &dst, const SrcView &src, const IdxView &idx, const ColView &col, size_t numCols)
ExecutionSpace::execution_space execution_space
static void permute(const DstView &dst, const SrcView &src, const DstIdxView &dst_idx, const SrcIdxView &src_idx, const DstColView &dst_col, const SrcColView &src_col, size_t numCols)
PermuteArrayMultiColumn(const DstView &dst_, const SrcView &src_, const DstIdxView &dst_idx_, const SrcIdxView &src_idx_, size_t numCols_)
Kokkos::View< const Sacado::UQ::PCE< SS > **, SP... > SrcView
static void unpack(const ExecutionSpace &execSpace, const DstView &dst, const SrcView &src, const IdxView &idx, const ColView &col, const Op &op, const size_t numCols)
Kokkos::View< const Sacado::UQ::PCE< SS > **, SP... > SrcView
Kokkos::View< Sacado::UQ::PCE< DS > **, DP... > DstView
KOKKOS_INLINE_FUNCTION void operator()(const size_type k, const size_type tidx) const
Kokkos::View< const Sacado::UQ::PCE< SS > *, SP... > SrcView
KOKKOS_INLINE_FUNCTION void operator()(const size_type k, const size_type tidx) const
Kokkos::View< Sacado::UQ::PCE< DS > *, DP... > DstView
Kokkos::View< const Sacado::UQ::PCE< SS > **, SP... > SrcView
Kokkos::DefaultExecutionSpace execution_space
expr1 expr1 expr1 expr2 expr1 expr1 c expr2 expr1 c fastAccessCoeff(j) - expr2.val(j)
execution_space::size_type size_type
static void pack(const DstView &dst, const SrcView &src, const IdxView &idx, size_t numCols)
KOKKOS_INLINE_FUNCTION void operator()(const size_type k) const
DstView::execution_space execution_space
Kokkos::View< Sacado::UQ::PCE< DS > *, DP... > DstView
DstView::execution_space execution_space
KOKKOS_INLINE_FUNCTION void operator()(const size_type k) const
PermuteArrayMultiColumnVariableStride(const DstView &dst_, const SrcView &src_, const DstIdxView &dst_idx_, const SrcIdxView &src_idx_, const DstColView &dst_col_, const SrcColView &src_col_, size_t numCols_)
ExecutionSpace::execution_space execution_space
static void permute(const DstView &dst, const SrcView &src, const DstIdxView &dst_idx, const SrcIdxView &src_idx, size_t numCols)
KOKKOS_INLINE_FUNCTION PCE< Storage > min(const typename PCE< Storage >::value_type &a, const PCE< Storage > &b)
Kokkos::View< Sacado::UQ::PCE< DS > **, DP... > DstView
static void unpack(const ExecutionSpace &execSpace, const DstView &dst, const SrcView &src, const IdxView &idx, const Op &op, const size_t numCols)
static void pack(const DstView &dst, const SrcView &src, const IdxView &idx, size_t col)
Kokkos::View< const Sacado::UQ::PCE< SS > *, SP... > SrcView
execution_space::size_type size_type
DstView::execution_space execution_space
KOKKOS_INLINE_FUNCTION constexpr std::enable_if< is_view_uq_pce< View< T, P... > >::value, unsigned >::type dimension_scalar(const View< T, P... > &view)
KOKKOS_INLINE_FUNCTION void operator()(const size_type k) const
KOKKOS_INLINE_FUNCTION void operator()(const size_type k) const
execution_space::size_type size_type
KOKKOS_INLINE_FUNCTION void operator()(const size_type k) const
KOKKOS_INLINE_FUNCTION void operator()(const size_type k, const size_type tidx) const
KOKKOS_INLINE_FUNCTION void operator()(const size_type k) const
Team-based parallel work configuration for Sacado::MP::Vector.
execution_space::size_type size_type
UnpackArrayMultiColumn(const ExecutionSpace &, const DstView &dst_, const SrcView &src_, const IdxView &idx_, const Op &op_, const size_t numCols_)
KOKKOS_INLINE_FUNCTION void operator()(const size_type k, const size_type tidx) const
Kokkos::View< const Sacado::UQ::PCE< SS > **, SP... > SrcView
PackArrayMultiColumnVariableStride(const DstView &dst_, const SrcView &src_, const IdxView &idx_, const ColView &col_, size_t numCols_)
DstView::execution_space execution_space
execution_space::size_type size_type
DstView::execution_space execution_space
KOKKOS_INLINE_FUNCTION void operator()(const size_type k, const size_type tidx) const
PackArraySingleColumn(const DstView &dst_, const SrcView &src_, const IdxView &idx_, size_t col_)
Kokkos::View< Sacado::UQ::PCE< DS > **, DP... > DstView
Kokkos::View< const Sacado::UQ::PCE< SS > **, SP... > SrcView
Kokkos::View< Sacado::UQ::PCE< DS > *, DP... > DstView
KOKKOS_INLINE_FUNCTION void operator()(const size_type k, const size_type tidx) const
UnpackArrayMultiColumnVariableStride(const ExecutionSpace &, const DstView &dst_, const SrcView &src_, const IdxView &idx_, const ColView &col_, const Op &op_, size_t numCols_)
execution_space::size_type size_type
execution_space::size_type size_type
PackArrayMultiColumn(const DstView &dst_, const SrcView &src_, const IdxView &idx_, size_t numCols_)