43 #ifndef IFPACK2_DETAILS_AMESOS2WRAPPER_DEF_HPP 44 #define IFPACK2_DETAILS_AMESOS2WRAPPER_DEF_HPP 46 #ifdef HAVE_IFPACK2_AMESOS2 48 #include "Ifpack2_LocalFilter.hpp" 49 #include "Trilinos_Details_LinearSolverFactory.hpp" 50 #include "Trilinos_Details_LinearSolver.hpp" 51 #include "Teuchos_TimeMonitor.hpp" 52 #include "Teuchos_TypeNameTraits.hpp" 66 template <
class MatrixType>
70 InitializeTime_ (0.0),
76 IsInitialized_ (false),
81 template <
class MatrixType>
85 template <
class MatrixType>
88 using Teuchos::ParameterList;
97 RCP<ParameterList> theList;
98 if (params.name () ==
"Amesos2") {
99 theList = rcp (
new ParameterList (params));
100 }
else if (params.isSublist (
"Amesos2")) {
102 ParameterList subpl = params.sublist (
"Amesos2");
103 theList = rcp (
new ParameterList (subpl));
104 theList->setName (
"Amesos2");
105 if (params.isParameter (
"Amesos2 solver name")) {
106 SolverName_ = params.get<std::string>(
"Amesos2 solver name");
111 TEUCHOS_TEST_FOR_EXCEPTION(
112 true, std::runtime_error,
"The ParameterList passed to Amesos2 must be " 113 "called \"Amesos2\".");
118 if (solver_.is_null ()) {
119 parameterList_ = theList;
125 solver_->setParameters(theList);
129 template <
class MatrixType>
130 Teuchos::RCP<const Teuchos::Comm<int> >
132 TEUCHOS_TEST_FOR_EXCEPTION(
133 A_.is_null (), std::runtime_error,
"Ifpack2::Amesos2Wrapper::getComm: " 134 "The matrix is null. Please call setMatrix() with a nonnull input " 135 "before calling this method.");
136 return A_->getComm ();
140 template <
class MatrixType>
141 Teuchos::RCP<const typename Amesos2Wrapper<MatrixType>::row_matrix_type>
147 template <
class MatrixType>
148 Teuchos::RCP<const typename Amesos2Wrapper<MatrixType>::map_type>
151 TEUCHOS_TEST_FOR_EXCEPTION(
152 A_.is_null (), std::runtime_error,
"Ifpack2::Amesos2Wrapper::getDomainMap: " 153 "The matrix is null. Please call setMatrix() with a nonnull input " 154 "before calling this method.");
155 return A_->getDomainMap ();
159 template <
class MatrixType>
160 Teuchos::RCP<const typename Amesos2Wrapper<MatrixType>::map_type>
163 TEUCHOS_TEST_FOR_EXCEPTION(
164 A_.is_null (), std::runtime_error,
"Ifpack2::Amesos2Wrapper::getRangeMap: " 165 "The matrix is null. Please call setMatrix() with a nonnull input " 166 "before calling this method.");
167 return A_->getRangeMap ();
171 template <
class MatrixType>
177 template <
class MatrixType>
179 return NumInitialize_;
183 template <
class MatrixType>
189 template <
class MatrixType>
195 template <
class MatrixType>
197 return InitializeTime_;
201 template<
class MatrixType>
207 template<
class MatrixType>
212 template<
class MatrixType>
219 IsInitialized_ =
false;
240 template<
class MatrixType>
241 Teuchos::RCP<const typename Amesos2Wrapper<MatrixType>::row_matrix_type>
246 using Teuchos::rcp_dynamic_cast;
247 using Teuchos::rcp_implicit_cast;
252 if (A->getRowMap ()->getComm ()->getSize () == 1 ||
253 A->getRowMap ()->isSameAs (* (A->getColMap ()))) {
260 RCP<const LocalFilter<row_matrix_type> > A_lf_r =
262 if (! A_lf_r.is_null ()) {
263 return rcp_implicit_cast<
const row_matrix_type> (A_lf_r);
269 return rcp (
new LocalFilter<row_matrix_type> (A));
274 template<
class MatrixType>
279 using Teuchos::rcp_const_cast;
280 using Teuchos::rcp_dynamic_cast;
282 using Teuchos::TimeMonitor;
283 using Teuchos::Array;
284 using Teuchos::ArrayView;
286 const std::string timerName (
"Ifpack2::Amesos2Wrapper::initialize");
287 RCP<Time> timer = TimeMonitor::lookupCounter (timerName);
288 if (timer.is_null ()) {
289 timer = TimeMonitor::getNewCounter (timerName);
293 TimeMonitor timeMon (*timer);
296 TEUCHOS_TEST_FOR_EXCEPTION(
297 A_.is_null (), std::runtime_error,
"Ifpack2::Amesos2Wrapper::initialize: " 298 "The matrix to precondition is null. Please call setMatrix() with a " 299 "nonnull input before calling this method.");
302 IsInitialized_ =
false;
305 RCP<const row_matrix_type> A_local = makeLocalFilter (A_);
306 TEUCHOS_TEST_FOR_EXCEPTION(
307 A_local.is_null (), std::logic_error,
"Ifpack2::AmesosWrapper::initialize: " 308 "makeLocalFilter returned null; it failed to compute A_local. " 309 "Please report this bug to the Ifpack2 developers.");
316 if (A_local_crs_.is_null ()) {
318 Array<size_t> entriesPerRow(numRows);
321 entriesPerRow[i] = A_local->getNumEntriesInLocalRow(i);
323 RCP<crs_matrix_type> A_local_crs_nc =
325 A_local->getColMap (),
328 Teuchos::Array<local_ordinal_type> indices(A_local->getNodeMaxNumRowEntries());
329 Teuchos::Array<scalar_type> values(A_local->getNodeMaxNumRowEntries());
332 size_t numEntries = 0;
333 A_local->getLocalRowCopy(i, indices(), values(), numEntries);
334 ArrayView<const local_ordinal_type> indicesInsert(indices.data(), numEntries);
335 ArrayView<const scalar_type> valuesInsert(values.data(), numEntries);
336 A_local_crs_nc->insertLocalValues(i, indicesInsert, valuesInsert);
338 A_local_crs_nc->fillComplete (A_local->getDomainMap (), A_local->getRangeMap ());
351 if (! Trilinos::Details::Impl::rememberRegisteredSomeLinearSolverFactory (
"Amesos2")) {
352 Amesos2::Details::registerLinearSolverFactory ();
355 solver_ = Trilinos::Details::getLinearSolver<MV, OP, typename MV::mag_type> (
"Amesos2", SolverName_);
356 TEUCHOS_TEST_FOR_EXCEPTION
357 (solver_.is_null (), std::runtime_error,
"Ifpack2::Details::" 358 "Amesos2Wrapper::initialize: Failed to create Amesos2 solver!");
360 solver_->setMatrix (A_local_crs_);
362 if (parameterList_ != Teuchos::null) {
363 setParameters (*parameterList_);
364 parameterList_ = Teuchos::null;
369 solver_->symbolic ();
372 IsInitialized_ =
true;
377 InitializeTime_ = timer->totalElapsedTime ();
380 template<
class MatrixType>
385 using Teuchos::TimeMonitor;
388 if (! isInitialized ()) {
392 const std::string timerName (
"Ifpack2::Details::Amesos2Wrapper::compute");
393 RCP<Time> timer = TimeMonitor::lookupCounter (timerName);
394 if (timer.is_null ()) {
395 timer = TimeMonitor::getNewCounter (timerName);
399 TimeMonitor timeMon (*timer);
408 ComputeTime_ = timer->totalElapsedTime ();
412 template <
class MatrixType>
414 apply (
const Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& X,
415 Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type>& Y,
416 Teuchos::ETransp mode,
420 using Teuchos::ArrayView;
423 using Teuchos::rcpFromRef;
425 using Teuchos::TimeMonitor;
430 const std::string timerName (
"Ifpack2::Amesos2Wrapper::apply");
431 RCP<Time> timer = TimeMonitor::lookupCounter (timerName);
432 if (timer.is_null ()) {
433 timer = TimeMonitor::getNewCounter (timerName);
437 TimeMonitor timeMon (*timer);
439 TEUCHOS_TEST_FOR_EXCEPTION(
440 ! isComputed (), std::runtime_error,
441 "Ifpack2::Amesos2Wrapper::apply: You must call compute() to compute the " 442 "incomplete factorization, before calling apply().");
444 TEUCHOS_TEST_FOR_EXCEPTION(
445 X.getNumVectors() != Y.getNumVectors(), std::runtime_error,
446 "Ifpack2::Amesos2Wrapper::apply: X and Y must have the same number of columns. " 447 "X has " << X.getNumVectors () <<
" columns, but Y has " 448 << Y.getNumVectors () <<
" columns.");
450 TEUCHOS_TEST_FOR_EXCEPTION(
451 mode != Teuchos::NO_TRANS, std::logic_error,
452 "Ifpack2::Amesos2Wrapper::apply: Solving with the transpose (mode == " 453 "Teuchos::TRANS) or conjugate transpose (Teuchos::CONJ_TRANS) is not " 459 RCP<MV> Y_temp = (alpha != STS::one () || beta != STS::zero ()) ?
460 rcp (
new MV (Y.getMap (), Y.getNumVectors ())) :
466 RCP<const MV> X_temp;
468 auto X_lcl_host = X.getLocalViewHost ();
469 auto Y_lcl_host = Y.getLocalViewHost ();
471 if (X_lcl_host.data () == Y_lcl_host.data ()) {
472 X_temp = rcp (
new MV (X, Teuchos::Copy));
474 X_temp = rcpFromRef (X);
479 RCP<const MV> X_local;
485 const bool multipleProcs = (A_->getRowMap ()->getComm ()->getSize () > 1) || (X.getMap ()->getComm ()->getSize () > 1);
490 X_local = X_temp->offsetView (A_local_crs_->getDomainMap (), 0);
491 Y_local = Y_temp->offsetViewNonConst (A_local_crs_->getRangeMap (), 0);
500 solver_->solve (*Y_local, *X_local);
502 if (alpha != STS::one () || beta != STS::zero ()) {
503 Y.update (alpha, *Y_temp, beta);
511 ApplyTime_ = timer->totalElapsedTime ();
515 template <
class MatrixType>
517 using Teuchos::TypeNameTraits;
518 std::ostringstream os;
523 os <<
"\"Ifpack2::Amesos2Wrapper\": {";
524 if (this->getObjectLabel () !=
"") {
525 os <<
"Label: \"" << this->getObjectLabel () <<
"\", ";
527 os <<
"Initialized: " << (isInitialized () ?
"true" :
"false")
528 <<
", Computed: " << (isComputed () ?
"true" :
"false");
530 if (A_local_crs_.is_null ()) {
531 os <<
", Matrix: null";
534 os <<
", Global matrix dimensions: [" 535 << A_local_crs_->getGlobalNumRows () <<
", " << A_local_crs_->getGlobalNumCols () <<
"]";
540 if (! solver_.is_null ()) {
541 Teuchos::Describable* d =
dynamic_cast<Teuchos::Describable*
> (solver_.getRawPtr ());
544 os << d->description ();
553 template <
class MatrixType>
557 const Teuchos::EVerbosityLevel verbLevel)
const 560 using Teuchos::OSTab;
562 using Teuchos::TypeNameTraits;
565 const Teuchos::EVerbosityLevel vl = (verbLevel == Teuchos::VERB_DEFAULT) ?
566 Teuchos::VERB_LOW : verbLevel;
570 if (vl > Teuchos::VERB_NONE) {
571 out <<
"\"Ifpack2::Amesos2Wrapper\":" << endl;
573 out <<
"MatrixType: \"" << TypeNameTraits<MatrixType>::name ()
576 if (this->getObjectLabel () !=
"") {
577 out <<
"Label: \"" << this->getObjectLabel () <<
"\"" << endl;
580 out <<
"Initialized: " << (isInitialized () ?
"true" :
"false") << endl;
581 out <<
"Computed: " << (isComputed () ?
"true" :
"false") << endl;
582 out <<
"Number of initialize calls: " << getNumInitialize () << endl;
583 out <<
"Number of compute calls: " << getNumCompute () << endl;
584 out <<
"Number of apply calls: " << getNumApply () << endl;
585 out <<
"Total time in seconds for initialize: " << getInitializeTime () << endl;
586 out <<
"Total time in seconds for compute: " << getComputeTime () << endl;
587 out <<
"Total time in seconds for apply: " << getApplyTime () << endl;
589 if (vl > Teuchos::VERB_LOW) {
590 out <<
"Matrix:" << endl;
591 A_local_crs_->describe (out, vl);
603 #define IFPACK2_DETAILS_AMESOS2WRAPPER_INSTANT(S,LO,GO,N) \ 604 template class Ifpack2::Details::Amesos2Wrapper< Tpetra::RowMatrix<S, LO, GO, N> >; 608 #define IFPACK2_DETAILS_AMESOS2WRAPPER_INSTANT(S,LO,GO,N) 610 #endif // HAVE_IFPACK2_AMESOS2 611 #endif // IFPACK2_DETAILS_AMESOS2WRAPPER_DEF_HPP double getApplyTime() const
The total time in seconds spent in successful calls to apply().
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:208
void setParameters(const Teuchos::ParameterList ¶ms)
Set parameters.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:86
bool hasTransposeApply() const
Whether this object's apply() method can apply the transpose (or conjugate transpose, if applicable).
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:172
Teuchos::RCP< const row_matrix_type > getMatrix() const
The input matrix; the matrix to be preconditioned.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:142
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, resulting in Y.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:414
void initialize()
Compute the preordering and symbolic factorization of the matrix.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:275
Ifpack2 implementation details.
Teuchos::RCP< const map_type > getDomainMap() const
Tpetra::Map representing the domain of this operator.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:149
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to the given output stream.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:556
std::string description() const
A one-line description of this object.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:516
MatrixType::scalar_type scalar_type
The type of the entries of the input MatrixType.
Definition: Ifpack2_Details_Amesos2Wrapper_decl.hpp:117
Tpetra::CrsMatrix< scalar_type, local_ordinal_type, global_ordinal_type, node_type > crs_matrix_type
Type of the Tpetra::CrsMatrix specialization that this class uses.
Definition: Ifpack2_Details_Amesos2Wrapper_decl.hpp:149
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:57
Wrapper class for direct solvers in Amesos2.
Definition: Ifpack2_Details_Amesos2Wrapper_decl.hpp:102
int getNumCompute() const
The total number of successful calls to compute().
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:184
virtual ~Amesos2Wrapper()
Destructor.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:82
double getInitializeTime() const
The total time in seconds spent in successful calls to initialize().
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:196
Teuchos::RCP< const map_type > getRangeMap() const
Tpetra::Map representing the range of this operator.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:161
MatrixType::local_ordinal_type local_ordinal_type
The type of local indices in the input MatrixType.
Definition: Ifpack2_Details_Amesos2Wrapper_decl.hpp:120
Amesos2Wrapper(const Teuchos::RCP< const row_matrix_type > &A)
Constructor.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:68
double getComputeTime() const
The total time in seconds spent in successful calls to compute().
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:202
int getNumApply() const
The total number of successful calls to apply().
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:190
void compute()
Compute the numeric factorization of the matrix.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:381
Access only local rows and columns of a sparse matrix.
Definition: Ifpack2_LocalFilter_decl.hpp:160
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:73
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
The input matrix's communicator.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:131
void registerLinearSolverFactory()
int getNumInitialize() const
The total number of successful calls to initialize().
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:178
virtual void setMatrix(const Teuchos::RCP< const row_matrix_type > &A)
Change the matrix to be preconditioned.
Definition: Ifpack2_Details_Amesos2Wrapper_def.hpp:213