20 #ifndef LIBMESH_TRILINOS_PRECONDITIONER_H 21 #define LIBMESH_TRILINOS_PRECONDITIONER_H 23 #include "libmesh/libmesh_config.h" 25 #ifdef LIBMESH_TRILINOS_HAVE_EPETRA 28 #include "libmesh/preconditioner.h" 29 #include "libmesh/libmesh_common.h" 30 #include "libmesh/reference_counted_object.h" 31 #include "libmesh/libmesh.h" 35 #include "libmesh/ignore_warnings.h" 36 #include "Epetra_Operator.h" 37 #include "Epetra_FECrsMatrix.h" 38 #include "Teuchos_ParameterList.hpp" 39 #include "libmesh/restore_warnings.h" 48 template <
typename T>
class SparseMatrix;
49 template <
typename T>
class NumericVector;
50 template <
typename T>
class ShellMatrix;
64 public Epetra_Operator
80 virtual void clear ()
override {}
82 virtual void init ()
override;
87 void set_params(Teuchos::ParameterList & list);
92 Epetra_FECrsMatrix *
mat() {
return _mat; }
125 virtual int Apply(
const Epetra_MultiVector & X, Epetra_MultiVector & Y)
const override;
126 virtual int ApplyInverse(
const Epetra_MultiVector & r, Epetra_MultiVector & z)
const override;
127 virtual double NormInf()
const override;
128 virtual const char *
Label()
const override;
131 virtual const Epetra_Comm &
Comm()
const override;
140 template <
typename T>
151 template <
typename T>
160 #endif // LIBMESH_TRILINOS_HAVE_EPETRA 161 #endif // LIBMESH_TRILINOS_PRECONDITIONER_H void set_params(Teuchos::ParameterList &list)
Stores a copy of the ParameterList list internally.
virtual void init() override
Initialize data structures if not done so already.
virtual ~TrilinosPreconditioner()
Destructor.
TrilinosPreconditioner(const libMesh::Parallel::Communicator &comm)
Constructor.
Teuchos::ParameterList _param_list
Parameter list to be used for building the preconditioner.
virtual bool HasNormInf() const override
Epetra_FECrsMatrix * mat()
virtual void apply(const NumericVector< T > &x, NumericVector< T > &y) override
Computes the preconditioned vector y based on input vector x.
Provides a uniform interface to vector storage schemes for different linear algebra libraries...
const Parallel::Communicator & comm() const
The libMesh namespace provides an interface to certain functionality in the library.
void set_preconditioner_type(const PreconditionerType &preconditioner_type)
Sets the Trilinos preconditioner to use based on the libMesh PreconditionerType.
virtual void clear() override
Release all memory and clear data structures.
This class provides an interface to the suite of preconditioners available from Trilinos.
This class provides a uniform interface for preconditioners.
Epetra_FECrsMatrix * _mat
Trilinos matrix that's been pulled out of the _matrix object.
void compute()
Compute the preconditioner.
virtual const Epetra_Map & OperatorDomainMap() const override
virtual bool UseTranspose() const override
virtual double NormInf() const override
PreconditionerType
Defines an enum for preconditioner types.
virtual int SetUseTranspose(bool UseTranspose) override
virtual const char * Label() const override
virtual const Epetra_Map & OperatorRangeMap() const override
virtual int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const override
Epetra_Operator * _prec
Trilinos preconditioner.
virtual const Epetra_Comm & Comm() const override
virtual int ApplyInverse(const Epetra_MultiVector &r, Epetra_MultiVector &z) const override
void ErrorVector unsigned int