* ParameterList) and is similar to the deal.II class SolverSelector.
*
* @ingroup TrilinosWrappers
- * @author Martin Kronbichler, 2008, 2009
+ * @author Martin Kronbichler, 2008, 2009; extension for full compatibility
+ * with LinearOperator class: Jean-Paul Pelteret, 2015
*/
class SolverBase
{
* Trilinos is chosen.
*/
void
- solve (Epetra_Operator &A,
+ solve (const Epetra_Operator &A,
VectorBase &x,
const VectorBase &b,
const PreconditionBase &preconditioner);
+ /**
+ * Solve the linear system <tt>Ax=b</tt> where both <tt>A</tt> and its
+ * @p precondtioner are an operator.
+ * This function can be used when both <tt>A</tt> and the @p preconditioner
+ * are LinearOperators derived from a TrilinosPayload.
+ * Depending on the information provided by derived classes and the object
+ * passed as a preconditioner, one of the linear solvers and preconditioners
+ * of Trilinos is chosen.
+ */
+ void
+ solve (const Epetra_Operator &A,
+ VectorBase &x,
+ const VectorBase &b,
+ const Epetra_Operator &preconditioner);
+
+ /**
+ * Solve the linear system <tt>Ax=b</tt> where <tt>A</tt> is an operator,
+ * and the vectors @p x and @p b are native Trilinos vector types.
+ * This function can be used when <tt>A</tt> is a LinearOperators derived
+ * from a TrilinosPayload.
+ * Depending on the information provided by derived classes and the object
+ * passed as a preconditioner, one of the linear solvers and preconditioners
+ * of Trilinos is chosen.
+ */
+ void
+ solve (const Epetra_Operator &A,
+ Epetra_MultiVector &x,
+ const Epetra_MultiVector &b,
+ const PreconditionBase &preconditioner);
+
+ /**
+ * Solve the linear system <tt>Ax=b</tt> where both <tt>A</tt> and its
+ * @p precondtioner are an operator, and the vectors @p x and @p b are
+ * native Trilinos vector types.
+ * This function can be used when both <tt>A</tt> and the @p preconditioner
+ * are LinearOperators derived from a TrilinosPayload.
+ * Depending on the information provided by derived classes and the object
+ * passed as a preconditioner, one of the linear solvers and preconditioners
+ * of Trilinos is chosen.
+ */
+ void
+ solve (const Epetra_Operator &A,
+ Epetra_MultiVector &x,
+ const Epetra_MultiVector &b,
+ const Epetra_Operator &preconditioner);
+
+
+
/**
* Solve the linear system <tt>Ax=b</tt>. Depending on the information
* provided by derived classes and the object passed as a preconditioner,
* The solve function is used to set properly the Epetra_LinearProblem,
* once it is done this function solves the linear problem.
*/
- void do_solve(const PreconditionBase &preconditioner);
+ template<typename Preconditioner>
+ void do_solve(const Preconditioner &preconditioner);
+
+ /**
+ * A function that sets the preconditioner that the solver will apply
+ */
+ template<typename Preconditioner>
+ void set_preconditioner (AztecOO &solver,
+ const Preconditioner &preconditioner);
/**
* A structure that collects the Trilinos sparse matrix, the right hand
+ // Note: "A" is set as a constant reference so that all patterns for ::solve
+ // can be used by the inverse_operator of LinearOperator
void
- SolverBase::solve (Epetra_Operator &A,
+ SolverBase::solve (const Epetra_Operator &A,
VectorBase &x,
const VectorBase &b,
const PreconditionBase &preconditioner)
// We need an Epetra_LinearProblem object to let the AztecOO solver know
// about the matrix and vectors.
linear_problem.reset
- (new Epetra_LinearProblem(&A,
+ (new Epetra_LinearProblem(const_cast<Epetra_Operator *>(&A),
&x.trilinos_vector(),
const_cast<Epetra_MultiVector *>(&b.trilinos_vector())));
+ // Note: "A" is set as a constant reference so that all patterns for ::solve
+ // can be used by the inverse_operator of LinearOperator
+ void
+ SolverBase::solve (const Epetra_Operator &A,
+ VectorBase &x,
+ const VectorBase &b,
+ const Epetra_Operator &preconditioner)
+ {
+ linear_problem.reset();
+
+ // We need an Epetra_LinearProblem object to let the AztecOO solver know
+ // about the matrix and vectors.
+ linear_problem.reset
+ (new Epetra_LinearProblem(const_cast<Epetra_Operator *>(&A),
+ &x.trilinos_vector(),
+ const_cast<Epetra_MultiVector *>(&b.trilinos_vector())));
+
+ do_solve(preconditioner);
+ }
+
+
+
+ // Note: "A" is set as a constant reference so that all patterns for ::solve
+ // can be used by the inverse_operator of LinearOperator
+ void
+ SolverBase::solve (const Epetra_Operator &A,
+ Epetra_MultiVector &x,
+ const Epetra_MultiVector &b,
+ const PreconditionBase &preconditioner)
+ {
+ linear_problem.reset();
+
+ // We need an Epetra_LinearProblem object to let the AztecOO solver know
+ // about the matrix and vectors.
+ linear_problem.reset
+ (new Epetra_LinearProblem(const_cast<Epetra_Operator *>(&A),
+ &x,
+ const_cast<Epetra_MultiVector *>(&b)));
+
+ do_solve(preconditioner);
+ }
+
+
+
+ // Note: "A" is set as a constant reference so that all patterns for ::solve
+ // can be used by the inverse_operator of LinearOperator
+ void
+ SolverBase::solve (const Epetra_Operator &A,
+ Epetra_MultiVector &x,
+ const Epetra_MultiVector &b,
+ const Epetra_Operator &preconditioner)
+ {
+ linear_problem.reset();
+
+ // We need an Epetra_LinearProblem object to let the AztecOO solver know
+ // about the matrix and vectors.
+ linear_problem.reset
+ (new Epetra_LinearProblem(const_cast<Epetra_Operator *>(&A),
+ &x,
+ const_cast<Epetra_MultiVector *>(&b)));
+
+ do_solve(preconditioner);
+ }
+
+
+
void
SolverBase::solve (const SparseMatrix &A,
dealii::Vector<double> &x,
}
-
+ template<typename Preconditioner>
void
- SolverBase::do_solve(const PreconditionBase &preconditioner)
+ SolverBase::do_solve(const Preconditioner &preconditioner)
{
int ierr;
Assert (false, ExcNotImplemented());
}
- // Introduce the preconditioner, if the identity preconditioner is used,
- // the precondioner is set to none, ...
- if (preconditioner.preconditioner.use_count()!=0)
- {
- ierr = solver.SetPrecOperator (const_cast<Epetra_Operator *>
- (preconditioner.preconditioner.get()));
- AssertThrow (ierr == 0, ExcTrilinosError(ierr));
- }
- else
- solver.SetAztecOption(AZ_precond,AZ_none);
+ // Set the preconditioner
+ set_preconditioner(solver, preconditioner);
// ... set some options, ...
solver.SetAztecOption (AZ_output, additional_data.output_solver_details ?
+ template<>
+ void
+ SolverBase::set_preconditioner(AztecOO &solver,
+ const PreconditionBase &preconditioner)
+ {
+ // Introduce the preconditioner, if the identity preconditioner is used,
+ // the precondioner is set to none, ...
+ if (preconditioner.preconditioner.use_count()!=0)
+ {
+ const int ierr = solver.SetPrecOperator (const_cast<Epetra_Operator *>
+ (preconditioner.preconditioner.get()));
+ AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+ }
+ else
+ solver.SetAztecOption(AZ_precond,AZ_none);
+ }
+
+
+ template<>
+ void
+ SolverBase::set_preconditioner(AztecOO &solver,
+ const Epetra_Operator &preconditioner)
+ {
+ const int ierr = solver.SetPrecOperator (const_cast<Epetra_Operator *>(&preconditioner));
+ AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+ }
/* ---------------------- SolverCG ------------------------ */
}
+
+// explicit instantiations
+// TODO: put these instantiations into generic file
+namespace TrilinosWrappers
+{
+ template void
+ SolverBase::do_solve(const PreconditionBase &preconditioner);
+
+ template void
+ SolverBase::do_solve(const Epetra_Operator &preconditioner);
+
+ template void
+ SolverBase::set_preconditioner(AztecOO &solver,
+ const PreconditionBase &preconditioner);
+
+ template void
+ SolverBase::set_preconditioner(AztecOO &solver,
+ const Epetra_Operator &preconditioner);
+}
+
DEAL_II_NAMESPACE_CLOSE
#endif // DEAL_II_WITH_PETSC