43 #ifndef IFPACK2_IDENTITY_SOLVER_DEF_HPP
44 #define IFPACK2_IDENTITY_SOLVER_DEF_HPP
46 #include "Ifpack2_IdentitySolver_decl.hpp"
50 template<
class MatrixType>
54 isInitialized_ (false),
65 template<
class MatrixType>
70 template<
class MatrixType>
75 template<
class MatrixType>
79 matrix_.is_null (), std::runtime_error,
"Ifpack2::IdentitySolver: "
80 "You must call setMatrix() with a nonnull input matrix "
81 "before you may call initialize() or compute().");
84 ! matrix_->getDomainMap ()->isCompatible (* (matrix_->getRangeMap ())),
85 std::invalid_argument,
86 "Ifpack2::IdentitySolver: The domain and range Maps "
87 "of the input matrix must be compatible.");
92 if (! matrix_->getDomainMap ()->isSameAs (* (matrix_->getRangeMap ()))) {
93 export_ =
Teuchos::rcp (
new export_type (matrix_->getDomainMap (),
94 matrix_->getRangeMap ()));
99 export_ = Teuchos::null;
102 isInitialized_ =
true;
106 template<
class MatrixType>
110 matrix_.is_null (), std::runtime_error,
"Ifpack2::IdentitySolver: "
111 "You must call setMatrix() with a nonnull input matrix "
112 "before you may call initialize() or compute().");
114 if (! isInitialized_) {
122 template<
class MatrixType>
124 apply (
const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
125 Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
136 ! isComputed (), std::runtime_error,
137 "Ifpack2::IdentitySolver::apply: If compute() has not yet been called, "
138 "or if you have changed the matrix via setMatrix(), "
139 "you must call compute() before you may call this method.");
144 if (export_.is_null ()) {
145 Y.update (alpha, X, beta);
148 if (alpha == STS::one () && beta == STS::zero ()) {
149 Y.doExport (X, *export_, Tpetra::REPLACE);
155 MV X_tmp (Y.getMap (), Y.getNumVectors ());
156 X_tmp.doExport (X, *export_, Tpetra::REPLACE);
157 Y.update (alpha, X_tmp, beta);
163 template <
class MatrixType>
165 return(numInitialize_);
168 template <
class MatrixType>
173 template <
class MatrixType>
178 template <
class MatrixType>
180 return(initializeTime_);
183 template<
class MatrixType>
185 return(computeTime_);
188 template<
class MatrixType>
193 template <
class MatrixType>
196 std::ostringstream os;
201 os <<
"\"Ifpack2::IdentitySolver\": {";
202 if (this->getObjectLabel () !=
"") {
203 os <<
"Label: \"" << this->getObjectLabel () <<
"\", ";
205 os <<
"Initialized: " << (isInitialized () ?
"true" :
"false") <<
", "
206 <<
"Computed: " << (isComputed () ?
"true" :
"false") <<
", ";
208 if (matrix_.is_null ()) {
209 os <<
"Matrix: null";
212 os <<
"Matrix: not null"
213 <<
", Global matrix dimensions: ["
214 << matrix_->getGlobalNumRows () <<
", "
215 << matrix_->getGlobalNumCols () <<
"]";
222 template <
class MatrixType>
234 out <<
"\"Ifpack2::IdentitySolver\":" << endl;
237 out <<
"numInitialize: " << numInitialize_ << endl;
238 out <<
"numCompute: " << numCompute_ << endl;
239 out <<
"numApply: " << numApply_ << endl;
243 template <
class MatrixType>
247 matrix_.is_null (), std::runtime_error,
"Ifpack2::IdentitySolver::getDomainMap: "
248 "The matrix is null. Please call setMatrix() with a nonnull input "
249 "before calling this method.");
250 return matrix_->getDomainMap ();
253 template <
class MatrixType>
257 matrix_.is_null (), std::runtime_error,
"Ifpack2::IdentitySolver::getRangeMap: "
258 "The matrix is null. Please call setMatrix() with a nonnull input "
259 "before calling this method.");
260 return matrix_->getRangeMap ();
263 template<
class MatrixType>
269 ! A.
is_null () && A->getComm ()->getSize () == 1 &&
270 A->getNodeNumRows () != A->getNodeNumCols (),
271 std::runtime_error,
"Ifpack2::IdentitySolver::setMatrix: If A's communicator only "
272 "contains one process, then A must be square. Instead, you provided a "
273 "matrix A with " << A->getNodeNumRows () <<
" rows and "
274 << A->getNodeNumCols () <<
" columns.");
279 isInitialized_ =
false;
281 export_ = Teuchos::null;
288 #define IFPACK2_IDENTITYSOLVER_INSTANT(S,LO,GO,N) \
289 template class Ifpack2::IdentitySolver< Tpetra::RowMatrix<S, LO, GO, N> >;
291 #endif // IFPACK2_IDENTITY_SOLVER_DEF_HPP
int getNumApply() const
Return the number of calls to apply().
Definition: Ifpack2_IdentitySolver_def.hpp:174
MatrixType::node_type node_type
Node type of the input matrix.
Definition: Ifpack2_IdentitySolver_decl.hpp:77
double getComputeTime() const
Return the time spent in compute().
Definition: Ifpack2_IdentitySolver_def.hpp:184
int getNumCompute() const
Return the number of calls to compute().
Definition: Ifpack2_IdentitySolver_def.hpp:169
virtual ~IdentitySolver()
Destructor.
Definition: Ifpack2_IdentitySolver_def.hpp:66
void apply(const Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &X, Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const
Apply the preconditioner to X, and put the result in Y.
Definition: Ifpack2_IdentitySolver_def.hpp:124
Teuchos::RCP< const map_type > getRangeMap() const
Return the Tpetra::Map object associated with the range of this operator.
Definition: Ifpack2_IdentitySolver_def.hpp:254
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void setParameters(const Teuchos::ParameterList ¶ms)
Set this object's parameters.
Definition: Ifpack2_IdentitySolver_def.hpp:71
void initialize()
Initialize.
Definition: Ifpack2_IdentitySolver_def.hpp:76
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< const map_type > getDomainMap() const
Return the Tpetra::Map object associated with the domain of this operator.
Definition: Ifpack2_IdentitySolver_def.hpp:244
IdentitySolver(const Teuchos::RCP< const row_matrix_type > &A)
Constructor: Takes the matrix to precondition.
Definition: Ifpack2_IdentitySolver_def.hpp:52
MatrixType::local_ordinal_type local_ordinal_type
Type of the local indices of the input matrix.
Definition: Ifpack2_IdentitySolver_decl.hpp:73
virtual void setMatrix(const Teuchos::RCP< const row_matrix_type > &A)
Set this preconditioner's matrix.
Definition: Ifpack2_IdentitySolver_def.hpp:265
std::string description() const
Return a simple one-line description of this object.
Definition: Ifpack2_IdentitySolver_def.hpp:194
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
Definition: Ifpack2_IdentitySolver_def.hpp:224
void compute()
Compute the preconditioner.
Definition: Ifpack2_IdentitySolver_def.hpp:107
int getNumInitialize() const
Return the number of calls to initialize().
Definition: Ifpack2_IdentitySolver_def.hpp:164
MatrixType::global_ordinal_type global_ordinal_type
Type of the global indices of the input matrix.
Definition: Ifpack2_IdentitySolver_decl.hpp:75
double getInitializeTime() const
Return the time spent in initialize().
Definition: Ifpack2_IdentitySolver_def.hpp:179
MatrixType::scalar_type scalar_type
Type of the entries of the input matrix.
Definition: Ifpack2_IdentitySolver_decl.hpp:71
double getApplyTime() const
Return the time spent in apply().
Definition: Ifpack2_IdentitySolver_def.hpp:189
static std::string name()