Tpetra parallel linear algebra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Tpetra_Details_getDiagCopyWithoutOffsets_decl.hpp
Go to the documentation of this file.
1 /*
2 // @HEADER
3 // ***********************************************************************
4 //
5 // Tpetra: Templated Linear Algebra Services Package
6 // Copyright (2008) Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39 //
40 // ************************************************************************
41 // @HEADER
42 */
43 
44 #ifndef TPETRA_DETAILS_GETDIAGCOPYWITHOUTOFFSETS_DECL_HPP
45 #define TPETRA_DETAILS_GETDIAGCOPYWITHOUTOFFSETS_DECL_HPP
46 
55 
56 #include "TpetraCore_config.h"
57 #include "Kokkos_Core.hpp"
58 #include "Kokkos_ArithTraits.hpp"
60 #include "Tpetra_RowMatrix_decl.hpp"
61 #include "Tpetra_Vector_decl.hpp"
62 #include <type_traits>
63 
64 namespace Tpetra {
65 namespace Details {
66 
79 template<class DiagType,
80  class LocalMapType,
81  class CrsMatrixType>
83  typedef typename LocalMapType::local_ordinal_type LO; // local ordinal type
84  typedef typename LocalMapType::global_ordinal_type GO; // global ordinal type
85  typedef typename CrsMatrixType::device_type device_type;
86  typedef typename CrsMatrixType::value_type scalar_type;
87  typedef typename CrsMatrixType::size_type offset_type;
88 
90  typedef LO value_type;
91 
99  CrsMatrixGetDiagCopyFunctor (const DiagType& D,
100  const LocalMapType& rowMap,
101  const LocalMapType& colMap,
102  const CrsMatrixType& A) :
103  D_ (D), rowMap_ (rowMap), colMap_ (colMap), A_ (A)
104  {}
105 
109  KOKKOS_FUNCTION void
110  operator () (const LO& lclRowInd, value_type& errCount) const
111  {
112  const LO INV = Tpetra::Details::OrdinalTraits<LO>::invalid ();
113  const scalar_type ZERO =
114  Kokkos::ArithTraits<scalar_type>::zero ();
115 
116  // If the row lacks a stored diagonal entry, then its value is zero.
117  D_(lclRowInd) = ZERO;
118  const GO gblInd = rowMap_.getGlobalElement (lclRowInd);
119  const LO lclColInd = colMap_.getLocalElement (gblInd);
120 
121  if (lclColInd != INV) {
122  auto curRow = A_.rowConst (lclRowInd);
123 
124  // FIXME (mfh 12 May 2016) Use binary search when the row is
125  // long enough. findRelOffset currently lives in KokkosKernels
126  // (in tpetra/kernels/src/Kokkos_Sparse_findRelOffset.hpp).
127  LO offset = 0;
128  const LO numEnt = curRow.length;
129  for ( ; offset < numEnt; ++offset) {
130  if (curRow.colidx(offset) == lclColInd) {
131  break;
132  }
133  }
134 
135  if (offset == numEnt) {
136  ++errCount;
137  }
138  else {
139  D_(lclRowInd) = curRow.value(offset);
140  }
141  }
142  }
143 
144 private:
146  DiagType D_;
148  LocalMapType rowMap_;
150  LocalMapType colMap_;
152  CrsMatrixType A_;
153 };
154 
155 
178 template<class DiagType,
179  class LocalMapType,
180  class CrsMatrixType>
181 static typename LocalMapType::local_ordinal_type
182 getDiagCopyWithoutOffsets (const DiagType& D,
183  const LocalMapType& rowMap,
184  const LocalMapType& colMap,
185  const CrsMatrixType& A)
186 {
187  static_assert (Kokkos::is_view<DiagType>::value,
188  "DiagType must be a Kokkos::View.");
189  static_assert (static_cast<int> (DiagType::rank) == 1,
190  "DiagType must be a 1-D Kokkos::View.");
191  static_assert (std::is_same<typename DiagType::value_type, typename DiagType::non_const_value_type>::value,
192  "DiagType must be a nonconst Kokkos::View.");
193  typedef typename LocalMapType::local_ordinal_type LO;
194  typedef typename CrsMatrixType::device_type device_type;
195  typedef typename device_type::execution_space execution_space;
196  typedef Kokkos::RangePolicy<execution_space, LO> policy_type;
197 
198  typedef Kokkos::View<typename DiagType::non_const_value_type*,
199  typename DiagType::array_layout,
200  typename DiagType::device_type,
201  Kokkos::MemoryUnmanaged> diag_type;
202  diag_type D_out = D;
204  functor (D_out, rowMap, colMap, A);
205  const LO numRows = static_cast<LO> (D.extent (0));
206  LO errCount = 0;
207  Kokkos::parallel_reduce (policy_type (0, numRows), functor, errCount);
208  return errCount;
209 }
210 
247 template<class SC, class LO, class GO, class NT>
248 LO
250  const ::Tpetra::RowMatrix<SC, LO, GO, NT>& A,
251  const bool debug =
252 #ifdef HAVE_TPETRA_DEBUG
253  true);
254 #else // ! HAVE_TPETRA_DEBUG
255  false);
256 #endif // HAVE_TPETRA_DEBUG
257 
258 } // namespace Details
259 } // namespace Tpetra
260 
261 #endif // TPETRA_DETAILS_GETDIAGCOPYWITHOUTOFFSETS_DECL_HPP
Import KokkosSparse::OrdinalTraits, a traits class for &quot;invalid&quot; (flag) values of integer types...
Functor that implements much of the one-argument overload of Tpetra::CrsMatrix::getLocalDiagCopy, for the case where the matrix is fill complete.
LO value_type
The result of the reduction; number of errors.
Declaration of the Tpetra::Vector class.
KOKKOS_FUNCTION void operator()(const LO &lclRowInd, value_type &errCount) const
Operator for Kokkos::parallel_for.
static LocalMapType::local_ordinal_type getDiagCopyWithoutOffsets(const DiagType &D, const LocalMapType &rowMap, const LocalMapType &colMap, const CrsMatrixType &A)
Given a locally indexed, local sparse matrix, and corresponding local row and column Maps...
LO getLocalDiagCopyWithoutOffsetsNotFillComplete(::Tpetra::Vector< SC, LO, GO, NT > &diag, const ::Tpetra::RowMatrix< SC, LO, GO, NT > &A, const bool debug=false)
Given a locally indexed, global sparse matrix, extract the matrix&#39;s diagonal entries into a Tpetra::V...
Replace old values with zero.
A distributed dense vector.
CrsMatrixGetDiagCopyFunctor(const DiagType &D, const LocalMapType &rowMap, const LocalMapType &colMap, const CrsMatrixType &A)
Constructor.