48 #ifndef __IFPACK2_CRSARRAYS_DECL_HPP__
49 #define __IFPACK2_CRSARRAYS_DECL_HPP__
51 #include <Tpetra_RowMatrix.hpp>
52 #include <Tpetra_CrsMatrix.hpp>
53 #include <Kokkos_DefaultNode.hpp>
54 #include <KokkosSparse_CrsMatrix.hpp>
55 #include <Ifpack2_LocalFilter.hpp>
56 #include <Ifpack2_ReorderFilter.hpp>
65 template<
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal,
typename Node>
68 typedef typename Node::device_type device_type;
69 typedef typename device_type::execution_space execution_space;
70 typedef Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> TRowMatrix;
71 typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> TCrsMatrix;
74 typedef KokkosSparse::CrsMatrix<Scalar, LocalOrdinal, execution_space> KCrsMatrix;
75 typedef Kokkos::View<LocalOrdinal*, execution_space> OrdinalArray;
76 typedef Kokkos::View<Scalar*, execution_space> ScalarArray;
77 typedef Kokkos::View<LocalOrdinal*, Kokkos::HostSpace> OrdinalArrayHost;
78 typedef Kokkos::View<Scalar*, Kokkos::HostSpace> ScalarArrayHost;
80 typedef Kokkos::Serial functor_space;
81 typedef Kokkos::RangePolicy<functor_space, int> RangePol;
86 CountEntriesFunctor(
const TRowMatrix* A_, OrdinalArrayHost& rowptrs_) : A(A_), rowptrs(rowptrs_) {}
87 KOKKOS_INLINE_FUNCTION
void operator()(
const int row)
const
89 rowptrs(row) = A->getNumEntriesInLocalRow(row);
92 OrdinalArrayHost& rowptrs;
98 GetIndicesFunctor(
const TRowMatrix* A_,
const OrdinalArrayHost& rowptrs_, OrdinalArrayHost& colinds_) :
99 A(A_), rowptrs(rowptrs_), colinds(colinds_) {}
100 KOKKOS_INLINE_FUNCTION
void operator()(
const int row)
const
102 LocalOrdinal offset = rowptrs(row);
103 size_t entries = rowptrs(row + 1) - offset;
106 auto valsView = valsArray();
107 auto indicesView = indicesArray();
108 A->getLocalRowCopy(row, indicesView, valsView, entries);
109 std::sort(indicesView.getRawPtr(), indicesView().getRawPtr() + entries);
110 for(LocalOrdinal i = 0; i < entries; i++)
112 colinds(offset + i) = indicesView[i];
116 const OrdinalArrayHost& rowptrs;
117 OrdinalArrayHost& colinds;
123 GetValuesFunctor(
const TRowMatrix* A_,
const OrdinalArrayHost& rowptrs_, ScalarArrayHost& vals_) :
124 A(A_), rowptrs(rowptrs_), vals(vals_) {}
125 KOKKOS_INLINE_FUNCTION
void operator()(
const int row)
const
127 LocalOrdinal offset = rowptrs(row);
128 std::cout <<
"Rowptr[" << row <<
"] = " << offset <<
'\n';
129 size_t entries = rowptrs(row + 1) - offset;
132 auto valsView = valsArray();
133 auto indicesView = indicesArray();
134 A->getLocalRowCopy(row, indicesView, valsView, entries);
135 Tpetra::sort2(indicesView.getRawPtr(), indicesView.getRawPtr() + entries, valsView.getRawPtr());
136 for(LocalOrdinal i = 0; i < entries; i++)
138 vals(offset + i) = valsView[i];
142 const OrdinalArrayHost& rowptrs;
143 ScalarArrayHost& vals;
150 static void getValues(
const TRowMatrix* A, ScalarArray& vals, OrdinalArrayHost& rowptrs)
152 auto Acrs =
dynamic_cast<const TCrsMatrix*
>(A);
155 getValuesCrs(Acrs, vals);
159 LocalOrdinal nrows = A->getNodeNumRows();
160 ScalarArrayHost valsHost(
"Values (host)", rowptrs[nrows]);
161 vals = ScalarArray(
"Values", rowptrs[nrows]);
162 GetValuesFunctor funct(A, rowptrs, valsHost);
163 Kokkos::parallel_for(RangePol(0, nrows), funct);
164 Kokkos::deep_copy(vals, valsHost);
172 static void getStructure(
const TRowMatrix* A, OrdinalArrayHost& rowptrsHost, OrdinalArray& rowptrs, OrdinalArray& colinds)
174 auto Acrs =
dynamic_cast<const TCrsMatrix*
>(A);
177 getStructureCrs(Acrs, rowptrs, colinds);
183 LocalOrdinal nrows = A->getNodeNumRows();
184 rowptrsHost = OrdinalArrayHost(
"RowPtrs (host)", nrows + 1);
185 rowptrs = OrdinalArray(
"RowPtrs", nrows + 1);
186 CountEntriesFunctor countFunct(A, rowptrsHost);
187 Kokkos::parallel_for(RangePol(0, nrows), countFunct);
190 LocalOrdinal accum = 0;
191 for(LocalOrdinal i = 0; i < nrows; i++)
193 LocalOrdinal temp = rowptrsHost[i];
194 rowptrsHost[i] = accum;
197 rowptrsHost[nrows] = accum;
200 OrdinalArrayHost colindsHost = OrdinalArrayHost(
"ColInds (host)", rowptrsHost[nrows]);
201 colinds = OrdinalArray(
"ColInds", rowptrsHost[nrows]);
202 GetIndicesFunctor indicesFunct(A, rowptrsHost, colindsHost);
203 Kokkos::parallel_for(RangePol(0, nrows), indicesFunct);
205 Kokkos::deep_copy(rowptrs, rowptrsHost);
206 Kokkos::deep_copy(colinds, colindsHost);
212 static void getValuesCrs(
const TCrsMatrix* A, ScalarArray& vals)
214 vals = A->getLocalMatrix().values;
218 static void getStructureCrs(
const TCrsMatrix* A, OrdinalArray& rowptrs, OrdinalArray& colinds)
221 auto rowmap = A->getLocalMatrix().graph.row_map;
223 rowptrs = OrdinalArray(
"RowPtrs", A->getNodeNumRows() + 1);
224 for(
size_t i = 0; i < rowmap.extent(0); i++)
226 rowptrs[i] = rowmap[i];
228 colinds = A->getLocalMatrix().graph.entries;
Wraps a Tpetra::RowMatrix in a filter that reorders local rows and columns.
Definition: Ifpack2_ReorderFilter_decl.hpp:69
Functor for counting matrix entries per row in parallel.
Definition: Ifpack2_Details_CrsArrays.hpp:84
Functor for getting matrix values in parallel.
Definition: Ifpack2_Details_CrsArrays.hpp:121
Functor for getting column indices in parallel.
Definition: Ifpack2_Details_CrsArrays.hpp:96
Access only local rows and columns of a sparse matrix.
Definition: Ifpack2_LocalFilter_decl.hpp:160