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" 55 template <
class MatrixType>
57 Hiptmair (
const Teuchos::RCP<const row_matrix_type>& A,
58 const Teuchos::RCP<const row_matrix_type>& PtAP,
59 const Teuchos::RCP<const row_matrix_type>& P) :
64 precType1_ (
"CHEBYSHEV"),
65 precType2_ (
"CHEBYSHEV"),
67 ZeroStartingSolution_ (true),
69 IsInitialized_ (false),
74 InitializeTime_ (0.0),
80 template <
class MatrixType>
83 template <
class MatrixType>
87 using Teuchos::ParameterList;
88 using Teuchos::Exceptions::InvalidParameterName;
89 using Teuchos::Exceptions::InvalidParameterType;
91 ParameterList params = plist;
98 std::string precType1 = precType1_;
99 std::string precType2 = precType2_;
100 std::string preOrPost = preOrPost_;
101 Teuchos::ParameterList precList1 = precList1_;
102 Teuchos::ParameterList precList2 = precList2_;
103 bool zeroStartingSolution = ZeroStartingSolution_;
105 precType1 = params.get(
"hiptmair: smoother type 1", precType1);
106 precType2 = params.get(
"hiptmair: smoother type 2", precType2);
107 precList1 = params.get(
"hiptmair: smoother list 1", precList1);
108 precList2 = params.get(
"hiptmair: smoother list 2", precList2);
109 preOrPost = params.get(
"hiptmair: pre or post", preOrPost);
110 zeroStartingSolution = params.get(
"hiptmair: zero starting solution",
111 zeroStartingSolution);
114 precType1_ = precType1;
115 precType2_ = precType2;
116 precList1_ = precList1;
117 precList2_ = precList2;
118 preOrPost_ = preOrPost;
119 ZeroStartingSolution_ = zeroStartingSolution;
123 template <
class MatrixType>
124 Teuchos::RCP<const Teuchos::Comm<int> >
126 TEUCHOS_TEST_FOR_EXCEPTION(
127 A_.is_null (), std::runtime_error,
"Ifpack2::Hiptmair::getComm: " 128 "The input matrix A is null. Please call setMatrix() with a nonnull " 129 "input matrix before calling this method.");
130 return A_->getComm ();
134 template <
class MatrixType>
135 Teuchos::RCP<const Tpetra::RowMatrix<typename MatrixType::scalar_type,typename MatrixType::local_ordinal_type,typename MatrixType::global_ordinal_type,typename MatrixType::node_type> >
141 template <
class MatrixType>
142 Teuchos::RCP<const Tpetra::Map<typename MatrixType::local_ordinal_type,typename MatrixType::global_ordinal_type,typename MatrixType::node_type> >
145 TEUCHOS_TEST_FOR_EXCEPTION(
146 A_.is_null (), std::runtime_error,
"Ifpack2::Hiptmair::getDomainMap: " 147 "The input matrix A is null. Please call setMatrix() with a nonnull " 148 "input matrix before calling this method.");
149 return A_->getDomainMap ();
153 template <
class MatrixType>
154 Teuchos::RCP<const Tpetra::Map<typename MatrixType::local_ordinal_type,typename MatrixType::global_ordinal_type,typename MatrixType::node_type> >
157 TEUCHOS_TEST_FOR_EXCEPTION(
158 A_.is_null (), std::runtime_error,
"Ifpack2::Hiptmair::getRangeMap: " 159 "The input matrix A is null. Please call setMatrix() with a nonnull " 160 "input matrix before calling this method.");
161 return A_->getRangeMap ();
165 template <
class MatrixType>
173 template <
class MatrixType>
175 return NumInitialize_;
179 template <
class MatrixType>
185 template <
class MatrixType>
191 template <
class MatrixType>
193 return InitializeTime_;
197 template <
class MatrixType>
203 template <
class MatrixType>
209 template <
class MatrixType>
212 using Teuchos::ParameterList;
216 TEUCHOS_TEST_FOR_EXCEPTION(
217 A_.is_null (), std::runtime_error,
"Ifpack2::Hiptmair::initialize: " 218 "The input matrix A is null. Please call setMatrix() with a nonnull " 219 "input matrix before calling this method.");
222 IsInitialized_ =
false;
225 Teuchos::Time timer (
"initialize");
227 Teuchos::TimeMonitor timeMon (timer);
231 ifpack2_prec1_=factory.
create(precType1_,A_);
232 ifpack2_prec1_->initialize();
233 ifpack2_prec1_->setParameters(precList1_);
235 ifpack2_prec2_=factory.
create(precType2_,PtAP_);
236 ifpack2_prec2_->initialize();
237 ifpack2_prec2_->setParameters(precList2_);
240 IsInitialized_ =
true;
242 InitializeTime_ += timer.totalElapsedTime ();
246 template <
class MatrixType>
249 TEUCHOS_TEST_FOR_EXCEPTION(
250 A_.is_null (), std::runtime_error,
"Ifpack2::Hiptmair::compute: " 251 "The input matrix A is null. Please call setMatrix() with a nonnull " 252 "input matrix before calling this method.");
259 Teuchos::Time timer (
"compute");
261 Teuchos::TimeMonitor timeMon (timer);
262 ifpack2_prec1_->compute();
263 ifpack2_prec2_->compute();
267 ComputeTime_ += timer.totalElapsedTime ();
271 template <
class MatrixType>
273 apply (
const Tpetra::MultiVector<
typename MatrixType::scalar_type,
274 typename MatrixType::local_ordinal_type,
275 typename MatrixType::global_ordinal_type,
276 typename MatrixType::node_type>& X,
277 Tpetra::MultiVector<
typename MatrixType::scalar_type,
278 typename MatrixType::local_ordinal_type,
279 typename MatrixType::global_ordinal_type,
280 typename MatrixType::node_type>& Y,
281 Teuchos::ETransp mode,
282 typename MatrixType::scalar_type alpha,
283 typename MatrixType::scalar_type beta)
const 287 using Teuchos::rcpFromRef;
290 TEUCHOS_TEST_FOR_EXCEPTION(
292 "Ifpack2::Hiptmair::apply: You must call compute() before you may call apply().");
293 TEUCHOS_TEST_FOR_EXCEPTION(
294 X.getNumVectors () != Y.getNumVectors (), std::invalid_argument,
295 "Ifpack2::Hiptmair::apply: The MultiVector inputs X and Y do not have the " 296 "same number of columns. X.getNumVectors() = " << X.getNumVectors ()
297 <<
" != Y.getNumVectors() = " << Y.getNumVectors () <<
".");
300 TEUCHOS_TEST_FOR_EXCEPTION(
301 alpha != STS::one (), std::logic_error,
302 "Ifpack2::Hiptmair::apply: alpha != 1 has not been implemented.");
303 TEUCHOS_TEST_FOR_EXCEPTION(
304 beta != STS::zero (), std::logic_error,
305 "Ifpack2::Hiptmair::apply: zero != 0 has not been implemented.");
306 TEUCHOS_TEST_FOR_EXCEPTION(
307 mode != Teuchos::NO_TRANS, std::logic_error,
308 "Ifpack2::Hiptmair::apply: mode != Teuchos::NO_TRANS has not been implemented.");
310 Teuchos::Time timer (
"apply");
312 Teuchos::TimeMonitor timeMon (timer);
318 auto X_lcl_host = X.template getLocalView<Kokkos::HostSpace> ();
319 auto Y_lcl_host = Y.template getLocalView<Kokkos::HostSpace> ();
320 if (X_lcl_host.ptr_on_device () == Y_lcl_host.ptr_on_device ()) {
321 Xcopy = rcp (
new MV (X, Teuchos::Copy));
323 Xcopy = rcpFromRef (X);
327 RCP<MV> Ycopy = rcpFromRef (Y);
328 if (ZeroStartingSolution_) {
329 Ycopy->putScalar (STS::zero ());
333 applyHiptmairSmoother (*Xcopy, *Ycopy);
337 ApplyTime_ += timer.totalElapsedTime ();
340 template <
class MatrixType>
343 typename MatrixType::local_ordinal_type,
344 typename MatrixType::global_ordinal_type,
345 typename MatrixType::node_type>& X,
346 Tpetra::MultiVector<
typename MatrixType::scalar_type,
347 typename MatrixType::local_ordinal_type,
348 typename MatrixType::global_ordinal_type,
349 typename MatrixType::node_type>& Y)
const 353 using Teuchos::rcpFromRef;
356 const scalar_type ZERO = STS::zero ();
357 const scalar_type ONE = STS::one ();
359 RCP<MV> res1 = rcp (
new MV (A_->getRowMap (), X.getNumVectors ()));
360 RCP<MV> vec1 = rcp (
new MV (A_->getRowMap (), X.getNumVectors ()));
361 RCP<MV> res2 = rcp (
new MV (PtAP_->getRowMap (), X.getNumVectors ()));
362 RCP<MV> vec2 = rcp (
new MV (PtAP_->getRowMap (), X.getNumVectors ()));
364 if (preOrPost_ ==
"pre" || preOrPost_ ==
"both") {
366 A_->apply (Y, *res1);
367 res1->update (ONE, X, -ONE);
368 vec1->putScalar (ZERO);
369 ifpack2_prec1_->apply (*res1, *vec1);
370 Y.update (ONE, *vec1, ONE);
374 A_->apply (Y, *res1);
375 res1->update (ONE, X, -ONE);
376 P_->apply (*res1, *res2, Teuchos::TRANS);
377 vec2->putScalar (ZERO);
378 ifpack2_prec2_->apply (*res2, *vec2);
379 P_->apply (*vec2, *vec1, Teuchos::NO_TRANS);
380 Y.update (ONE,*vec1,ONE);
382 if (preOrPost_ ==
"post" || preOrPost_ ==
"both") {
384 A_->apply (Y, *res1);
385 res1->update (ONE, X, -ONE);
386 vec1->putScalar (ZERO);
387 ifpack2_prec1_->apply (*res1, *vec1);
388 Y.update (ONE, *vec1, ONE);
392 template <
class MatrixType>
395 std::ostringstream os;
400 os <<
"\"Ifpack2::Hiptmair\": {";
401 if (this->getObjectLabel () !=
"") {
402 os <<
"Label: \"" << this->getObjectLabel () <<
"\", ";
404 os <<
"Initialized: " << (
isInitialized () ?
"true" :
"false") <<
", " 405 <<
"Computed: " << (
isComputed () ?
"true" :
"false") <<
", ";
408 os <<
"Matrix: null";
411 os <<
"Matrix: not null" 412 <<
", Global matrix dimensions: [" 413 << A_->getGlobalNumRows () <<
", " << A_->getGlobalNumCols () <<
"]";
421 template <
class MatrixType>
424 const Teuchos::EVerbosityLevel verbLevel)
const 428 using Teuchos::VERB_DEFAULT;
429 using Teuchos::VERB_NONE;
430 using Teuchos::VERB_LOW;
431 using Teuchos::VERB_MEDIUM;
432 using Teuchos::VERB_HIGH;
433 using Teuchos::VERB_EXTREME;
435 const Teuchos::EVerbosityLevel vl =
436 (verbLevel == VERB_DEFAULT) ? VERB_LOW : verbLevel;
438 if (vl != VERB_NONE) {
440 Teuchos::OSTab tab0 (out);
441 out <<
"\"Ifpack2::Hiptmair\":";
443 Teuchos::OSTab tab1 (out);
444 if (this->getObjectLabel () !=
"") {
445 out <<
"Label: " << this->getObjectLabel () << endl;
447 out <<
"Initialized: " << (
isInitialized () ?
"true" :
"false") << endl
448 <<
"Computed: " << (
isComputed () ?
"true" :
"false") << endl
449 <<
"Global number of rows: " << A_->getGlobalNumRows () << endl
450 <<
"Global number of columns: " << A_->getGlobalNumCols () << endl
453 out <<
" null" << endl;
455 A_->describe (out, vl);
462 #define IFPACK2_HIPTMAIR_INSTANT(S,LO,GO,N) \ 463 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:210
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:68
MatrixType::local_ordinal_type local_ordinal_type
The type of local indices in the input MatrixType.
Definition: Ifpack2_Hiptmair_decl.hpp:84
MatrixType::scalar_type scalar_type
The type of the entries of the input MatrixType.
Definition: Ifpack2_Hiptmair_decl.hpp:81
double getComputeTime() const
Returns the time spent in Compute().
Definition: Ifpack2_Hiptmair_def.hpp:198
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
Returns the operator's communicator.
Definition: Ifpack2_Hiptmair_def.hpp:125
bool hasTransposeApply() const
Whether this object's apply() method can apply the transpose (or conjugate transpose, if applicable).
Definition: Ifpack2_Hiptmair_def.hpp:166
bool isInitialized() const
Return true if initialize() completed successfully, else false.
Definition: Ifpack2_Hiptmair_decl.hpp:140
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:423
double getApplyTime() const
Returns the time spent in apply().
Definition: Ifpack2_Hiptmair_def.hpp:204
void setParameters(const Teuchos::ParameterList ¶ms)
Set the preconditioner's parameters.
Definition: Ifpack2_Hiptmair_def.hpp:84
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:136
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:155
MatrixType::global_ordinal_type global_ordinal_type
The type of global indices in the input MatrixType.
Definition: Ifpack2_Hiptmair_decl.hpp:87
int getNumApply() const
Returns the number of calls to apply().
Definition: Ifpack2_Hiptmair_def.hpp:186
Wrapper for Hiptmair smoothers.
Definition: Ifpack2_Hiptmair_decl.hpp:70
Hiptmair(const Teuchos::RCP< const row_matrix_type > &A, const Teuchos::RCP< const row_matrix_type > &PtAP, const Teuchos::RCP< const row_matrix_type > &P)
Constructor that takes 3 Tpetra matrices.
Definition: Ifpack2_Hiptmair_def.hpp:57
int getNumCompute() const
Returns the number of calls to Compute().
Definition: Ifpack2_Hiptmair_def.hpp:180
int getNumInitialize() const
Returns the number of calls to Initialize().
Definition: Ifpack2_Hiptmair_def.hpp:174
MatrixType::node_type node_type
The Node type used by the input MatrixType.
Definition: Ifpack2_Hiptmair_decl.hpp:90
bool isComputed() const
Return true if compute() completed successfully, else false.
Definition: Ifpack2_Hiptmair_decl.hpp:148
void compute()
Do any initialization that depends on the input matrix's values.
Definition: Ifpack2_Hiptmair_def.hpp:247
virtual ~Hiptmair()
Destructor.
Definition: Ifpack2_Hiptmair_def.hpp:81
std::string description() const
Return a simple one-line description of this object.
Definition: Ifpack2_Hiptmair_def.hpp:393
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:143
double getInitializeTime() const
Returns the time spent in Initialize().
Definition: Ifpack2_Hiptmair_def.hpp:192
"Factory" for creating single-level preconditioners.
Definition: Ifpack2_Details_OneLevelFactory_decl.hpp:110
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:72
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:273