20 #include "libmesh/libmesh_common.h" 22 #ifdef LIBMESH_TRILINOS_HAVE_AZTECOO 26 #include "libmesh/libmesh_logging.h" 27 #include "libmesh/enum_to_string.h" 28 #include "libmesh/trilinos_aztec_linear_solver.h" 29 #include "libmesh/trilinos_epetra_matrix.h" 30 #include "libmesh/trilinos_epetra_vector.h" 31 #include "libmesh/enum_preconditioner_type.h" 32 #include "libmesh/enum_solver_type.h" 33 #include "libmesh/enum_convergence_flags.h" 60 this->_solver_type =
GMRES;
62 if (this->n_processors() == 1)
79 _linear_solver =
new AztecOO();
83 switch(this->_preconditioner_type)
86 _linear_solver->SetAztecOption(AZ_precond,AZ_dom_decomp);
87 _linear_solver->SetAztecOption(AZ_subdomain_solve,AZ_ilu);
91 _linear_solver->SetAztecOption(AZ_precond,AZ_Jacobi);
95 _linear_solver->SetAztecOption(AZ_precond,AZ_dom_decomp);
96 _linear_solver->SetAztecOption(AZ_subdomain_solve,AZ_icc);
100 _linear_solver->SetAztecOption(AZ_precond,AZ_dom_decomp);
101 _linear_solver->SetAztecOption(AZ_subdomain_solve,AZ_lu);
105 _linear_solver->SetAztecOption(AZ_precond,AZ_dom_decomp);
106 _linear_solver->SetAztecOption(AZ_subdomain_solve,AZ_ilu);
114 template <
typename T>
115 std::pair<unsigned int, Real>
120 const std::optional<double> tol,
121 const std::optional<unsigned int> m_its)
123 LOG_SCOPE(
"solve()",
"AztecLinearSolver");
128 EpetraVector<T> * solution = cast_ptr<EpetraVector<T> *>(&solution_in);
139 const int max_its = this->get_int_solver_setting(
"max_its", m_its);
140 const double rel_tol = this->get_real_solver_setting(
"rel_tol", tol);
142 _linear_solver->SetAztecOption(AZ_max_iter,max_its);
143 _linear_solver->SetAztecParam(AZ_tol,rel_tol);
145 Epetra_FECrsMatrix * emat = matrix->
mat();
146 Epetra_Vector * esol = solution->vec();
147 Epetra_Vector * erhs = rhs->vec();
149 _linear_solver->Iterate(emat, esol, erhs, max_its, rel_tol);
152 return std::make_pair(_linear_solver->NumIters(), _linear_solver->TrueResidual());
157 template <
typename T>
158 std::pair<unsigned int, Real>
162 const std::optional<double>,
163 const std::optional<unsigned int>)
165 libmesh_not_implemented();
170 template <
typename T>
171 std::pair<unsigned int, Real>
176 const std::optional<double>,
177 const std::optional<unsigned int>)
179 libmesh_not_implemented();
184 template <
typename T>
187 libmesh_not_implemented();
193 template <
typename T>
196 return _linear_solver->TrueResidual();
201 template <
typename T>
204 const double *
status = _linear_solver->GetAztecStatus();
206 switch (static_cast<int>(
status[AZ_why]))
212 libMesh::out <<
"AztecOO failed to converge within maximum iterations.\n";
215 libMesh::out <<
"AztecOO failed to support a user-requested parameter.\n";
218 libMesh::out <<
"AztecOO encountered numerical breakdown.\n";
221 libMesh::out <<
"AztecOO encountered numerical loss of precision.\n";
224 libMesh::out <<
"AztecOO encountered an ill-conditioned GMRES Hessian.\n";
227 libMesh::out <<
"AztecOO reported an unrecognized condition.\n";
234 template <
typename T>
237 const double *
status = _linear_solver->GetAztecStatus();
239 switch (static_cast<int>(
status[AZ_why]))
253 template <
typename T>
256 switch (this->_solver_type)
259 _linear_solver->SetAztecOption(AZ_solver, AZ_cg);
return;
262 _linear_solver->SetAztecOption(AZ_solver, AZ_cgs);
return;
265 _linear_solver->SetAztecOption(AZ_solver, AZ_tfqmr);
return;
268 _linear_solver->SetAztecOption(AZ_solver, AZ_bicgstab);
return;
271 _linear_solver->SetAztecOption(AZ_solver, AZ_gmres);
return;
276 <<
"Continuing with AztecOO defaults" << std::endl;
288 #endif // #ifdef LIBMESH_TRILINOS_HAVE_AZTECOO
This class provides a nice interface to the Trilinos Epetra_Vector object.
void get_residual_history(std::vector< double > &hist)
Fills the input vector with the sequence of residual norms from the latest iterative solve...
LinearConvergenceReason
Linear solver convergence flags (taken from the PETSc flags).
Provides a uniform interface to vector storage schemes for different linear algebra libraries...
The libMesh namespace provides an interface to certain functionality in the library.
This base class can be inherited from to provide interfaces to linear solvers from different packages...
processor_id_type n_processors() const
bool _is_initialized
Flag that tells if init() has been called.
virtual void clear() override
Release all memory and clear data structures.
void init(triangulateio &t)
Initializes the fields of t to nullptr/0 as necessary.
AztecLinearSolver(const libMesh::Parallel::Communicator &comm)
Constructor.
virtual std::pair< unsigned int, Real > solve(SparseMatrix< T > &matrix_in, NumericVector< T > &solution_in, NumericVector< T > &rhs_in, const std::optional< double > tol=std::nullopt, const std::optional< unsigned int > m_its=std::nullopt) override
Call the Aztec solver.
virtual void init(const char *name=nullptr) override
Initialize data structures if not done so already.
virtual void close() override
Calls the SparseMatrix's internal assembly routines, ensuring that the values are consistent across p...
std::string enum_to_string(const T e)
This class provides an interface to AztecOO iterative solvers that is compatible with the libMesh Lin...
void set_solver_type()
Tells AztecOO to use the user-specified solver stored in _solver_type.
virtual void print_converged_reason() const override
Prints a useful message about why the latest linear solve con(di)verged.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
PreconditionerType _preconditioner_type
Enum stating with type of preconditioner to use.
virtual LinearConvergenceReason get_converged_reason() const override
bool initialized()
Checks that library initialization has been done.
Epetra_FECrsMatrix * mat()
Generic shell matrix, i.e.
This class provides a nice interface to the Epetra data structures for parallel, sparse matrices...
Real get_initial_residual()