43 #ifndef IFPACK2_HIPTMAIR_DEF_HPP 44 #define IFPACK2_HIPTMAIR_DEF_HPP 46 #include "Ifpack2_Details_OneLevelFactory.hpp" 47 #include "Ifpack2_Parameters.hpp" 48 #include "Teuchos_TimeMonitor.hpp" 49 #include "Tpetra_MultiVector.hpp" 56 template <
class MatrixType>
58 Hiptmair (
const Teuchos::RCP<const row_matrix_type>& A,
59 const Teuchos::RCP<const row_matrix_type>& PtAP,
60 const Teuchos::RCP<const row_matrix_type>& P) :
65 precType1_ (
"CHEBYSHEV"),
66 precType2_ (
"CHEBYSHEV"),
68 ZeroStartingSolution_ (true),
70 IsInitialized_ (false),
75 InitializeTime_ (0.0),
83 template <
class MatrixType>
85 Hiptmair (
const Teuchos::RCP<const row_matrix_type>& A):
90 precType1_ (
"CHEBYSHEV"),
91 precType2_ (
"CHEBYSHEV"),
93 ZeroStartingSolution_ (true),
95 IsInitialized_ (false),
100 InitializeTime_ (0.0),
106 template <
class MatrixType>
109 template <
class MatrixType>
114 using Teuchos::ParameterList;
115 using Teuchos::Exceptions::InvalidParameterName;
116 using Teuchos::Exceptions::InvalidParameterType;
118 ParameterList params = plist;
125 std::string precType1 = precType1_;
126 std::string precType2 = precType2_;
127 std::string preOrPost = preOrPost_;
128 Teuchos::ParameterList precList1 = precList1_;
129 Teuchos::ParameterList precList2 = precList2_;
130 bool zeroStartingSolution = ZeroStartingSolution_;
132 precType1 = params.get(
"hiptmair: smoother type 1", precType1);
133 precType2 = params.get(
"hiptmair: smoother type 2", precType2);
134 precList1 = params.get(
"hiptmair: smoother list 1", precList1);
135 precList2 = params.get(
"hiptmair: smoother list 2", precList2);
136 preOrPost = params.get(
"hiptmair: pre or post", preOrPost);
137 zeroStartingSolution = params.get(
"hiptmair: zero starting solution",
138 zeroStartingSolution);
143 PtAP_ = params.get<RCP<row_matrix_type> >(
"PtAP");
145 P_ = params.get<RCP<row_matrix_type> >(
"P");
149 precType1_ = precType1;
150 precType2_ = precType2;
151 precList1_ = precList1;
152 precList2_ = precList2;
153 preOrPost_ = preOrPost;
154 ZeroStartingSolution_ = zeroStartingSolution;
158 template <
class MatrixType>
159 Teuchos::RCP<const Teuchos::Comm<int> >
161 TEUCHOS_TEST_FOR_EXCEPTION(
162 A_.is_null (), std::runtime_error,
"Ifpack2::Hiptmair::getComm: " 163 "The input matrix A is null. Please call setMatrix() with a nonnull " 164 "input matrix before calling this method.");
165 return A_->getComm ();
169 template <
class MatrixType>
170 Teuchos::RCP<const Tpetra::RowMatrix<typename MatrixType::scalar_type,typename MatrixType::local_ordinal_type,typename MatrixType::global_ordinal_type,typename MatrixType::node_type> >
176 template <
class MatrixType>
177 Teuchos::RCP<const Tpetra::Map<typename MatrixType::local_ordinal_type,typename MatrixType::global_ordinal_type,typename MatrixType::node_type> >
180 TEUCHOS_TEST_FOR_EXCEPTION(
181 A_.is_null (), std::runtime_error,
"Ifpack2::Hiptmair::getDomainMap: " 182 "The input matrix A is null. Please call setMatrix() with a nonnull " 183 "input matrix before calling this method.");
184 return A_->getDomainMap ();
188 template <
class MatrixType>
189 Teuchos::RCP<const Tpetra::Map<typename MatrixType::local_ordinal_type,typename MatrixType::global_ordinal_type,typename MatrixType::node_type> >
192 TEUCHOS_TEST_FOR_EXCEPTION(
193 A_.is_null (), std::runtime_error,
"Ifpack2::Hiptmair::getRangeMap: " 194 "The input matrix A is null. Please call setMatrix() with a nonnull " 195 "input matrix before calling this method.");
196 return A_->getRangeMap ();
200 template <
class MatrixType>
208 template <
class MatrixType>
210 return NumInitialize_;
214 template <
class MatrixType>
220 template <
class MatrixType>
226 template <
class MatrixType>
228 return InitializeTime_;
232 template <
class MatrixType>
238 template <
class MatrixType>
244 template <
class MatrixType>
247 using Teuchos::ParameterList;
251 TEUCHOS_TEST_FOR_EXCEPTION(
252 A_.is_null (), std::runtime_error,
"Ifpack2::Hiptmair::initialize: " 253 "The input matrix A is null. Please call setMatrix() with a nonnull " 254 "input matrix before calling this method.");
257 IsInitialized_ =
false;
260 Teuchos::Time timer (
"initialize");
261 double startTime = timer.wallTime();
263 Teuchos::TimeMonitor timeMon (timer);
267 ifpack2_prec1_=factory.
create(precType1_,A_);
268 ifpack2_prec1_->initialize();
269 ifpack2_prec1_->setParameters(precList1_);
271 ifpack2_prec2_=factory.
create(precType2_,PtAP_);
272 ifpack2_prec2_->initialize();
273 ifpack2_prec2_->setParameters(precList2_);
276 IsInitialized_ =
true;
278 InitializeTime_ += (timer.wallTime() - startTime);
282 template <
class MatrixType>
285 TEUCHOS_TEST_FOR_EXCEPTION(
286 A_.is_null (), std::runtime_error,
"Ifpack2::Hiptmair::compute: " 287 "The input matrix A is null. Please call setMatrix() with a nonnull " 288 "input matrix before calling this method.");
291 if (! isInitialized ()) {
295 Teuchos::Time timer (
"compute");
296 double startTime = timer.wallTime();
298 Teuchos::TimeMonitor timeMon (timer);
299 ifpack2_prec1_->compute();
300 ifpack2_prec2_->compute();
304 ComputeTime_ += (timer.wallTime() - startTime);
308 template <
class MatrixType>
310 apply (
const Tpetra::MultiVector<
typename MatrixType::scalar_type,
311 typename MatrixType::local_ordinal_type,
312 typename MatrixType::global_ordinal_type,
313 typename MatrixType::node_type>& X,
314 Tpetra::MultiVector<
typename MatrixType::scalar_type,
315 typename MatrixType::local_ordinal_type,
316 typename MatrixType::global_ordinal_type,
317 typename MatrixType::node_type>& Y,
318 Teuchos::ETransp mode,
319 typename MatrixType::scalar_type alpha,
320 typename MatrixType::scalar_type beta)
const 324 using Teuchos::rcpFromRef;
327 TEUCHOS_TEST_FOR_EXCEPTION(
328 ! isComputed (), std::runtime_error,
329 "Ifpack2::Hiptmair::apply: You must call compute() before you may call apply().");
330 TEUCHOS_TEST_FOR_EXCEPTION(
331 X.getNumVectors () != Y.getNumVectors (), std::invalid_argument,
332 "Ifpack2::Hiptmair::apply: The MultiVector inputs X and Y do not have the " 333 "same number of columns. X.getNumVectors() = " << X.getNumVectors ()
334 <<
" != Y.getNumVectors() = " << Y.getNumVectors () <<
".");
337 TEUCHOS_TEST_FOR_EXCEPTION(
338 alpha != STS::one (), std::logic_error,
339 "Ifpack2::Hiptmair::apply: alpha != 1 has not been implemented.");
340 TEUCHOS_TEST_FOR_EXCEPTION(
341 beta != STS::zero (), std::logic_error,
342 "Ifpack2::Hiptmair::apply: zero != 0 has not been implemented.");
343 TEUCHOS_TEST_FOR_EXCEPTION(
344 mode != Teuchos::NO_TRANS, std::logic_error,
345 "Ifpack2::Hiptmair::apply: mode != Teuchos::NO_TRANS has not been implemented.");
347 Teuchos::Time timer (
"apply");
348 double startTime = timer.wallTime();
350 Teuchos::TimeMonitor timeMon (timer);
357 Xcopy = rcp (
new MV (X, Teuchos::Copy));
359 Xcopy = rcpFromRef (X);
363 RCP<MV> Ycopy = rcpFromRef (Y);
364 if (ZeroStartingSolution_) {
365 Ycopy->putScalar (STS::zero ());
369 applyHiptmairSmoother (*Xcopy, *Ycopy);
373 ApplyTime_ += (timer.wallTime() - startTime);
376 template <
class MatrixType>
379 typename MatrixType::local_ordinal_type,
380 typename MatrixType::global_ordinal_type,
381 typename MatrixType::node_type>& X,
382 Tpetra::MultiVector<
typename MatrixType::scalar_type,
383 typename MatrixType::local_ordinal_type,
384 typename MatrixType::global_ordinal_type,
385 typename MatrixType::node_type>& Y)
const 389 using Teuchos::rcpFromRef;
390 typedef Tpetra::MultiVector<scalar_type, local_ordinal_type,
391 global_ordinal_type, node_type> MV;
392 const scalar_type ZERO = STS::zero ();
393 const scalar_type ONE = STS::one ();
395 RCP<MV> res1 = rcp (
new MV (A_->getRowMap (), X.getNumVectors ()));
396 RCP<MV> vec1 = rcp (
new MV (A_->getRowMap (), X.getNumVectors ()));
397 RCP<MV> res2 = rcp (
new MV (PtAP_->getRowMap (), X.getNumVectors ()));
398 RCP<MV> vec2 = rcp (
new MV (PtAP_->getRowMap (), X.getNumVectors ()));
400 if (preOrPost_ ==
"pre" || preOrPost_ ==
"both") {
402 A_->apply (Y, *res1);
403 res1->update (ONE, X, -ONE);
404 vec1->putScalar (ZERO);
405 ifpack2_prec1_->apply (*res1, *vec1);
406 Y.update (ONE, *vec1, ONE);
410 A_->apply (Y, *res1);
411 res1->update (ONE, X, -ONE);
412 P_->apply (*res1, *res2, Teuchos::TRANS);
413 vec2->putScalar (ZERO);
414 ifpack2_prec2_->apply (*res2, *vec2);
415 P_->apply (*vec2, *vec1, Teuchos::NO_TRANS);
416 Y.update (ONE,*vec1,ONE);
418 if (preOrPost_ ==
"post" || preOrPost_ ==
"both") {
420 A_->apply (Y, *res1);
421 res1->update (ONE, X, -ONE);
422 vec1->putScalar (ZERO);
423 ifpack2_prec1_->apply (*res1, *vec1);
424 Y.update (ONE, *vec1, ONE);
428 template <
class MatrixType>
431 std::ostringstream os;
436 os <<
"\"Ifpack2::Hiptmair\": {";
437 if (this->getObjectLabel () !=
"") {
438 os <<
"Label: \"" << this->getObjectLabel () <<
"\", ";
440 os <<
"Initialized: " << (isInitialized () ?
"true" :
"false") <<
", " 441 <<
"Computed: " << (isComputed () ?
"true" :
"false") <<
", ";
444 os <<
"Matrix: null";
447 os <<
"Matrix: not null" 448 <<
", Global matrix dimensions: [" 449 << A_->getGlobalNumRows () <<
", " << A_->getGlobalNumCols () <<
"]";
457 template <
class MatrixType>
460 const Teuchos::EVerbosityLevel verbLevel)
const 464 using Teuchos::VERB_DEFAULT;
465 using Teuchos::VERB_NONE;
466 using Teuchos::VERB_LOW;
467 using Teuchos::VERB_MEDIUM;
468 using Teuchos::VERB_HIGH;
469 using Teuchos::VERB_EXTREME;
471 const Teuchos::EVerbosityLevel vl =
472 (verbLevel == VERB_DEFAULT) ? VERB_LOW : verbLevel;
474 if (vl != VERB_NONE) {
476 Teuchos::OSTab tab0 (out);
477 out <<
"\"Ifpack2::Hiptmair\":";
479 Teuchos::OSTab tab1 (out);
480 if (this->getObjectLabel () !=
"") {
481 out <<
"Label: " << this->getObjectLabel () << endl;
483 out <<
"Initialized: " << (isInitialized () ?
"true" :
"false") << endl
484 <<
"Computed: " << (isComputed () ?
"true" :
"false") << endl
485 <<
"Global number of rows: " << A_->getGlobalNumRows () << endl
486 <<
"Global number of columns: " << A_->getGlobalNumCols () << endl
489 out <<
" null" << endl;
491 A_->describe (out, vl);
498 #define IFPACK2_HIPTMAIR_INSTANT(S,LO,GO,N) \ 499 template class Ifpack2::Hiptmair< Tpetra::RowMatrix<S, LO, GO, N> >; void initialize()
Do any initialization that depends on the input matrix's structure.
Definition: Ifpack2_Hiptmair_def.hpp:245
int getNumInitialize() const
Returns the number of calls to Initialize().
Definition: Ifpack2_Hiptmair_def.hpp:209
double getComputeTime() const
Returns the time spent in Compute().
Definition: Ifpack2_Hiptmair_def.hpp:233
MatrixType::local_ordinal_type local_ordinal_type
The type of local indices in the input MatrixType.
Definition: Ifpack2_Hiptmair_decl.hpp:85
Hiptmair(const Teuchos::RCP< const row_matrix_type > &A)
Constructor that takes 1 Tpetra matrix (assumes we'll get the rest off the parameter list) ...
Definition: Ifpack2_Hiptmair_def.hpp:85
double getApplyTime() const
Returns the time spent in apply().
Definition: Ifpack2_Hiptmair_def.hpp:239
MatrixType::scalar_type scalar_type
The type of the entries of the input MatrixType.
Definition: Ifpack2_Hiptmair_decl.hpp:82
Teuchos::RCP< const Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > > getRangeMap() const
Tpetra::Map representing the range of this operator.
Definition: Ifpack2_Hiptmair_def.hpp:190
int getNumCompute() const
Returns the number of calls to Compute().
Definition: Ifpack2_Hiptmair_def.hpp:215
void setParameters(const Teuchos::ParameterList ¶ms)
Set the preconditioner's parameters.
Definition: Ifpack2_Hiptmair_def.hpp:110
Teuchos::RCP< const Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > > getDomainMap() const
Tpetra::Map representing the domain of this operator.
Definition: Ifpack2_Hiptmair_def.hpp:178
Teuchos::RCP< const Tpetra::RowMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > > getMatrix() const
Returns a reference to the matrix to be preconditioned.
Definition: Ifpack2_Hiptmair_def.hpp:171
MatrixType::global_ordinal_type global_ordinal_type
The type of global indices in the input MatrixType.
Definition: Ifpack2_Hiptmair_decl.hpp:88
Wrapper for Hiptmair smoothers.
Definition: Ifpack2_Hiptmair_decl.hpp:71
MatrixType::node_type node_type
The Node type used by the input MatrixType.
Definition: Ifpack2_Hiptmair_decl.hpp:91
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_Hiptmair_def.hpp:459
void compute()
Do any initialization that depends on the input matrix's values.
Definition: Ifpack2_Hiptmair_def.hpp:283
virtual ~Hiptmair()
Destructor.
Definition: Ifpack2_Hiptmair_def.hpp:107
double getInitializeTime() const
Returns the time spent in Initialize().
Definition: Ifpack2_Hiptmair_def.hpp:227
std::string description() const
Return a simple one-line description of this object.
Definition: Ifpack2_Hiptmair_def.hpp:429
bool hasTransposeApply() const
Whether this object's apply() method can apply the transpose (or conjugate transpose, if applicable).
Definition: Ifpack2_Hiptmair_def.hpp:201
int getNumApply() const
Returns the number of calls to apply().
Definition: Ifpack2_Hiptmair_def.hpp:221
"Factory" for creating single-level preconditioners.
Definition: Ifpack2_Details_OneLevelFactory_decl.hpp:124
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:73
Teuchos::RCP< prec_type > create(const std::string &precType, const Teuchos::RCP< const row_matrix_type > &matrix) const
Create an instance of Preconditioner given the string name of the preconditioner type.
Definition: Ifpack2_Details_OneLevelFactory_def.hpp:84
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
Returns the operator's communicator.
Definition: Ifpack2_Hiptmair_def.hpp:160
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, putting the result in Y.
Definition: Ifpack2_Hiptmair_def.hpp:310