10 #ifndef TPETRA_LEFTANDORRIGHTSCALECRSMATRIX_DEF_HPP
11 #define TPETRA_LEFTANDORRIGHTSCALECRSMATRIX_DEF_HPP
20 #include "Tpetra_CrsMatrix.hpp"
21 #include "Tpetra_Vector.hpp"
25 #include "Teuchos_TestForException.hpp"
29 template <
class SC,
class LO,
class GO,
class NT>
32 #
if KOKKOS_VERSION >= 40799
33 const typename KokkosKernels::ArithTraits<SC>::mag_type*,
35 const typename Kokkos::ArithTraits<SC>::mag_type*,
37 typename NT::device_type>& rowScalingFactors,
39 #
if KOKKOS_VERSION >= 40799
40 const typename KokkosKernels::ArithTraits<SC>::mag_type*,
42 const typename Kokkos::ArithTraits<SC>::mag_type*,
44 typename NT::device_type>& colScalingFactors,
46 const bool rightScale,
47 const bool assumeSymmetric,
49 if (!leftScale && !rightScale) {
54 if (!A_fillComplete_on_input) {
62 static_cast<LO
>(A.
getRowMap()->getLocalNumElements());
63 TEUCHOS_TEST_FOR_EXCEPTION(A_lcl.numRows() != lclNumRows, std::invalid_argument,
64 "leftAndOrRightScaleCrsMatrix: Local matrix is not valid. "
65 "This means that A was not created with a local matrix, "
66 "and that fillComplete has never yet been called on A before. "
67 "Please call fillComplete on A at least once first "
68 "before calling this method.");
73 const bool divide = scaling == SCALING_DIVIDE;
87 if (A_fillComplete_on_input) {
88 Teuchos::RCP<Teuchos::ParameterList> params = Teuchos::parameterList();
89 params->set(
"No Nonlocal Changes",
true);
94 template <
class SC,
class LO,
class GO,
class NT>
97 #
if KOKKOS_VERSION >= 40799
98 typename KokkosKernels::ArithTraits<SC>::mag_type,
100 typename Kokkos::ArithTraits<SC>::mag_type,
102 LO, GO, NT>& rowScalingFactors,
104 #
if KOKKOS_VERSION >= 40799
105 typename KokkosKernels::ArithTraits<SC>::mag_type,
107 typename Kokkos::ArithTraits<SC>::mag_type,
109 LO, GO, NT>& colScalingFactors,
110 const bool leftScale,
111 const bool rightScale,
112 const bool assumeSymmetric,
114 using device_type =
typename NT::device_type;
115 using dev_memory_space =
typename device_type::memory_space;
116 #if KOKKOS_VERSION >= 40799
117 using mag_type =
typename KokkosKernels::ArithTraits<SC>::mag_type;
119 using mag_type =
typename Kokkos::ArithTraits<SC>::mag_type;
121 const char prefix[] =
"leftAndOrRightScaleCrsMatrix: ";
124 Kokkos::View<const mag_type*, device_type> row_lcl;
125 Kokkos::View<const mag_type*, device_type> col_lcl;
128 const bool same = rowScalingFactors.getMap()->isSameAs(*(A.
getRowMap()));
129 TEUCHOS_TEST_FOR_EXCEPTION(!same, std::invalid_argument, prefix <<
"rowScalingFactors's Map "
130 "must be the same as the CrsMatrix's row Map. If you see this "
131 "message, it's likely that you are using a range Map Vector and that "
132 "the CrsMatrix's row Map is overlapping.");
137 auto row_lcl_2d = rowScalingFactors.template getLocalView<dev_memory_space>(Access::ReadOnly);
138 row_lcl = Kokkos::subview(row_lcl_2d, Kokkos::ALL(), 0);
142 const bool same = colScalingFactors.getMap()->isSameAs(*(A.
getColMap()));
143 TEUCHOS_TEST_FOR_EXCEPTION(!same, std::invalid_argument, prefix <<
"colScalingFactors's Map "
144 "must be the same as the CrsMatrix's column Map. If you see this "
145 "message, it's likely that you are using a domain Map Vector.");
150 auto col_lcl_2d = colScalingFactors.template getLocalView<dev_memory_space>(Access::ReadOnly);
151 col_lcl = Kokkos::subview(col_lcl_2d, Kokkos::ALL(), 0);
155 assumeSymmetric, scaling);
166 #if KOKKOS_VERSION >= 40799
167 #define TPETRA_LEFTANDORRIGHTSCALECRSMATRIX_INSTANT(SC, LO, GO, NT) \
169 leftAndOrRightScaleCrsMatrix( \
170 Tpetra::CrsMatrix<SC, LO, GO, NT>& A, \
171 const Kokkos::View< \
172 const KokkosKernels::ArithTraits<SC>::mag_type*, \
173 NT::device_type>& rowScalingFactors, \
174 const Kokkos::View< \
175 const KokkosKernels::ArithTraits<SC>::mag_type*, \
176 NT::device_type>& colScalingFactors, \
177 const bool leftScale, \
178 const bool rightScale, \
179 const bool assumeSymmetric, \
180 const EScaling scaling); \
183 leftAndOrRightScaleCrsMatrix( \
184 Tpetra::CrsMatrix<SC, LO, GO, NT>& A, \
185 const Tpetra::Vector<KokkosKernels::ArithTraits<SC>::mag_type, LO, GO, NT>& rowScalingFactors, \
186 const Tpetra::Vector<KokkosKernels::ArithTraits<SC>::mag_type, LO, GO, NT>& colScalingFactors, \
187 const bool leftScale, \
188 const bool rightScale, \
189 const bool assumeSymmetric, \
190 const EScaling scaling);
193 #define TPETRA_LEFTANDORRIGHTSCALECRSMATRIX_INSTANT(SC, LO, GO, NT) \
195 leftAndOrRightScaleCrsMatrix( \
196 Tpetra::CrsMatrix<SC, LO, GO, NT>& A, \
197 const Kokkos::View< \
198 const Kokkos::ArithTraits<SC>::mag_type*, \
199 NT::device_type>& rowScalingFactors, \
200 const Kokkos::View< \
201 const Kokkos::ArithTraits<SC>::mag_type*, \
202 NT::device_type>& colScalingFactors, \
203 const bool leftScale, \
204 const bool rightScale, \
205 const bool assumeSymmetric, \
206 const EScaling scaling); \
209 leftAndOrRightScaleCrsMatrix( \
210 Tpetra::CrsMatrix<SC, LO, GO, NT>& A, \
211 const Tpetra::Vector<Kokkos::ArithTraits<SC>::mag_type, LO, GO, NT>& rowScalingFactors, \
212 const Tpetra::Vector<Kokkos::ArithTraits<SC>::mag_type, LO, GO, NT>& colScalingFactors, \
213 const bool leftScale, \
214 const bool rightScale, \
215 const bool assumeSymmetric, \
216 const EScaling scaling);
220 #endif // TPETRA_LEFTANDORRIGHTSCALECRSMATRIX_DEF_HPP
Sparse matrix that presents a row-oriented interface that lets users read or modify entries...
Teuchos::RCP< const map_type > getRangeMap() const override
The range Map of this matrix.
void leftAndOrRightScaleCrsMatrix(Tpetra::CrsMatrix< SC, LO, GO, NT > &A, const Kokkos::View< const typename Kokkos::ArithTraits< SC >::mag_type *, typename NT::device_type > &rowScalingFactors, const Kokkos::View< const typename Kokkos::ArithTraits< SC >::mag_type *, typename NT::device_type > &colScalingFactors, const bool leftScale, const bool rightScale, const bool assumeSymmetric, const EScaling scaling)
Left-scale and/or right-scale (in that order) the entries of the input Tpetra::CrsMatrix A...
EScaling
Whether "scaling" a matrix means multiplying or dividing its entries.
void resumeFill(const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Resume operations that may change the values or structure of the matrix.
static bool debug()
Whether Tpetra is in debug mode.
void leftScaleLocalCrsMatrix(const LocalSparseMatrixType &A_lcl, const ScalingFactorsViewType &scalingFactors, const bool assumeSymmetric, const bool divide=true)
Left-scale a KokkosSparse::CrsMatrix.
bool isFillComplete() const override
Whether the matrix is fill complete.
TPETRA_DETAILS_ALWAYS_INLINE local_matrix_device_type getLocalMatrixDevice() const
The local sparse matrix.
void fillComplete(const Teuchos::RCP< const map_type > &domainMap, const Teuchos::RCP< const map_type > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Tell the matrix that you are done changing its structure or values, and that you are ready to do comp...
Declaration and definition of Tpetra::Details::leftScaleLocalCrsMatrix.
Teuchos::RCP< const map_type > getDomainMap() const override
The domain Map of this matrix.
A distributed dense vector.
Teuchos::RCP< const map_type > getColMap() const override
The Map that describes the column distribution in this matrix.
Declaration and definition of Tpetra::Details::rightScaleLocalCrsMatrix.
void rightScaleLocalCrsMatrix(const LocalSparseMatrixType &A_lcl, const ScalingFactorsViewType &scalingFactors, const bool assumeSymmetric, const bool divide=true)
Right-scale a KokkosSparse::CrsMatrix.
Declaration of Tpetra::Details::Behavior, a class that describes Tpetra's behavior.
Teuchos::RCP< const map_type > getRowMap() const override
The Map that describes the row distribution in this matrix.