42#ifndef TPETRA_DETAILS_COMPUTEOFFSETS_HPP
43#define TPETRA_DETAILS_COMPUTEOFFSETS_HPP
52#include "TpetraCore_config.h"
75template<
class OffsetType,
78class ComputeOffsetsFromCounts {
80 static_assert (std::is_integral<OffsetType>::value,
81 "OffsetType must be a built-in integer.");
82 static_assert (std::is_integral<CountType>::value,
83 "CountType must be a built-in integer.");
84 static_assert (std::is_integral<SizeType>::value,
85 "SizeType must be a built-in integer.");
87 using offsets_view_type =
88 Kokkos::View<OffsetType*, Kokkos::AnonymousSpace>;
89 using counts_view_type =
90 Kokkos::View<const CountType*, Kokkos::AnonymousSpace>;
97 ComputeOffsetsFromCounts (
const offsets_view_type& offsets,
98 const counts_view_type& counts) :
101 size_ (counts.extent (0))
105 KOKKOS_INLINE_FUNCTION
void
106 operator () (
const SizeType i, OffsetType& update,
107 const bool finalPass)
const
109 const auto curVal = (i < size_) ? counts_[i] : OffsetType ();
111 offsets_[i] = update;
113 update += (i < size_) ? curVal : OffsetType ();
116 template<
class ExecutionSpace>
118 run (
const ExecutionSpace& execSpace,
119 const offsets_view_type& offsets,
120 const counts_view_type& counts)
122 const SizeType numCounts (counts.extent (0));
123 using range_type = Kokkos::RangePolicy<ExecutionSpace, SizeType>;
124 range_type range (execSpace, 0, numCounts + SizeType (1));
126 ComputeOffsetsFromCounts<OffsetType, CountType, SizeType>;
127 functor_type functor (offsets, counts);
128 OffsetType total (0);
129 const char funcName[] =
"Tpetra::Details::computeOffsetsFromCounts";
130 Kokkos::parallel_scan (funcName, range, functor, total);
136 offsets_view_type offsets_;
138 counts_view_type counts_;
152template<
class OffsetType,
155class ComputeOffsetsFromConstantCount {
157 static_assert (std::is_integral<OffsetType>::value,
158 "OffsetType must be a built-in integer.");
159 static_assert (std::is_integral<CountType>::value,
160 "CountType must be a built-in integer.");
161 static_assert (std::is_integral<SizeType>::value,
162 "SizeType must be a built-in integer.");
164 using offsets_view_type =
165 Kokkos::View<OffsetType*, Kokkos::AnonymousSpace>;
172 ComputeOffsetsFromConstantCount (
const offsets_view_type& offsets,
173 const CountType count) :
179 KOKKOS_INLINE_FUNCTION
void
180 operator () (
const SizeType i)
const
182 offsets_[i] = count_*i;
185 template<
class ExecutionSpace>
187 run (
const ExecutionSpace& execSpace,
188 const offsets_view_type& offsets,
189 const CountType count)
191 const SizeType numOffsets (offsets.extent (0));
192 using range_type = Kokkos::RangePolicy<ExecutionSpace, SizeType>;
193 range_type range (execSpace, 0, numOffsets);
195 ComputeOffsetsFromConstantCount<OffsetType, CountType, SizeType>;
196 functor_type functor (offsets, count);
197 const OffsetType total = numOffsets*count;
198 const char funcName[] =
199 "Tpetra::Details::computeOffsetsFromConstantCount";
200 Kokkos::parallel_for (funcName, range, functor);
206 offsets_view_type offsets_;
238template<
class ExecutionSpace,
239 class OffsetsViewType,
240 class CountsViewType,
241 class SizeType =
typename OffsetsViewType::size_type>
242typename OffsetsViewType::non_const_value_type
247 static_assert (Kokkos::is_execution_space<ExecutionSpace>::value,
248 "ExecutionSpace must be a Kokkos execution space.");
249 static_assert (Kokkos::is_view<OffsetsViewType>::value,
250 "OffsetsViewType (the type of ptr) must be a Kokkos::View.");
251 static_assert (Kokkos::is_view<CountsViewType>::value,
252 "CountsViewType (the type of counts) must be a Kokkos::View.");
253 static_assert (std::is_same<
typename OffsetsViewType::value_type,
254 typename OffsetsViewType::non_const_value_type>::value,
255 "OffsetsViewType (the type of ptr) must be a nonconst Kokkos::View.");
256 static_assert (
static_cast<int> (OffsetsViewType::rank) == 1,
257 "OffsetsViewType (the type of ptr) must be a rank-1 Kokkos::View.");
258 static_assert (
static_cast<int> (CountsViewType::rank) == 1,
259 "CountsViewType (the type of counts) must be a rank-1 Kokkos::View.");
261 using offset_type =
typename OffsetsViewType::non_const_value_type;
262 static_assert (std::is_integral<offset_type>::value,
263 "The entries of ptr must be built-in integers.");
264 using count_type =
typename CountsViewType::non_const_value_type;
265 static_assert (std::is_integral<count_type>::value,
266 "The entries of counts must be built-in integers.");
267 static_assert (std::is_integral<SizeType>::value,
268 "SizeType must be a built-in integer type.");
270 const char funcName[] =
"Tpetra::Details::computeOffsetsFromCounts";
274 offset_type
total (0);
279 ": counts.size() = " <<
numCounts <<
" >= ptr.size() = " <<
282 using Kokkos::AnonymousSpace;
292 typename offsets_device_type::memory_space;
295 Kokkos::SpaceAccessibility<
305 using Kokkos::view_alloc;
306 using Kokkos::WithoutInitializing;
324 class SizeType =
typename OffsetsViewType::size_type>
325typename OffsetsViewType::non_const_value_type
329 using execution_space =
typename OffsetsViewType::execution_space;
357 class SizeType =
typename OffsetsViewType::size_type>
358typename OffsetsViewType::non_const_value_type
362 static_assert (Kokkos::is_view<OffsetsViewType>::value,
363 "ptr must be a Kokkos::View.");
364 static_assert (std::is_same<
typename OffsetsViewType::value_type,
365 typename OffsetsViewType::non_const_value_type>::value,
366 "ptr must be a nonconst Kokkos::View.");
367 static_assert (
static_cast<int> (OffsetsViewType::rank) == 1,
368 "ptr must be a rank-1 Kokkos::View.");
370 using offset_type =
typename OffsetsViewType::non_const_value_type;
371 static_assert (std::is_integral<offset_type>::value,
372 "The type of each entry of ptr must be a "
373 "built-in integer.");
374 static_assert (std::is_integral<CountType>::value,
375 "CountType must be a built-in integer.");
376 static_assert (std::is_integral<SizeType>::value,
377 "SizeType must be a built-in integer.");
379 using device_type =
typename OffsetsViewType::device_type;
380 using execution_space =
typename device_type::execution_space;
382 offset_type
total (0);
383 if (
ptr.extent (0) != 0) {
Declaration and definition of Tpetra::Details::getEntryOnHost.
Struct that holds views of the contents of a CrsMatrix.
Implementation details of Tpetra.
OffsetsViewType::non_const_value_type computeOffsetsFromCounts(const ExecutionSpace &execSpace, const OffsetsViewType &ptr, const CountsViewType &counts)
Compute offsets from counts.
OffsetsViewType::non_const_value_type computeOffsetsFromConstantCount(const OffsetsViewType &ptr, const CountType count)
Compute offsets from a constant count.
Namespace Tpetra contains the class and methods constituting the Tpetra library.