48#ifndef __IFPACK2_CRSARRAYS_DECL_HPP__
49#define __IFPACK2_CRSARRAYS_DECL_HPP__
51#include <Tpetra_RowMatrix.hpp>
52#include <Tpetra_CrsMatrix.hpp>
53#include <KokkosCompat_DefaultNode.hpp>
54#include <KokkosSparse_CrsMatrix.hpp>
55#include <Ifpack2_LocalFilter.hpp>
56#include <Ifpack2_ReorderFilter.hpp>
65template<
typename Scalar,
typename ImplScalar,
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<ImplScalar, LocalOrdinal, execution_space> KCrsMatrix;
75 typedef Kokkos::View<LocalOrdinal*, execution_space> OrdinalArray;
76 typedef Kokkos::View<ImplScalar*, execution_space> ScalarArray;
77 typedef typename OrdinalArray::HostMirror OrdinalArrayHost;
79 typedef Kokkos::Serial functor_space;
80 typedef Kokkos::RangePolicy<functor_space, int> RangePol;
86 static void getValues(
const TRowMatrix* A, ScalarArray& vals, OrdinalArrayHost& rowptrs)
88 auto Acrs =
dynamic_cast<const TCrsMatrix*
>(A);
91 getValuesCrs(Acrs, vals);
94 using range_type = Kokkos::pair<int, int>;
95 using local_inds_host_view_type =
typename TRowMatrix::nonconst_local_inds_host_view_type;
96 using values_host_view_type =
typename TRowMatrix::nonconst_values_host_view_type;
97 using scalar_type =
typename values_host_view_type::value_type;
99 LocalOrdinal nrows = A->getLocalNumRows();
100 size_t nnz = A->getLocalNumEntries();
101 size_t maxNnz = A->getLocalMaxNumRowEntries();
103 vals = ScalarArray(
"Values", nnz);
104 auto valsHost = Kokkos::create_mirror(vals);
105 local_inds_host_view_type lclColInds (
"lclColinds", maxNnz);
108 for(LocalOrdinal i = 0; i < nrows; i++) {
109 size_t NumEntries = A->getNumEntriesInLocalRow(i);
110 auto constLclValues = Kokkos::subview (valsHost, range_type (nnz, nnz+NumEntries));
111 values_host_view_type lclValues (
const_cast<scalar_type*
>(constLclValues.data()), NumEntries);
113 A->getLocalRowCopy (i, lclColInds, lclValues, NumEntries);
116 Kokkos::deep_copy(vals, valsHost);
124 static void getStructure(
const TRowMatrix* A, OrdinalArrayHost& rowptrsHost, OrdinalArray& rowptrs, OrdinalArray& colinds)
126 auto Acrs =
dynamic_cast<const TCrsMatrix*
>(A);
129 getStructureCrs(Acrs, rowptrsHost, rowptrs, colinds);
135 LocalOrdinal nrows = A->getLocalNumRows();
136 rowptrsHost = OrdinalArrayHost(
"RowPtrs (host)", nrows + 1);
138 using range_type = Kokkos::pair<int, int>;
139 using values_host_view_type =
typename TRowMatrix::nonconst_values_host_view_type;
140 using local_inds_host_view_type =
typename TRowMatrix::nonconst_local_inds_host_view_type;
141 using local_ind_type =
typename local_inds_host_view_type::value_type;
142 size_t nnz = A->getLocalNumEntries();
143 size_t maxNnz = A->getLocalMaxNumRowEntries();
145 colinds = OrdinalArray(
"ColInds", nnz);
146 auto colindsHost = Kokkos::create_mirror(colinds);
147 values_host_view_type lclValues (
"lclValues", maxNnz);
150 rowptrsHost[0] = nnz;
151 for(LocalOrdinal i = 0; i < nrows; i++) {
152 size_t NumEntries = A->getNumEntriesInLocalRow(i);
153 auto constLclValues = Kokkos::subview (colindsHost, range_type (nnz, nnz+NumEntries));
154 local_inds_host_view_type lclColInds (
const_cast<local_ind_type*
>(constLclValues.data()), NumEntries);
155 A->getLocalRowCopy (i, lclColInds, lclValues, NumEntries);
158 rowptrsHost[i+1] = nnz;
161 rowptrs = OrdinalArray(
"RowPtrs", nrows + 1);
162 Kokkos::deep_copy(rowptrs, rowptrsHost);
163 Kokkos::deep_copy(colinds, colindsHost);
169 static void getValuesCrs(
const TCrsMatrix* A, ScalarArray& values_)
171 auto localA = A->getLocalMatrixDevice();
172 auto values = localA.values;
173 auto nnz = values.extent(0);
174 values_ = ScalarArray(
"Values", nnz );
175 Kokkos::deep_copy(values_, values);
179 static void getStructureCrs(
const TCrsMatrix* A, OrdinalArrayHost& rowptrsHost_, OrdinalArray& rowptrs_, OrdinalArray& colinds_)
182 auto localA = A->getLocalMatrixDevice();
183 auto rowptrs = localA.graph.row_map;
184 auto colinds = localA.graph.entries;
185 auto numRows = A->getLocalNumRows();
186 auto nnz = colinds.extent(0);
188 rowptrs_ = OrdinalArray(
"RowPtrs", numRows + 1);
189 colinds_ = OrdinalArray(
"ColInds", nnz );
190 Kokkos::deep_copy(rowptrs_, rowptrs);
191 Kokkos::deep_copy(colinds_, colinds);
193 rowptrsHost_ = Kokkos::create_mirror(rowptrs_);
194 Kokkos::deep_copy(rowptrsHost_, rowptrs_);
Access only local rows and columns of a sparse matrix.
Definition Ifpack2_LocalFilter_decl.hpp:163
Wraps a Tpetra::RowMatrix in a filter that reorders local rows and columns.
Definition Ifpack2_ReorderFilter_decl.hpp:70
Ifpack2 implementation details.
Preconditioners and smoothers for Tpetra sparse matrices.
Definition Ifpack2_AdditiveSchwarz_decl.hpp:74