43 #ifndef IFPACK2_DIAGONAL_DEF_HPP 44 #define IFPACK2_DIAGONAL_DEF_HPP 46 #include "Ifpack2_Diagonal_decl.hpp" 47 #include "Tpetra_CrsMatrix.hpp" 51 template<
class MatrixType>
52 Diagonal<MatrixType>::Diagonal (
const Teuchos::RCP<const row_matrix_type>& A) :
54 initializeTime_ (0.0),
60 isInitialized_ (false),
64 template<
class MatrixType>
65 Diagonal<MatrixType>::Diagonal (
const Teuchos::RCP<const crs_matrix_type>& A) :
67 initializeTime_ (0.0),
73 isInitialized_ (false),
77 template<
class MatrixType>
78 Diagonal<MatrixType>::Diagonal (
const Teuchos::RCP<const vector_type>& diag) :
79 userInverseDiag_ (diag),
81 initializeTime_ (0.0),
87 isInitialized_ (false),
91 template<
class MatrixType>
92 Diagonal<MatrixType>::~Diagonal ()
95 template<
class MatrixType>
96 Teuchos::RCP<const typename Diagonal<MatrixType>::map_type>
97 Diagonal<MatrixType>::getDomainMap ()
const 99 if (matrix_.is_null ()) {
100 if (userInverseDiag_.is_null ()) {
101 TEUCHOS_TEST_FOR_EXCEPTION(
102 true, std::runtime_error,
"Ifpack2::Diagonal::getDomainMap: " 103 "The input matrix A is null, and you did not provide a vector of " 104 "inverse diagonal entries. Please call setMatrix() with a nonnull " 105 "input matrix before calling this method.");
107 return userInverseDiag_->getMap ();
110 return matrix_->getDomainMap ();
114 template<
class MatrixType>
115 Teuchos::RCP<const typename Diagonal<MatrixType>::map_type>
116 Diagonal<MatrixType>::getRangeMap ()
const 118 if (matrix_.is_null ()) {
119 if (userInverseDiag_.is_null ()) {
120 TEUCHOS_TEST_FOR_EXCEPTION(
121 true, std::runtime_error,
"Ifpack2::Diagonal::getRangeMap: " 122 "The input matrix A is null, and you did not provide a vector of " 123 "inverse diagonal entries. Please call setMatrix() with a nonnull " 124 "input matrix before calling this method.");
126 return userInverseDiag_->getMap ();
129 return matrix_->getRangeMap ();
133 template<
class MatrixType>
134 void Diagonal<MatrixType>::
135 setParameters (
const Teuchos::ParameterList& )
138 template<
class MatrixType>
139 void Diagonal<MatrixType>::reset ()
141 inverseDiag_ = Teuchos::null;
142 offsets_ = Teuchos::null;
143 isInitialized_ =
false;
147 template<
class MatrixType>
148 void Diagonal<MatrixType>::
149 setMatrix (
const Teuchos::RCP<const row_matrix_type>& A)
151 if (A.getRawPtr () != matrix_.getRawPtr ()) {
157 template<
class MatrixType>
158 void Diagonal<MatrixType>::initialize ()
162 TEUCHOS_TEST_FOR_EXCEPTION(
163 matrix_.is_null () && userInverseDiag_.is_null (), std::runtime_error,
164 "Ifpack2::Diagonal::initialize: The matrix to precondition is null, " 165 "and you did not provide a Tpetra::Vector of diagonal entries. " 166 "Please call setMatrix() with a nonnull input before calling this method.");
174 if (! matrix_.is_null ()) {
179 Teuchos::RCP<const crs_matrix_type> A_crs =
180 Teuchos::rcp_dynamic_cast<
const crs_matrix_type> (matrix_);
181 if (A_crs.is_null ()) {
182 offsets_ = Teuchos::null;
185 A_crs->getLocalDiagOffsets (offsets_);
189 isInitialized_ =
true;
193 template<
class MatrixType>
194 void Diagonal<MatrixType>::compute ()
198 TEUCHOS_TEST_FOR_EXCEPTION(
199 matrix_.is_null () && userInverseDiag_.is_null (), std::runtime_error,
200 "Ifpack2::Diagonal::compute: The matrix to precondition is null, " 201 "and you did not provide a Tpetra::Vector of diagonal entries. " 202 "Please call setMatrix() with a nonnull input before calling this method.");
204 if (! isInitialized_) {
214 if (matrix_.is_null ()) {
215 inverseDiag_ = userInverseDiag_;
218 Teuchos::RCP<vector_type> tmpVec (
new vector_type (matrix_->getRowMap ()));
219 Teuchos::RCP<const crs_matrix_type> A_crs =
220 Teuchos::rcp_dynamic_cast<
const crs_matrix_type> (matrix_);
221 if (A_crs.is_null ()) {
223 matrix_->getLocalDiagCopy (*tmpVec);
228 A_crs->getLocalDiagCopy (*tmpVec, offsets_ ());
230 tmpVec->reciprocal (*tmpVec);
231 inverseDiag_ = tmpVec;
238 template<
class MatrixType>
239 void Diagonal<MatrixType>::
240 apply (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
241 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
244 scalar_type beta)
const 246 TEUCHOS_TEST_FOR_EXCEPTION(
247 !
isComputed (), std::runtime_error,
"Ifpack2::Diagonal::apply: You " 248 "must first call compute() before you may call apply(). Once you have " 249 "called compute(), you need not call it again unless the values in the " 250 "matrix have changed, or unless you have called setMatrix().");
257 Y.elementWiseMultiply (alpha, *inverseDiag_, X, beta);
261 template <
class MatrixType>
262 int Diagonal<MatrixType>::getNumInitialize()
const {
263 return numInitialize_;
266 template <
class MatrixType>
267 int Diagonal<MatrixType>::getNumCompute()
const {
271 template <
class MatrixType>
272 int Diagonal<MatrixType>::getNumApply()
const {
276 template <
class MatrixType>
277 double Diagonal<MatrixType>::getInitializeTime()
const {
278 return initializeTime_;
281 template<
class MatrixType>
282 double Diagonal<MatrixType>::getComputeTime()
const {
286 template<
class MatrixType>
287 double Diagonal<MatrixType>::getApplyTime()
const {
291 template <
class MatrixType>
292 std::string Diagonal<MatrixType>::description ()
const 294 std::ostringstream out;
299 out <<
"\"Ifpack2::Diagonal\": " 301 if (this->getObjectLabel () !=
"") {
302 out <<
"Label: \"" << this->getObjectLabel () <<
"\", ";
304 if (matrix_.is_null ()) {
305 out <<
"Matrix: null";
308 out <<
"Matrix: not null" 309 <<
", Global matrix dimensions: [" 310 << matrix_->getGlobalNumRows () <<
", " 311 << matrix_->getGlobalNumCols () <<
"]";
318 template <
class MatrixType>
319 void Diagonal<MatrixType>::
320 describe (Teuchos::FancyOStream &out,
321 const Teuchos::EVerbosityLevel verbLevel)
const 325 const Teuchos::EVerbosityLevel vl =
326 (verbLevel == Teuchos::VERB_DEFAULT) ? Teuchos::VERB_LOW : verbLevel;
327 if (vl != Teuchos::VERB_NONE) {
328 Teuchos::OSTab tab0 (out);
329 out <<
"\"Ifpack2::Diagonal\":";
330 Teuchos::OSTab tab1 (out);
331 out <<
"Template parameter: " 332 << Teuchos::TypeNameTraits<MatrixType>::name () << endl;
333 if (this->getObjectLabel () !=
"") {
334 out <<
"Label: \"" << this->getObjectLabel () <<
"\", ";
336 out <<
"Number of initialize calls: " << numInitialize_ << endl
337 <<
"Number of compute calls: " << numCompute_ << endl
338 <<
"Number of apply calls: " << numApply_ << endl;
344 #define IFPACK2_DIAGONAL_INSTANT(S,LO,GO,N) \ 345 template class Ifpack2::Diagonal< Tpetra::RowMatrix<S, LO, GO, N> >; virtual bool isComputed() const
Returns true if partitions have been computed successfully.
Definition: Ifpack2_OverlappingPartitioner_def.hpp:357
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:72