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"
61template<
class SC,
class LO,
class GO,
class NT>
65 const typename Kokkos::ArithTraits<SC>::mag_type*,
68 const typename Kokkos::ArithTraits<SC>::mag_type*,
71 const bool rightScale,
72 const bool assumeSymmetric,
75 if (! leftScale && ! rightScale) {
86 auto A_lcl =
A.getLocalMatrixDevice ();
88 static_cast<LO
> (
A.getRowMap ()->getLocalNumElements ());
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.");
116 Teuchos::RCP<Teuchos::ParameterList>
params = Teuchos::parameterList ();
117 params->set (
"No Nonlocal Changes",
true);
118 A.fillComplete (
A.getDomainMap (),
A.getRangeMap (),
params);
122template<
class SC,
class LO,
class GO,
class NT>
126 typename Kokkos::ArithTraits<SC>::mag_type,
129 typename Kokkos::ArithTraits<SC>::mag_type,
131 const bool leftScale,
132 const bool rightScale,
133 const bool assumeSymmetric,
136 using device_type =
typename NT::device_type;
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;
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.");
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.");
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);
Declaration of Tpetra::Details::Behavior, a class that describes Tpetra's behavior.
Declaration and definition of Tpetra::Details::leftScaleLocalCrsMatrix.
Declaration and definition of Tpetra::Details::rightScaleLocalCrsMatrix.
Struct that holds views of the contents of a CrsMatrix.
static bool debug()
Whether Tpetra is in debug mode.
A distributed dense vector.
void leftScaleLocalCrsMatrix(const LocalSparseMatrixType &A_lcl, const ScalingFactorsViewType &scalingFactors, const bool assumeSymmetric, const bool divide=true)
Left-scale a KokkosSparse::CrsMatrix.
void rightScaleLocalCrsMatrix(const LocalSparseMatrixType &A_lcl, const ScalingFactorsViewType &scalingFactors, const bool assumeSymmetric, const bool divide=true)
Right-scale a KokkosSparse::CrsMatrix.
Namespace Tpetra contains the class and methods constituting the Tpetra library.
EScaling
Whether "scaling" a matrix means multiplying or dividing its entries.
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.