42 #ifndef TPETRA_LEFTANDORRIGHTSCALECRSMATRIX_DEF_HPP
43 #define TPETRA_LEFTANDORRIGHTSCALECRSMATRIX_DEF_HPP
52 #include "Tpetra_CrsMatrix.hpp"
53 #include "Tpetra_Vector.hpp"
57 #include "Teuchos_TestForException.hpp"
61 template<
class SC,
class LO,
class GO,
class NT>
65 const typename Kokkos::ArithTraits<SC>::mag_type*,
66 typename NT::device_type>& rowScalingFactors,
68 const typename Kokkos::ArithTraits<SC>::mag_type*,
69 typename NT::device_type>& colScalingFactors,
71 const bool rightScale,
72 const bool assumeSymmetric,
75 if (! leftScale && ! rightScale) {
80 if (! A_fillComplete_on_input) {
88 static_cast<LO
> (A.
getRowMap ()->getLocalNumElements ());
89 TEUCHOS_TEST_FOR_EXCEPTION
90 (A_lcl.numRows () != lclNumRows, std::invalid_argument,
91 "leftAndOrRightScaleCrsMatrix: Local matrix is not valid. "
92 "This means that A was not created with a local matrix, "
93 "and that fillComplete has never yet been called on A before. "
94 "Please call fillComplete on A at least once first "
95 "before calling this method.");
101 const bool divide = scaling == SCALING_DIVIDE;
115 if (A_fillComplete_on_input) {
116 Teuchos::RCP<Teuchos::ParameterList> params = Teuchos::parameterList ();
117 params->set (
"No Nonlocal Changes",
true);
122 template<
class SC,
class LO,
class GO,
class NT>
126 typename Kokkos::ArithTraits<SC>::mag_type,
127 LO, GO, NT>& rowScalingFactors,
129 typename Kokkos::ArithTraits<SC>::mag_type,
130 LO, GO, NT>& colScalingFactors,
131 const bool leftScale,
132 const bool rightScale,
133 const bool assumeSymmetric,
136 using device_type =
typename NT::device_type;
137 using dev_memory_space =
typename device_type::memory_space;
138 using mag_type =
typename Kokkos::ArithTraits<SC>::mag_type;
139 const char prefix[] =
"leftAndOrRightScaleCrsMatrix: ";
142 Kokkos::View<const mag_type*, device_type> row_lcl;
143 Kokkos::View<const mag_type*, device_type> col_lcl;
146 const bool same = rowScalingFactors.getMap ()->isSameAs (* (A.
getRowMap ()));
147 TEUCHOS_TEST_FOR_EXCEPTION
148 (! same, std::invalid_argument, prefix <<
"rowScalingFactors's Map "
149 "must be the same as the CrsMatrix's row Map. If you see this "
150 "message, it's likely that you are using a range Map Vector and that "
151 "the CrsMatrix's row Map is overlapping.");
156 auto row_lcl_2d = rowScalingFactors.template getLocalView<dev_memory_space> (Access::ReadOnly);
157 row_lcl = Kokkos::subview (row_lcl_2d, Kokkos::ALL (), 0);
161 const bool same = colScalingFactors.getMap ()->isSameAs (* (A.
getColMap ()));
162 TEUCHOS_TEST_FOR_EXCEPTION
163 (! same, std::invalid_argument, prefix <<
"colScalingFactors's Map "
164 "must be the same as the CrsMatrix's column Map. If you see this "
165 "message, it's likely that you are using a domain Map Vector.");
170 auto col_lcl_2d = colScalingFactors.template getLocalView<dev_memory_space> (Access::ReadOnly);
171 col_lcl = Kokkos::subview (col_lcl_2d, Kokkos::ALL (), 0);
175 assumeSymmetric, scaling);
186 #define TPETRA_LEFTANDORRIGHTSCALECRSMATRIX_INSTANT(SC,LO,GO,NT) \
188 leftAndOrRightScaleCrsMatrix ( \
189 Tpetra::CrsMatrix<SC, LO, GO, NT>& A, \
190 const Kokkos::View< \
191 const Kokkos::ArithTraits<SC>::mag_type*, \
192 NT::device_type>& rowScalingFactors, \
193 const Kokkos::View< \
194 const Kokkos::ArithTraits<SC>::mag_type*, \
195 NT::device_type>& colScalingFactors, \
196 const bool leftScale, \
197 const bool rightScale, \
198 const bool assumeSymmetric, \
199 const EScaling scaling); \
202 leftAndOrRightScaleCrsMatrix ( \
203 Tpetra::CrsMatrix<SC, LO, GO, NT>& A, \
204 const Tpetra::Vector<Kokkos::ArithTraits<SC>::mag_type, LO, GO, NT>& rowScalingFactors, \
205 const Tpetra::Vector<Kokkos::ArithTraits<SC>::mag_type, LO, GO, NT>& colScalingFactors, \
206 const bool leftScale, \
207 const bool rightScale, \
208 const bool assumeSymmetric, \
209 const EScaling scaling);
211 #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.
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.