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<
69 Kokkos::View<Sacado::UQ::PCE<DS>*,DP...>,
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<
120 Kokkos::View<Sacado::UQ::PCE<DS>*,DP...>,
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<
179 Kokkos::View<Sacado::UQ::PCE<DS>*,DP...>,
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 DS,
typename ... DP,
239 typename SS,
typename ... SP,
240 typename IdxView,
typename Op>
241 struct UnpackArrayMultiColumn<
242 Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
243 Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>,
247 typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>
DstView;
248 typedef Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>
SrcView;
252 static const unsigned BlockSize = 32;
265 dst(dst_), src(src_), idx(idx_), op(op_), numCols(numCols_) {}
267 KOKKOS_INLINE_FUNCTION
270 const size_t offset = k*numCols;
271 for (
size_t j = 0;
j < numCols; ++
j)
272 op( dst(localRow,
j), src(offset+
j) );
275 KOKKOS_INLINE_FUNCTION
278 const size_t offset = k*numCols;
279 for (
size_t j = 0;
j < numCols; ++
j)
281 op( dst(localRow,
j).fastAccessCoeff(i),
282 src(offset+
j).fastAccessCoeff(i) );
291 Kokkos::parallel_for(
293 UnpackArrayMultiColumn(dst,src,idx,op,numCols) );
295 Kokkos::parallel_for( idx.size(),
296 UnpackArrayMultiColumn(dst,src,idx,op,numCols) );
300 template <
typename DS,
typename ... DP,
301 typename SS,
typename ... SP,
302 typename IdxView,
typename ColView,
typename Op>
303 struct UnpackArrayMultiColumnVariableStride<
304 Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
305 Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>,
310 typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>
DstView;
311 typedef Kokkos::View<const Sacado::UQ::PCE<SS>*,SP...>
SrcView;
315 static const unsigned BlockSize = 32;
330 dst(dst_), src(src_), idx(idx_), col(col_), op(op_), numCols(numCols_) {}
332 KOKKOS_INLINE_FUNCTION
335 const size_t offset = k*numCols;
336 for (
size_t j = 0;
j < numCols; ++
j)
337 op( dst(localRow,col(
j)), src(offset+
j) );
340 KOKKOS_INLINE_FUNCTION
343 const size_t offset = k*numCols;
344 for (
size_t j = 0;
j < numCols; ++
j)
346 op( dst(localRow,col(
j)).fastAccessCoeff(i),
347 src(offset+
j).fastAccessCoeff(i) );
357 Kokkos::parallel_for(
359 UnpackArrayMultiColumnVariableStride(dst,src,idx,col,op,numCols) );
361 Kokkos::parallel_for( idx.size(),
362 UnpackArrayMultiColumnVariableStride(
363 dst,src,idx,col,op,numCols) );
367 template <
typename DS,
typename ... DP,
368 typename SS,
typename ... SP,
369 typename DstIdxView,
typename SrcIdxView>
370 struct PermuteArrayMultiColumn<
371 Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
372 Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
376 typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>
DstView;
377 typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>
SrcView;
381 static const unsigned BlockSize = 32;
391 const DstIdxView& dst_idx_,
392 const SrcIdxView& src_idx_,
394 dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
397 KOKKOS_INLINE_FUNCTION
401 for (
size_t j = 0;
j < numCols; ++
j)
402 dst(toRow,
j) = src(fromRow,
j);
405 KOKKOS_INLINE_FUNCTION
409 for (
size_t j = 0;
j < numCols; ++
j)
411 dst(toRow,
j).fastAccessCoeff(i) =
412 src(fromRow,
j).fastAccessCoeff(i);
417 const DstIdxView& dst_idx,
418 const SrcIdxView& src_idx,
422 Kokkos::parallel_for(
424 PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
426 Kokkos::parallel_for(
427 n, PermuteArrayMultiColumn(dst,src,dst_idx,src_idx,numCols) );
431 template <
typename DS,
typename ... DP,
432 typename SS,
typename ... SP,
433 typename DstIdxView,
typename SrcIdxView,
434 typename DstColView,
typename SrcColView>
435 struct PermuteArrayMultiColumnVariableStride<
436 Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>,
437 Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>,
438 DstIdxView, SrcIdxView,
439 DstColView, SrcColView >
441 typedef Kokkos::View<Sacado::UQ::PCE<DS>**,DP...>
DstView;
442 typedef Kokkos::View<const Sacado::UQ::PCE<SS>**,SP...>
SrcView;
446 static const unsigned BlockSize = 32;
458 const DstIdxView& dst_idx_,
459 const SrcIdxView& src_idx_,
460 const DstColView& dst_col_,
461 const SrcColView& src_col_,
463 dst(dst_), src(src_), dst_idx(dst_idx_), src_idx(src_idx_),
464 dst_col(dst_col_), src_col(src_col_),
467 KOKKOS_INLINE_FUNCTION
471 for (
size_t j = 0;
j < numCols; ++
j)
472 dst(toRow, dst_col(
j)) = src(fromRow, src_col(
j));
475 KOKKOS_INLINE_FUNCTION
479 for (
size_t j = 0;
j < numCols; ++
j)
481 dst(toRow, dst_col(
j)).fastAccessCoeff(i) =
482 src(fromRow, src_col(
j)).fastAccessCoeff(i);
487 const DstIdxView& dst_idx,
488 const SrcIdxView& src_idx,
489 const DstColView& dst_col,
490 const SrcColView& src_col,
494 Kokkos::parallel_for(
496 PermuteArrayMultiColumnVariableStride(
497 dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
499 Kokkos::parallel_for(
500 n, PermuteArrayMultiColumnVariableStride(
501 dst,src,dst_idx,src_idx,dst_col,src_col,numCols) );
509 #endif // TPETRA_KOKKOS_REFACTOR_DETAILS_MULTI_VECTOR_DIST_OBJECT_KERNELS_UQ_PCE_HPP
static void pack(const DstView &dst, const SrcView &src, const IdxView &idx, const ColView &col, size_t numCols)
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
Kokkos::View< const Sacado::UQ::PCE< SS > *, SP... > SrcView
Kokkos::View< const Sacado::UQ::PCE< SS > **, SP... > SrcView
KOKKOS_INLINE_FUNCTION void operator()(const size_type k, const size_type tidx) const
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
KOKKOS_INLINE_FUNCTION void operator()(const size_type k, const size_type tidx) const
static void unpack(const DstView &dst, const SrcView &src, const IdxView &idx, const ColView &col, const Op &op, size_t numCols)
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 size_type tidx) const
KOKKOS_INLINE_FUNCTION void operator()(const size_type k) const
UnpackArrayMultiColumnVariableStride(const DstView &dst_, const SrcView &src_, const IdxView &idx_, const ColView &col_, const Op &op_, size_t numCols_)
DstView::execution_space execution_space
Kokkos::View< Sacado::UQ::PCE< DS > *, DP... > DstView
DstView::execution_space execution_space
UnpackArrayMultiColumn(const DstView &dst_, const SrcView &src_, const IdxView &idx_, const Op &op_, size_t numCols_)
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_)
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 pack(const DstView &dst, const SrcView &src, const IdxView &idx, size_t col)
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.
static void unpack(const DstView &dst, const SrcView &src, const IdxView &idx, const Op &op, size_t numCols)
execution_space::size_type size_type
KOKKOS_INLINE_FUNCTION void operator()(const size_type k) const
DstView::execution_space execution_space
DstView::execution_space execution_space
Kokkos::View< Sacado::UQ::PCE< DS > **, DP... > DstView
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< 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
Kokkos::View< Sacado::UQ::PCE< DS > **, DP... > DstView
execution_space::size_type size_type
DstView::execution_space execution_space
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
KOKKOS_INLINE_FUNCTION void operator()(const size_type k) const
execution_space::size_type size_type
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_)