42#ifndef TPETRA_DETAILS_GETNUMDIAGS_HPP
43#define TPETRA_DETAILS_GETNUMDIAGS_HPP
52#include "Tpetra_CrsGraph.hpp"
53#include "Teuchos_CommHelpers.hpp"
65 template<
class LocalGraphType,
class LocalMapType>
71 G_ (
G), rowMap_ (rowMap), colMap_ (colMap)
76 using result_type =
typename LocalMapType::local_ordinal_type;
81 result_type& diagCount)
const
83 using LO =
typename LocalMapType::local_ordinal_type;
84 using LOT = typename ::Tpetra::Details::OrdinalTraits<LO>;
92 const LO
lclDiagCol = colMap_.getLocalElement (rowMap_.getGlobalElement (
lclRow));
113 template<
class LO,
class GO,
class NT>
114 typename ::Tpetra::CrsGraph<LO, GO, NT>::local_ordinal_type
115 countLocalNumDiagsInFillCompleteGraph (const ::Tpetra::CrsGraph<LO, GO, NT>&
G)
118 using local_map_type =
typename crs_graph_type::map_type::local_map_type;
119 using local_graph_device_type =
typename crs_graph_type::local_graph_device_type;
121 using execution_space =
typename crs_graph_type::device_type::execution_space;
122 using policy_type = Kokkos::RangePolicy<execution_space, LO>;
124 const auto rowMap =
G.getRowMap ();
125 const auto colMap =
G.getColMap ();
126 if (rowMap.get () ==
nullptr || colMap.get () ==
nullptr) {
131 functor_type f (G.getLocalGraphDevice (), rowMap->getLocalMap (), colMap->getLocalMap ());
132 Kokkos::parallel_reduce (policy_type (0, G.getLocalNumRows ()), f, lclNumDiags);
146 template<
class MapType>
147 typename MapType::local_ordinal_type
152 return colMap.getLocalElement (rowMap.getGlobalElement (
lclRow));
156 template<
class LO,
class GO,
class NT>
157 typename ::Tpetra::RowGraph<LO, GO, NT>::local_ordinal_type
160 using LOT = typename ::Tpetra::Details::OrdinalTraits<LO>;
162 const auto rowMap =
G.getRowMap ();
163 const auto colMap =
G.getColMap ();
164 if (rowMap.get () ==
nullptr || colMap.get () ==
nullptr) {
169 (!
G.supportsRowViews (), std::logic_error,
"Not implemented!");
171 typename ::Tpetra::RowGraph<LO, GO, NT>::local_inds_host_view_type
173 const LO
lclNumRows =
static_cast<LO
> (
G.getLocalNumRows ());
180 const LO
lclDiagCol = colMap->getLocalElement (rowMap->getGlobalElement (
lclRow));
200 template<
class LO,
class GO,
class NT>
201 typename ::Tpetra::RowGraph<LO, GO, NT>::local_ordinal_type
204 using LOT = typename ::Tpetra::Details::OrdinalTraits<LO>;
206 const auto rowMap =
G.getRowMap ();
207 const auto colMap =
G.getColMap ();
208 if (rowMap.get () ==
nullptr || colMap.get () ==
nullptr) {
212 using inds_type = typename ::Tpetra::RowGraph<LO,GO,NT>::nonconst_local_inds_host_view_type;
214 const LO
lclNumRows =
static_cast<LO
> (
G.getLocalNumRows ());
226 colMap->getLocalElement (rowMap->getGlobalElement (
lclRow));
246 template<
class LO,
class GO,
class NT>
247 typename ::Tpetra::RowGraph<LO, GO, NT>::local_ordinal_type
250 const auto rowMap =
G.getRowMap ();
251 if (rowMap.get () ==
nullptr) {
255 typename ::Tpetra::RowGraph<LO,GO,NT>::global_inds_host_view_type
257 const LO
lclNumRows =
static_cast<LO
> (
G.getLocalNumRows ());
281 template<
class LO,
class GO,
class NT>
282 typename ::Tpetra::RowGraph<LO, GO, NT>::local_ordinal_type
285 using gids_type = typename ::Tpetra::RowGraph<LO,GO,NT>::nonconst_global_inds_host_view_type ;
286 const auto rowMap =
G.getRowMap ();
287 if (rowMap.get () ==
nullptr) {
292 const LO
lclNumRows =
static_cast<LO
> (
G.getLocalNumRows ());
328 template<
class MatrixType>
330 static typename MatrixType::local_ordinal_type
333 using LO =
typename MatrixType::local_ordinal_type;
334 using GO =
typename MatrixType::global_ordinal_type;
335 using NT =
typename MatrixType::node_type;
338 auto G =
A.getGraph ();
339 if (
G.get () ==
nullptr) {
349 template<
class LO,
class GO,
class NT>
352 getLocalNumDiags (const ::Tpetra::RowGraph<LO, GO, NT>&
G)
356 const crs_graph_type*
G_crs =
dynamic_cast<const crs_graph_type*
> (&
G);
357 if (
G_crs !=
nullptr &&
G_crs->isFillComplete ()) {
358 return countLocalNumDiagsInFillCompleteGraph (*
G_crs);
361 if (
G.isLocallyIndexed ()) {
362 if (
G.supportsRowViews ()) {
369 else if (
G.isGloballyIndexed ()) {
370 if (
G.supportsRowViews ()) {
385 template<
class LO,
class GO,
class NT>
388 getLocalNumDiags (const ::Tpetra::CrsGraph<LO, GO, NT>&
G)
398template<
class CrsGraphType>
399typename CrsGraphType::local_ordinal_type
407template<
class CrsGraphType>
408typename CrsGraphType::global_ordinal_type
411 using GO =
typename CrsGraphType::global_ordinal_type;
413 const auto map =
G.getRowMap ();
414 if (
map.get () ==
nullptr) {
418 const auto comm =
map->getComm ();
419 if (comm.get () ==
nullptr) {
425 using Teuchos::REDUCE_SUM;
426 using Teuchos::reduceAll;
427 using Teuchos::outArg;
Import KokkosSparse::OrdinalTraits, a traits class for "invalid" (flag) values of integer types,...
typename::Tpetra::RowGraph< LO, GO, NT >::local_ordinal_type countLocalNumDiagsInNonFillCompleteGloballyIndexedGraphWithoutRowViews(const ::Tpetra::RowGraph< LO, GO, NT > &G)
Return local number of diagonal entries.
typename::Tpetra::RowGraph< LO, GO, NT >::local_ordinal_type countLocalNumDiagsInNonFillCompleteGloballyIndexedGraphWithRowViews(const ::Tpetra::RowGraph< LO, GO, NT > &G)
Return local number of diagonal entries.
MapType::local_ordinal_type getLocalDiagonalColumnIndex(const typename MapType::local_ordinal_type lclRow, const MapType &rowMap, const MapType &colMap)
Local columm index of diagonal entry.
typename::Tpetra::RowGraph< LO, GO, NT >::local_ordinal_type countLocalNumDiagsInNonFillCompleteLocallyIndexedGraphWithoutRowViews(const ::Tpetra::RowGraph< LO, GO, NT > &G)
Return local number of diagonal entries.
typename::Tpetra::RowGraph< LO, GO, NT >::local_ordinal_type countLocalNumDiagsInNonFillCompleteLocallyIndexedGraphWithRowViews(const ::Tpetra::RowGraph< LO, GO, NT > &G)
Return local number of diagonal entries.
A distributed graph accessed by rows (adjacency lists) and stored sparsely.
Struct that holds views of the contents of a CrsMatrix.
Kokkos::parallel_reduce functor for counting the local number of diagonal entries in a sparse graph.
KOKKOS_INLINE_FUNCTION void operator()(const typename LocalMapType::local_ordinal_type lclRow, result_type &diagCount) const
Reduction function: result is (diagonal count, error count).
An abstract interface for graphs accessed by rows.
Implementation details of Tpetra.
CrsGraphType::global_ordinal_type getGlobalNumDiags(const CrsGraphType &G)
Number of populated diagonal entries in the given sparse graph, over all processes in the graph's (MP...
CrsGraphType::local_ordinal_type getLocalNumDiags(const CrsGraphType &G)
Number of populated diagonal entries in the given sparse graph, on the calling (MPI) process.
Namespace Tpetra contains the class and methods constituting the Tpetra library.
Implementation of Tpetra::Details::getLocalNumDiags (see below).