Tpetra parallel linear algebra  Version of the Day
Tpetra_Details_computeOffsets.hpp
Go to the documentation of this file.
1 /*
2 // @HEADER
3 // ***********************************************************************
4 //
5 // Tpetra: Templated Linear Algebra Services Package
6 // Copyright (2008) Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39 //
40 // ************************************************************************
41 // @HEADER
42 */
43 
44 #ifndef TPETRA_DETAILS_COMPUTEOFFSETS_HPP
45 #define TPETRA_DETAILS_COMPUTEOFFSETS_HPP
46 
52 
53 #include "TpetraCore_config.h"
54 #include "Kokkos_Core.hpp"
55 #include <limits>
56 #include <type_traits>
57 
58 namespace Tpetra {
59 namespace Details {
60 
61 //
62 // Implementation details for computeOffsetsFromCounts (see below).
63 // Users should skip over this anonymous namespace.
64 //
65 namespace { // (anonymous)
66 
84 template<class OffsetsViewType,
85  class CountsViewType,
86  class SizeType = typename OffsetsViewType::size_type>
87 class ComputeOffsetsFromCounts {
88 public:
89  static_assert (Kokkos::Impl::is_view<OffsetsViewType>::value,
90  "OffsetsViewType (the type of ptr) must be a Kokkos::View.");
91  static_assert (Kokkos::Impl::is_view<CountsViewType>::value,
92  "CountsViewType (the type of counts) must be a Kokkos::View.");
93  static_assert (std::is_same<typename OffsetsViewType::value_type,
94  typename OffsetsViewType::non_const_value_type>::value,
95  "OffsetsViewType (the type of ptr) must be a nonconst Kokkos::View.");
96  static_assert (static_cast<int> (OffsetsViewType::rank) == 1,
97  "OffsetsViewType (the type of ptr) must be a rank-1 Kokkos::View.");
98  static_assert (static_cast<int> (CountsViewType::rank) == 1,
99  "CountsViewType (the type of counts) must be a rank-1 Kokkos::View.");
100  static_assert (std::is_integral<typename OffsetsViewType::non_const_value_type>::value,
101  "The entries of ptr must be built-in integers.");
102  static_assert (std::is_integral<typename CountsViewType::non_const_value_type>::value,
103  "The entries of counts must be built-in integers.");
104  static_assert (std::is_integral<SizeType>::value,
105  "SizeType must be a built-in integer type.");
106 
107  typedef OffsetsViewType offsets_view_type;
108  typedef typename CountsViewType::const_type counts_view_type;
109  typedef SizeType size_type;
110  typedef typename OffsetsViewType::non_const_value_type value_type;
111 
117  ComputeOffsetsFromCounts (const offsets_view_type& offsets,
118  const counts_view_type& counts) :
119  offsets_ (offsets),
120  counts_ (counts),
121  size_ (counts.dimension_0 ())
122  {}
123 
125  KOKKOS_INLINE_FUNCTION void init (value_type& dst) const
126  {
127  dst = 0;
128  }
129 
131  KOKKOS_INLINE_FUNCTION void
132  join (volatile value_type& dst,
133  const volatile value_type& src) const
134  {
135  dst += src;
136  }
137 
139  KOKKOS_INLINE_FUNCTION void
140  operator () (const size_type& i, value_type& update, const bool final) const
141  {
142  if (final) {
143  offsets_[i] = update;
144  }
145  if (i < size_) {
146  update += counts_[i];
147  }
148  }
149 
150 private:
152  offsets_view_type offsets_;
154  counts_view_type counts_;
156  size_type size_;
157 };
158 
176 template<class OffsetsViewType,
177  class CountType,
178  class SizeType = typename OffsetsViewType::size_type>
179 class ComputeOffsetsFromConstantCount {
180 public:
181  static_assert (Kokkos::Impl::is_view<OffsetsViewType>::value,
182  "OffsetsViewType (the type of ptr) must be a Kokkos::View.");
183  static_assert (std::is_same<typename OffsetsViewType::value_type,
184  typename OffsetsViewType::non_const_value_type>::value,
185  "OffsetsViewType (the type of ptr) must be a nonconst Kokkos::View.");
186  static_assert (static_cast<int> (OffsetsViewType::rank) == 1,
187  "OffsetsViewType (the type of ptr) must be a rank-1 Kokkos::View.");
188  static_assert (std::is_integral<typename OffsetsViewType::non_const_value_type>::value,
189  "The entries of ptr must be built-in integers.");
190  static_assert (std::is_integral<CountType>::value,
191  "CountType must be a built-in integer type.");
192  static_assert (std::is_integral<SizeType>::value,
193  "SizeType must be a built-in integer type.");
194 
195  typedef OffsetsViewType offsets_view_type;
196  typedef CountType count_type;
197  typedef SizeType size_type;
198  typedef typename offsets_view_type::non_const_value_type value_type;
199 
205  ComputeOffsetsFromConstantCount (const offsets_view_type& offsets,
206  const count_type count) :
207  offsets_ (offsets),
208  count_ (count),
209  size_ (offsets_.dimension_0 () == 0 ?
210  static_cast<size_type> (0) :
211  static_cast<size_type> (offsets_.dimension_0 () - 1))
212  {}
213 
215  KOKKOS_INLINE_FUNCTION void init (value_type& dst) const
216  {
217  dst = 0;
218  }
219 
221  KOKKOS_INLINE_FUNCTION void
222  join (volatile value_type& dst,
223  const volatile value_type& src) const
224  {
225  dst += src;
226  }
227 
229  KOKKOS_INLINE_FUNCTION void
230  operator () (const size_type& i, value_type& update, const bool final) const
231  {
232  if (final) {
233  offsets_[i] = update;
234  }
235  if (i < size_) {
236  update += count_;
237  }
238  }
239 
240 private:
242  offsets_view_type offsets_;
244  count_type count_;
246  size_type size_;
247 };
248 
249 } // namespace (anonymous)
250 
274 template<class OffsetsViewType,
275  class CountsViewType,
276  class SizeType = typename OffsetsViewType::size_type>
277 typename OffsetsViewType::non_const_value_type
278 computeOffsetsFromCounts (const OffsetsViewType& ptr,
279  const CountsViewType& counts)
280 {
281  static_assert (Kokkos::Impl::is_view<OffsetsViewType>::value,
282  "OffsetsViewType (the type of ptr) must be a Kokkos::View.");
283  static_assert (Kokkos::Impl::is_view<CountsViewType>::value,
284  "CountsViewType (the type of counts) must be a Kokkos::View.");
285  static_assert (std::is_same<typename OffsetsViewType::value_type,
286  typename OffsetsViewType::non_const_value_type>::value,
287  "OffsetsViewType (the type of ptr) must be a nonconst Kokkos::View.");
288  static_assert (static_cast<int> (OffsetsViewType::rank) == 1,
289  "OffsetsViewType (the type of ptr) must be a rank-1 Kokkos::View.");
290  static_assert (static_cast<int> (CountsViewType::rank) == 1,
291  "CountsViewType (the type of counts) must be a rank-1 Kokkos::View.");
292  static_assert (std::is_integral<typename OffsetsViewType::non_const_value_type>::value,
293  "The entries of ptr must be built-in integers.");
294  static_assert (std::is_integral<typename CountsViewType::non_const_value_type>::value,
295  "The entries of counts must be built-in integers.");
296  static_assert (std::is_integral<SizeType>::value,
297  "SizeType must be a built-in integer type.");
298 
299  typedef typename CountsViewType::non_const_value_type count_type;
300  typedef typename OffsetsViewType::non_const_value_type offset_type;
301  typedef typename OffsetsViewType::device_type device_type;
302  typedef typename device_type::execution_space execution_space;
303  typedef typename device_type::memory_space memory_space;
304 
305  const auto numOffsets = ptr.size ();
306  const auto numCounts = counts.size ();
307  if (numOffsets != 0) {
308  TEUCHOS_TEST_FOR_EXCEPTION
309  (numCounts >= numOffsets, std::invalid_argument,
310  "computeOffsetsFromCounts: counts.dimension_0() = " << numCounts
311  << " >= ptr.dimension_0() = " << numOffsets << ".");
312 
313  Kokkos::RangePolicy<execution_space, SizeType> range (0, numCounts+1);
314  try {
315  // We always want to run in the offsets' execution space, since
316  // that is the output argument. (This gives us first touch, if
317  // applicable, and in general improves locality.) However, we
318  // need to make sure that we can access counts from this
319  // execution space. If we can't, we need to make a temporary
320  // "device" copy of counts in offsets' memory space.
321 
322  // The first template parameter needs to be a memory space.
323  constexpr bool countsAccessibleFromOffsets =
324  Kokkos::Impl::VerifyExecutionCanAccessMemorySpace<memory_space,
325  typename CountsViewType::memory_space>::value;
326  if (countsAccessibleFromOffsets) {
327  typedef ComputeOffsetsFromCounts<OffsetsViewType, CountsViewType,
328  SizeType> functor_type;
329  // offsets' execution space can access counts
330  functor_type functor (ptr, counts);
331  Kokkos::parallel_scan (range, functor);
332  }
333  else {
334  // Make a temporary copy of counts in offsets' execution
335  // space. Use the same array layout as the original, so we
336  // can deep copy.
337  typedef Kokkos::View<count_type*, typename CountsViewType::array_layout,
338  device_type> dev_counts_type;
339  dev_counts_type counts_d ("counts_d", numCounts);
340  Kokkos::deep_copy (counts_d, counts);
341 
342  typedef ComputeOffsetsFromCounts<OffsetsViewType, dev_counts_type,
343  SizeType> functor_type;
344  functor_type functor (ptr, counts_d);
345  Kokkos::parallel_scan (range, functor);
346  }
347  }
348  catch (std::exception& e) {
349  TEUCHOS_TEST_FOR_EXCEPTION
350  (true, std::runtime_error, "computeOffsetsFromCounts: parallel_scan "
351  "(with device_type Kokkos::Device<" <<
352  typeid (execution_space).name () << ", " <<
353  typeid (memory_space).name () << ">) threw an std::exception: "
354  << e.what ());
355  }
356  catch (...) {
357  TEUCHOS_TEST_FOR_EXCEPTION
358  (true, std::runtime_error, "Kokkos::parallel_scan threw an "
359  "exception not a subclass of std::exception");
360  }
361 
362  // Get the sum of all the entries of counts from the last entry of
363  // ptr. The second branch of this 'if' always works, but we save
364  // a little time by specializing for non-CUDA execution spaces.
365  if (Kokkos::Impl::VerifyExecutionCanAccessMemorySpace<Kokkos::HostSpace,
366  memory_space>::value) {
367  return ptr[numCounts];
368  }
369  else {
370  auto ptr_last = Kokkos::subview (ptr, numCounts);
371  auto ptr_last_h = Kokkos::create_mirror_view (ptr_last);
372  Kokkos::deep_copy (ptr_last_h, ptr_last);
373  return ptr_last_h ();
374  }
375  }
376  else {
377  return static_cast<offset_type> (0);
378  }
379 }
380 
403 template<class OffsetsViewType,
404  class CountType,
405  class SizeType = typename OffsetsViewType::size_type>
406 typename OffsetsViewType::non_const_value_type
407 computeOffsetsFromConstantCount (const OffsetsViewType& ptr,
408  const CountType& count)
409 {
410  static_assert (Kokkos::Impl::is_view<OffsetsViewType>::value,
411  "OffsetsViewType (the type of ptr) must be a Kokkos::View.");
412  static_assert (std::is_same<typename OffsetsViewType::value_type,
413  typename OffsetsViewType::non_const_value_type>::value,
414  "OffsetsViewType (the type of ptr) must be a nonconst Kokkos::View.");
415  static_assert (static_cast<int> (OffsetsViewType::rank) == 1,
416  "OffsetsViewType (the type of ptr) must be a rank-1 Kokkos::View.");
417  static_assert (std::is_integral<typename OffsetsViewType::non_const_value_type>::value,
418  "The entries of ptr must be built-in integers.");
419  static_assert (std::is_integral<CountType>::value,
420  "CountType must be a built-in integer type.");
421  static_assert (std::is_integral<SizeType>::value,
422  "SizeType must be a built-in integer type.");
423 
424  typedef typename std::decay<CountType>::type count_type;
425  typedef typename OffsetsViewType::non_const_value_type offset_type;
426  typedef typename OffsetsViewType::device_type device_type;
427  typedef typename device_type::execution_space execution_space;
428  typedef typename device_type::memory_space memory_space;
429 
430  const auto numOffsets = ptr.size ();
431  if (numOffsets != 0) {
432  ComputeOffsetsFromConstantCount<OffsetsViewType, count_type,
433  SizeType> functor (ptr, count);
434  Kokkos::RangePolicy<execution_space, SizeType> range (0, numOffsets);
435  try {
436 #ifdef KOKKOS_HAVE_CUDA
437  constexpr bool isCuda = std::is_same<execution_space, Kokkos::Cuda>::value;
438 #else
439  constexpr bool isCuda = false;
440 #endif // KOKKOS_HAVE_CUDA
441  // FIXME (mfh 26 Jun 2016) Kokkos::Cuda doesn't like explicit
442  // use of RangePolicy here; it reports a build error. Not sure
443  // why. Using an integer range works just fine.
444  if (isCuda) {
445  Kokkos::parallel_scan (numOffsets, functor);
446  }
447  else {
448  Kokkos::parallel_scan (range, functor);
449  }
450  }
451  catch (std::exception& e) {
452  TEUCHOS_TEST_FOR_EXCEPTION
453  (true, std::runtime_error, "computeOffsetsFromConstantCount: "
454  "parallel_scan (with device_type Kokkos::Device<" <<
455  typeid (execution_space).name () << ", " <<
456  typeid (memory_space).name () << ">) threw an std::exception: "
457  << e.what ());
458  }
459  catch (...) {
460  TEUCHOS_TEST_FOR_EXCEPTION
461  (true, std::runtime_error, "Kokkos::parallel_scan threw an "
462  "exception not a subclass of std::exception");
463  }
464 
465  // Get the sum of all the entries of counts from the last entry of
466  // ptr. The second branch of this 'if' always works, but we save
467  // a little time by specializing for non-CUDA execution spaces.
468  if (Kokkos::Impl::VerifyExecutionCanAccessMemorySpace<Kokkos::HostSpace,
469  memory_space>::value) {
470  return ptr[numOffsets - 1];
471  }
472  else {
473  auto ptr_last = Kokkos::subview (ptr, numOffsets - 1);
474  auto ptr_last_h = Kokkos::create_mirror_view (ptr_last);
475  Kokkos::deep_copy (ptr_last_h, ptr_last);
476  return ptr_last_h ();
477  }
478  }
479  else {
480  return static_cast<offset_type> (0);
481  }
482 }
483 
484 } // namespace Details
485 } // namespace Tpetra
486 
487 #endif // TPETRA_DETAILS_COMPUTEOFFSETS_HPP
Namespace Tpetra contains the class and methods constituting the Tpetra library.
void deep_copy(MultiVector< DS, DL, DG, DN, dstClassic > &dst, const MultiVector< SS, SL, SG, SN, srcClassic > &src)
Copy the contents of the MultiVector src into dst.
Implementation details of Tpetra.
OffsetsViewType::non_const_value_type computeOffsetsFromCounts(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.