Tpetra parallel linear algebra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Tpetra_Details_determineLocalTriangularStructure.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_DETERMINELOCALTRIANGULARSTRUCTURE_HPP
11 #define TPETRA_DETAILS_DETERMINELOCALTRIANGULARSTRUCTURE_HPP
12 
19 
20 #include "Kokkos_Core.hpp"
22 
23 namespace Tpetra {
24 namespace Details {
25 
29 template <class LO>
39 };
40 
41 namespace Impl {
55 template <class LocalGraphType, class LocalMapType>
57  public:
58  // Result can't be more than the number of local rows, so
59  // local_ordinal_type is appropriate.
60  using result_type =
62 
72  DetermineLocalTriangularStructure(const LocalGraphType& G,
73  const LocalMapType& rowMap,
74  const LocalMapType& colMap,
75  const bool ignoreMapsForTriangularStructure)
76  : G_(G)
77  , rowMap_(rowMap)
78  , colMap_(colMap)
79  , ignoreMapsForTriangularStructure_(ignoreMapsForTriangularStructure) {}
80 
82  KOKKOS_INLINE_FUNCTION void init(result_type& dst) const {
83  dst.diagCount = 0;
84  dst.maxNumRowEnt = 0;
85  dst.couldBeLowerTriangular = true; // well, we don't know yet, do we?
86  dst.couldBeUpperTriangular = true; // ditto
87  }
88 
89  KOKKOS_INLINE_FUNCTION void
90  join(result_type& dst,
91  const result_type& src) const {
92  dst.diagCount += src.diagCount;
93  dst.maxNumRowEnt = (src.maxNumRowEnt > dst.maxNumRowEnt) ? src.maxNumRowEnt : dst.maxNumRowEnt;
94  dst.couldBeLowerTriangular &= src.couldBeLowerTriangular;
95  dst.couldBeUpperTriangular &= src.couldBeUpperTriangular;
96  }
97 
99  KOKKOS_INLINE_FUNCTION void
100  operator()(const typename LocalMapType::local_ordinal_type lclRow,
101  result_type& result) const {
102  using LO = typename LocalMapType::local_ordinal_type;
103  using GO = typename LocalMapType::global_ordinal_type;
104  using LOT = typename ::Tpetra::Details::OrdinalTraits<LO>;
105 
106  auto G_row = G_.rowConst(lclRow);
107  const LO numEnt = G_row.length;
108  if (numEnt != 0) {
109  result.maxNumRowEnt = (numEnt > result.maxNumRowEnt) ? numEnt : result.maxNumRowEnt;
110  // Use global row and column indices to find the diagonal
111  // entry. Caller promises that local row index is in the row
112  // Map on the calling process.
113  const GO gblDiagCol = rowMap_.getGlobalElement(lclRow);
114  const LO lclDiagCol = colMap_.getLocalElement(gblDiagCol);
115  // If it's not in the column Map, then there's no diagonal entry.
116  if (lclDiagCol != LOT::invalid()) {
117  // TODO (mfh 25 Apr 2018) Use findRelOffset to optimize for
118  // the sorted case, but note that it requires operator[].
119  bool foundDiag = false; // don't count duplicates
120 
121  if (ignoreMapsForTriangularStructure_) {
122  for (LO k = 0; k < numEnt && !foundDiag; ++k) {
123  const LO lclCol = G_row(k);
124  if (lclCol == lclDiagCol) {
125  foundDiag = true;
126  }
127  }
128  // mfh 30 Apr 2018: See GitHub Issue #2658. Per
129  // current Tpetra::CrsGraph::computeLocalConstants
130  // behavior, assume that local column indices are
131  // sorted in each row.
132  if (numEnt > LO(0)) {
133  const LO smallestLclCol = G_row(0);
134  const LO largestLclCol = G_row(numEnt - 1); // could be same
135 
136  if (smallestLclCol < lclRow) {
137  result.couldBeUpperTriangular = false;
138  }
139  if (lclRow < largestLclCol) {
140  result.couldBeLowerTriangular = false;
141  }
142  }
143  } else {
144  for (LO k = 0; k < numEnt &&
145  ((!foundDiag) ||
146  result.couldBeLowerTriangular ||
147  result.couldBeUpperTriangular);
148  ++k) {
149  const LO lclCol = G_row(k);
150  if (lclCol == lclDiagCol) {
151  foundDiag = true;
152  } else {
153  const GO gblCol = colMap_.getGlobalElement(lclCol);
154  if (gblCol < gblDiagCol) {
155  result.couldBeUpperTriangular = false;
156  }
157  if (gblDiagCol < gblCol) {
158  result.couldBeLowerTriangular = false;
159  }
160  }
161  } // for each entry in lclRow
162  } // if-else ignoreMapsForTriangularStructure
163 
164  if (foundDiag) {
165  ++(result.diagCount);
166  }
167  }
168  }
169  }
170 
171  private:
172  LocalGraphType G_;
173  LocalMapType rowMap_;
174  LocalMapType colMap_;
175  bool ignoreMapsForTriangularStructure_;
176 };
177 
178 } // namespace Impl
179 
197 template <class LocalGraphType, class LocalMapType>
198 LocalTriangularStructureResult<typename LocalMapType::local_ordinal_type>
199 determineLocalTriangularStructure(const LocalGraphType& G,
200  const LocalMapType& rowMap,
201  const LocalMapType& colMap,
202  const bool ignoreMapsForTriangularStructure) {
203  using LO = typename LocalMapType::local_ordinal_type;
204  using execution_space = typename LocalGraphType::device_type::execution_space;
205  using range_type = Kokkos::RangePolicy<execution_space, LO>;
206  using functor_type =
208 
209  LocalTriangularStructureResult<LO> result{0, 0, true, true};
210  Kokkos::parallel_reduce("Tpetra::Details::determineLocalTriangularStructure",
211  range_type(0, G.numRows()),
212  functor_type(G, rowMap, colMap,
213  ignoreMapsForTriangularStructure),
214  result);
215  return result;
216 }
217 
218 } // namespace Details
219 } // namespace Tpetra
220 
221 #endif // TPETRA_DETAILS_DETERMINELOCALTRIANGULARSTRUCTURE_HPP
Import KokkosSparse::OrdinalTraits, a traits class for &quot;invalid&quot; (flag) values of integer types...
bool couldBeLowerTriangular
Whether the graph is locally structurally lower triangular.
KOKKOS_INLINE_FUNCTION void init(result_type &dst) const
Set the initial value of the reduction result.
KOKKOS_INLINE_FUNCTION void operator()(const typename LocalMapType::local_ordinal_type lclRow, result_type &result) const
Reduction operator: result is (diagonal count, error count).
DetermineLocalTriangularStructure(const LocalGraphType &G, const LocalMapType &rowMap, const LocalMapType &colMap, const bool ignoreMapsForTriangularStructure)
Constructor.
LocalTriangularStructureResult< typename LocalMapType::local_ordinal_type > determineLocalTriangularStructure(const LocalGraphType &G, const LocalMapType &rowMap, const LocalMapType &colMap, const bool ignoreMapsForTriangularStructure)
Count the local number of diagonal entries in a local sparse graph, and determine whether the local p...
bool couldBeUpperTriangular
Whether the graph is locally structurally upper triangular.
Implementation of Tpetra::Details::determineLocalTriangularStructure (which see below).