40 #ifndef TPETRA_DETAILS_LOCALROWOFFSETS_DEF_HPP
41 #define TPETRA_DETAILS_LOCALROWOFFSETS_DEF_HPP
47 #include "TpetraCore_config.h"
48 #include "Tpetra_CrsGraph.hpp"
49 #include "Tpetra_RowGraph.hpp"
52 #include "Kokkos_Core.hpp"
54 #ifdef TPETRA_ENABLE_DEPRECATED_CODE
61 template <
class LO,
class GO,
class NT>
62 std::pair<typename LocalRowOffsetsResult<NT>::offsets_type,
size_t>
63 localRowCounts (
const RowGraph<LO, GO, NT>& G)
65 using result_type = LocalRowOffsetsResult<NT>;
66 using offsets_type =
typename result_type::offsets_type;
67 using offset_type =
typename result_type::offset_type;
69 const LO lclNumRows (G.getLocalNumRows ());
70 offsets_type entPerRow;
71 if (lclNumRows != 0) {
72 using Kokkos::view_alloc;
73 using Kokkos::WithoutInitializing;
75 offsets_type (view_alloc (
"entPerRow", WithoutInitializing),
78 using host = Kokkos::DefaultHostExecutionSpace;
79 auto entPerRow_h = Kokkos::create_mirror_view (host (), entPerRow);
86 for (LO i = 0; i < lclNumRows; ++i) {
87 const size_t lclNumEnt = G.getNumEntriesInLocalRow (i);
88 entPerRow_h[i] = offset_type (lclNumEnt);
89 maxNumEnt = maxNumEnt < lclNumEnt ? lclNumEnt : maxNumEnt;
92 using execution_space =
typename NT::execution_space;
94 return {entPerRow, maxNumEnt};
97 template <
class LO,
class GO,
class NT>
98 LocalRowOffsetsResult<NT>
99 localRowOffsetsFromRowGraph (
const RowGraph<LO, GO, NT>& G)
101 using result_type = LocalRowOffsetsResult<NT>;
102 using offsets_type =
typename result_type::offsets_type;
103 using offset_type =
typename result_type::offset_type;
105 offsets_type entPerRow;
106 size_t maxNumEnt = 0;
108 auto result = localRowCounts (G);
109 entPerRow = result.first;
110 maxNumEnt = result.second;
113 const LO lclNumRows (G.getLocalNumRows ());
116 if (lclNumRows != 0) {
117 using Kokkos::view_alloc;
118 using Kokkos::WithoutInitializing;
119 ptr = offsets_type (view_alloc (
"ptr", WithoutInitializing),
124 return {ptr, nnz, maxNumEnt};
127 template <
class LO,
class GO,
class NT>
128 LocalRowOffsetsResult<NT>
129 localRowOffsetsFromFillCompleteCrsGraph (
const CrsGraph<LO, GO, NT>& G)
131 using Kokkos::view_alloc;
132 using Kokkos::WithoutInitializing;
133 using result_type = LocalRowOffsetsResult<NT>;
134 using offsets_type =
typename result_type::offsets_type;
135 using offset_type =
typename result_type::offset_type;
137 auto G_lcl = G.getLocalGraphDevice ();
138 offsets_type ptr (view_alloc (
"ptr", WithoutInitializing),
139 G_lcl.row_map.extent (0));
143 const offset_type nnz = G.getLocalNumEntries ();
144 const size_t maxNumEnt = G.getLocalMaxNumRowEntries ();
145 return {ptr, nnz, maxNumEnt};
150 template <
class LO,
class GO,
class NT>
151 LocalRowOffsetsResult<NT>
152 localRowOffsets (
const RowGraph<LO, GO, NT>& G)
154 if (G.isFillComplete ()) {
155 using crs_graph_type = CrsGraph<LO, GO, NT>;
156 const crs_graph_type* G_crs =
157 dynamic_cast<const crs_graph_type*
> (&G);
158 if (G_crs !=
nullptr) {
159 return Impl::localRowOffsetsFromFillCompleteCrsGraph (*G_crs);
162 return Impl::localRowOffsetsFromRowGraph (G);
174 #define TPETRA_DETAILS_LOCALROWOFFSETS_INSTANT(LO, GO, NT) \
175 namespace Details { \
178 template std::pair<LocalRowOffsetsResult<NT>::offsets_type, size_t> \
179 localRowCounts (const RowGraph<LO, GO, NT>& G); \
181 template LocalRowOffsetsResult<NT> \
182 localRowOffsetsFromRowGraph (const RowGraph<LO, GO, NT>& G); \
184 template LocalRowOffsetsResult<NT> \
185 localRowOffsetsFromFillCompleteCrsGraph (const CrsGraph<LO, GO, NT>& G); \
189 template LocalRowOffsetsResult<NT> \
190 localRowOffsets (const RowGraph<LO, GO, NT>& A); \
193 #endif // TPETRA_ENABLE_DEPRECATED_CODE
195 #endif // TPETRA_DETAILS_LOCALROWOFFSETS_DEF_HPP
void deep_copy(MultiVector< DS, DL, DG, DN > &dst, const MultiVector< SS, SL, SG, SN > &src)
Copy the contents of the MultiVector src into dst.
Declare and define the functions Tpetra::Details::computeOffsetsFromCounts and Tpetra::computeOffsets...
OffsetsViewType::non_const_value_type computeOffsetsFromCounts(const ExecutionSpace &execSpace, const OffsetsViewType &ptr, const CountsViewType &counts)
Compute offsets from counts.
Declaration and definition of Tpetra::Details::getEntryOnHost.