Tpetra parallel linear algebra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Tpetra_Details_EquilibrationInfo.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Tpetra: Templated Linear Algebra Services Package
5 // Copyright (2008) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ************************************************************************
40 // @HEADER
41 
42 #ifndef TPETRA_DETAILS_EQUILIBRATIONINFO_HPP
43 #define TPETRA_DETAILS_EQUILIBRATIONINFO_HPP
44 
47 
48 #include "TpetraCore_config.h"
49 #include "Kokkos_ArithTraits.hpp"
50 #include "Kokkos_Core.hpp"
51 
52 namespace Tpetra {
53 namespace Details {
54 
78 template<class ScalarType, class DeviceType>
80  using val_type = typename Kokkos::ArithTraits<ScalarType>::val_type;
81  using mag_type = typename Kokkos::ArithTraits<val_type>::mag_type;
82  using device_type = typename DeviceType::device_type;
83  using host_device_type = typename Kokkos::View<mag_type*, device_type>::HostMirror::device_type;
85 
87  foundInf (false),
88  foundNan (false),
89  foundZeroDiag (false),
90  foundZeroRowNorm (false)
91  {}
92 
93  EquilibrationInfo (const std::size_t lclNumRows,
94  const std::size_t lclNumCols,
95  const bool assumeSymmetric_) :
96  rowNorms (Kokkos::View<mag_type*, device_type> ("rowNorms", lclNumRows)),
97  rowDiagonalEntries (Kokkos::View<val_type*, device_type> ("rowDiagonalEntries", lclNumRows)),
98  colNorms (Kokkos::View<mag_type*, device_type> ("colNorms", lclNumCols)),
99  colDiagonalEntries (Kokkos::View<val_type*, device_type> ("colDiagonalEntries",
100  assumeSymmetric_ ?
101  std::size_t (0) :
102  lclNumCols)),
103  rowScaledColNorms (Kokkos::View<mag_type*, device_type> ("rowScaledColNorms",
104  assumeSymmetric_ ?
105  std::size_t (0) :
106  lclNumCols)),
107  assumeSymmetric (assumeSymmetric_),
108  foundInf (false),
109  foundNan (false),
110  foundZeroDiag (false),
111  foundZeroRowNorm (false)
112  {}
113 
114  EquilibrationInfo (const Kokkos::View<mag_type*, device_type>& rowNorms_,
115  const Kokkos::View<val_type*, device_type>& rowDiagonalEntries_,
116  const Kokkos::View<mag_type*, device_type>& colNorms_,
117  const Kokkos::View<val_type*, device_type>& colDiagonalEntries_,
118  const Kokkos::View<mag_type*, device_type>& rowScaledColNorms_,
119  const bool assumeSymmetric_,
120  const bool foundInf_,
121  const bool foundNan_,
122  const bool foundZeroDiag_,
123  const bool foundZeroRowNorm_) :
124  rowNorms (rowNorms_),
125  rowDiagonalEntries (rowDiagonalEntries_),
126  colNorms (colNorms_),
127  colDiagonalEntries (colDiagonalEntries_),
128  rowScaledColNorms (rowScaledColNorms_),
129  assumeSymmetric (assumeSymmetric_),
130  foundInf (foundInf_),
131  foundNan (foundNan_),
132  foundZeroDiag (foundZeroDiag_),
133  foundZeroRowNorm (foundZeroRowNorm_)
134  {}
135 
137  template<class SrcDeviceType>
138  void
140  {
141  using execution_space = typename device_type::execution_space;
142  // DEEP_COPY REVIEW - DEVICE-TO-DEVICE
143  Kokkos::deep_copy (execution_space(), rowNorms, src.rowNorms);
144  // DEEP_COPY REVIEW - DEVICE-TO-DEVICE
145  Kokkos::deep_copy (execution_space(), rowDiagonalEntries, src.rowDiagonalEntries);
146  // DEEP_COPY REVIEW - DEVICE-TO-DEVICE
147  Kokkos::deep_copy (execution_space(), colNorms, src.colNorms);
148  if (src.colDiagonalEntries.extent (0) == 0) {
150  Kokkos::View<val_type*, device_type> ("colDiagonalEntries", 0);
151  }
152  else {
153  // DEEP_COPY REVIEW - DEVICE-TO-DEVICE
154  Kokkos::deep_copy (execution_space(), colDiagonalEntries, src.colDiagonalEntries);
155  }
156  if (src.rowScaledColNorms.extent (0) == 0) {
158  Kokkos::View<mag_type*, device_type> ("rowScaledColNorms", 0);
159  }
160  else {
161  // DEEP_COPY REVIEW - DEVICE-TO-DEVICE
162  Kokkos::deep_copy (execution_space(), rowScaledColNorms, src.rowScaledColNorms);
163  }
164 
166  foundInf = src.foundInf;
167  foundNan = src.foundNan;
170  }
171 
173  createMirrorView ()
174  {
175  auto rowNorms_h = Kokkos::create_mirror_view (rowNorms);
176  auto rowDiagonalEntries_h = Kokkos::create_mirror_view (rowDiagonalEntries);
177  auto colNorms_h = Kokkos::create_mirror_view (colNorms);
178  auto colDiagonalEntries_h = Kokkos::create_mirror_view (colDiagonalEntries);
179  auto rowScaledColNorms_h = Kokkos::create_mirror_view (rowScaledColNorms);
180 
181  return HostMirror {rowNorms_h, rowDiagonalEntries_h, colNorms_h,
182  colDiagonalEntries_h, rowScaledColNorms_h, assumeSymmetric,
184  }
185 
186  // We call a row a "diagonally dominant row" if the absolute value
187  // of the diagonal entry is >= the sum of the absolute values of the
188  // off-diagonal entries. The row norm is the sum of those two
189  // things, so this means diagAbsVal >= rowNorm - diagAbsVal. Ditto
190  // for a column.
191 
193  Kokkos::View<mag_type*, device_type> rowNorms;
194 
196  Kokkos::View<val_type*, device_type> rowDiagonalEntries;
197 
203  Kokkos::View<mag_type*, device_type> colNorms;
204 
208  Kokkos::View<val_type*, device_type> colDiagonalEntries;
209 
220  Kokkos::View<mag_type*, device_type> rowScaledColNorms;
221 
227 
229  bool foundInf;
230 
232  bool foundNan;
233 
236 
239 };
240 
241 } // namespace Details
242 } // namespace Tpetra
243 
244 #endif // TPETRA_DETAILS_EQUILIBRATIONINFO_HPP
Kokkos::View< mag_type *, device_type > colNorms
One-norms of the matrix&#39;s columns, distributed via the column Map.
Kokkos::View< mag_type *, device_type > rowNorms
One-norms of the matrix&#39;s rows, distributed via the row Map.
bool foundZeroDiag
Found a zero diagonal entry somewhere in the matrix.
bool assumeSymmetric
Whether to assume that the matrix is (globally) symmetric.
Kokkos::View< val_type *, device_type > colDiagonalEntries
Diagonal entries of the matrix, distributed via the column Map.
Struct storing results of Tpetra::computeRowAndColumnOneNorms.
void deep_copy(MultiVector< DS, DL, DG, DN > &dst, const MultiVector< SS, SL, SG, SN > &src)
Copy the contents of the MultiVector src into dst.
Kokkos::View< val_type *, device_type > rowDiagonalEntries
Diagonal entries of the matrix, distributed via the row Map.
Kokkos::View< mag_type *, device_type > rowScaledColNorms
One-norms of the matrix&#39;s columns, after the matrix&#39;s rows have been scaled by rowNorms.
void assign(const EquilibrationInfo< ScalarType, SrcDeviceType > &src)
Deep-copy src into *this.
bool foundNan
Found a NaN somewhere in the matrix.
bool foundInf
Found an Inf somewhere in the matrix.
bool foundZeroRowNorm
At least one row of the matrix has a zero norm.