Tpetra parallel linear algebra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Tpetra_Details_getGraphDiagOffsets_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Tpetra: Templated Linear Algebra Services Package
4 //
5 // Copyright 2008 NTESS and the Tpetra contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef TPETRA_DETAILS_GETGRAPHDIAGOFFSETS_DECL_HPP
11 #define TPETRA_DETAILS_GETGRAPHDIAGOFFSETS_DECL_HPP
12 
17 
18 #include "TpetraCore_config.h"
19 #include "Kokkos_Core.hpp"
20 #include "KokkosSparse_StaticCrsGraph.hpp"
22 #include <type_traits>
23 
24 namespace Tpetra {
25 namespace Details {
26 namespace Impl {
27 
41 template<class LO,
42  class GO,
43  class DeviceType,
44  class DiagOffsetType = size_t>
46 public:
47  typedef typename DeviceType::device_type device_type;
48  typedef DiagOffsetType diag_offset_type;
49  typedef ::Kokkos::View<diag_offset_type*,
50  device_type,
51  ::Kokkos::MemoryUnmanaged> diag_offsets_type;
52  using local_graph_device_type = ::KokkosSparse::StaticCrsGraph<LO,
53  ::Kokkos::LayoutLeft,
54  device_type, void, size_t>;
55  typedef ::Tpetra::Details::LocalMap<LO, GO, device_type> local_map_type;
56  typedef ::Kokkos::View<const typename local_graph_device_type::size_type*,
57  ::Kokkos::LayoutLeft,
58  device_type,
59  ::Kokkos::MemoryUnmanaged> row_offsets_type;
60  // This is unmanaged for performance, because we need to take
61  // subviews inside the functor.
62  typedef ::Kokkos::View<const LO*,
63  ::Kokkos::LayoutLeft,
64  device_type,
65  ::Kokkos::MemoryUnmanaged> lcl_col_inds_type;
66 
68  GetGraphDiagOffsets (const diag_offsets_type& diagOffsets,
69  const local_map_type& lclRowMap,
70  const local_map_type& lclColMap,
71  const row_offsets_type& ptr,
72  const lcl_col_inds_type& ind,
73  const bool isSorted);
74 
76  KOKKOS_FUNCTION void operator() (const LO& lclRowInd) const;
77 
78 private:
79  diag_offsets_type diagOffsets_;
80  local_map_type lclRowMap_;
81  local_map_type lclColMap_;
82  row_offsets_type ptr_;
83  lcl_col_inds_type ind_;
84  bool isSorted_;
85 };
86 
87 } // namespace Impl
88 
89 template<class DiagOffsetsType,
90  class LclMapType,
91  class RowOffsetsType,
92  class LclColIndsType>
93 void
94 getGraphDiagOffsets (const DiagOffsetsType& diagOffsets,
95  const LclMapType& lclRowMap,
96  const LclMapType& lclColMap,
97  const RowOffsetsType& ptr,
98  const LclColIndsType& ind,
99  const bool isSorted)
100 {
101  static_assert (Kokkos::is_view<DiagOffsetsType>::value,
102  "DiagOffsetsType (the type of diagOffsets) must be a Kokkos::View.");
103  static_assert (Kokkos::is_view<RowOffsetsType>::value,
104  "RowOffsetsType (the type of ptr) must be a Kokkos::View.");
105  static_assert (Kokkos::is_view<LclColIndsType>::value,
106  "LclColIndsType (the type of ind) must be a Kokkos::View.");
107  static_assert (static_cast<int> (DiagOffsetsType::rank) == 1,
108  "DiagOffsetsType (the type of diagOffsets) must be a rank-1 Kokkos::View.");
109  static_assert (static_cast<int> (RowOffsetsType::rank) == 1,
110  "RowOffsetsType (the type of ptr) must be a rank-1 Kokkos::View.");
111  static_assert (static_cast<int> (LclColIndsType::rank) == 1,
112  "LclColIndsType (the type of ind) must be a rank-1 Kokkos::View.");
113  typedef typename DiagOffsetsType::non_const_value_type diag_offset_type;
114  static_assert (std::is_same<typename DiagOffsetsType::value_type, diag_offset_type>::value,
115  "DiagOffsetsType (the type of diagOffsets) must be a nonconst "
116  "Kokkos::View, since it is the output argument of this function.");
117  static_assert (std::is_integral<diag_offset_type>::value,
118  "The type of each entry of diagOffsets must be an integer.");
119  typedef typename LclColIndsType::non_const_value_type local_ordinal_type;
120  static_assert (std::is_integral<local_ordinal_type>::value,
121  "The type of each entry of ind (the array of column indices) "
122  "must be an integer.");
123  static_assert (sizeof (diag_offset_type) >= sizeof (local_ordinal_type),
124  "Diagonal offset type must be big enough to count the number "
125  "of column indices, since the diagonal entry (if it exists) "
126  "may be anywhere in a row.");
127  typedef typename RowOffsetsType::non_const_value_type row_offset_type;
128  static_assert (std::is_integral<row_offset_type>::value,
129  "The type of each entry of ptr must be an integer.");
130 
131  typedef typename LclMapType::local_ordinal_type LO;
132  typedef typename LclMapType::global_ordinal_type GO;
133  typedef typename LclMapType::device_type DT;
134 
136  // The functor's constructor runs the functor.
137  impl_type impl (diagOffsets, lclRowMap, lclColMap, ptr, ind, isSorted);
138 }
139 
140 } // namespace Details
141 } // namespace Tpetra
142 
143 #endif // TPETRA_DETAILS_GETGRAPHDIAGOFFSETS_DECL_HPP
KOKKOS_FUNCTION void operator()(const LO &lclRowInd) const
Kokkos::parallel_for loop body.
Declaration and definition of the Tpetra::Map class, an implementation detail of Tpetra::Map.
Implementation detail of Tpetra::Details::getGraphDiagOffsets, which in turn is an implementation det...
GetGraphDiagOffsets(const diag_offsets_type &diagOffsets, const local_map_type &lclRowMap, const local_map_type &lclColMap, const row_offsets_type &ptr, const lcl_col_inds_type &ind, const bool isSorted)
Constructor; also runs the functor.