17 #ifndef __deal2__trilinos_precondition_h
18 #define __deal2__trilinos_precondition_h
21 #include <deal.II/base/config.h>
23 #ifdef DEAL_II_WITH_TRILINOS
25 # include <deal.II/base/subscriptor.h>
26 # include <deal.II/base/std_cxx1x/shared_ptr.h>
28 # include <deal.II/lac/trilinos_vector_base.h>
29 # include <deal.II/lac/parallel_vector.h>
31 # ifdef DEAL_II_WITH_MPI
32 # include <Epetra_MpiComm.h>
34 # include <Epetra_SerialComm.h>
36 # include <Epetra_Map.h>
38 # include <Teuchos_ParameterList.hpp>
39 # include <Epetra_Operator.h>
40 # include <Epetra_Vector.h>
43 class Ifpack_Preconditioner;
44 class Ifpack_Chebyshev;
47 class MultiLevelPreconditioner;
51 DEAL_II_NAMESPACE_OPEN
55 template <
typename number>
class Vector;
62 namespace TrilinosWrappers
141 const ::Vector<double> &src)
const;
151 const ::Vector<double> &src)
const;
160 const ::parallel::distributed::Vector<double> &src)
const;
169 const ::parallel::distributed::Vector<double> &src)
const;
176 <<
"The sparse matrix the preconditioner is based on "
177 <<
"uses a map that is not compatible to the one in vector "
179 <<
". Check preconditioner and matrix setup.");
182 friend class PreconditionStokes;
199 #ifdef DEAL_II_WITH_MPI
385 const unsigned int overlap = 0);
513 const unsigned int overlap = 0);
666 const unsigned int overlap = 0);
834 const unsigned int overlap = 0);
1000 const unsigned int overlap = 0);
1304 const bool w_cyle =
false,
1306 const std::vector<std::vector<bool> > &
constant_modes = std::vector<std::vector<bool> > (1),
1445 const Teuchos::ParameterList &ml_parameters);
1457 template <
typename number>
1458 void initialize (const ::SparseMatrix<number> &deal_ii_sparse_matrix,
1460 const double drop_tolerance = 1e-13,
1461 const ::SparsityPattern *use_this_sparsity = 0);
1540 const ::Vector<double> &src)
const;
1548 const ::Vector<double> &src)
const;
1556 const ::parallel::distributed::Vector<double> &src)
const;
1564 const ::parallel::distributed::Vector<double> &src)
const;
1580 ExcNonMatchingMaps(
"dst"));
1582 ExcNonMatchingMaps(
"src"));
1592 const VectorBase &src)
const
1595 ExcNonMatchingMaps(
"dst"));
1597 ExcNonMatchingMaps(
"src"));
1600 const int ierr =
preconditioner->ApplyInverse (src.trilinos_vector(),
1601 dst.trilinos_vector());
1633 const ::Vector<double> &src)
const
1637 AssertDimension (static_cast<TrilinosWrappers::types::int_type>(src.size()),
1639 Epetra_Vector tril_dst (View,
preconditioner->OperatorDomainMap(),
1641 Epetra_Vector tril_src (View,
preconditioner->OperatorRangeMap(),
1642 const_cast<double *
>(src.begin()));
1644 const int ierr =
preconditioner->ApplyInverse (tril_src, tril_dst);
1651 const ::Vector<double> &src)
const
1655 AssertDimension (static_cast<TrilinosWrappers::types::int_type>(src.size()),
1657 Epetra_Vector tril_dst (View,
preconditioner->OperatorDomainMap(),
1659 Epetra_Vector tril_src (View,
preconditioner->OperatorRangeMap(),
1660 const_cast<double *
>(src.begin()));
1663 const int ierr =
preconditioner->ApplyInverse (tril_src, tril_dst);
1679 Epetra_Vector tril_dst (View,
preconditioner->OperatorDomainMap(),
1681 Epetra_Vector tril_src (View,
preconditioner->OperatorRangeMap(),
1682 const_cast<double *
>(src.
begin()));
1684 const int ierr =
preconditioner->ApplyInverse (tril_src, tril_dst);
1697 Epetra_Vector tril_dst (View,
preconditioner->OperatorDomainMap(),
1699 Epetra_Vector tril_src (View,
preconditioner->OperatorRangeMap(),
1700 const_cast<double *
>(src.
begin()));
1703 const int ierr =
preconditioner->ApplyInverse (tril_src, tril_dst);
1716 DEAL_II_NAMESPACE_CLOSE
1718 #endif // DEAL_II_WITH_TRILINOS
AdditionalData(const unsigned int degree=1, const double max_eigenvalue=10., const double eigenvalue_ratio=30., const double min_eigenvalue=1., const double min_diagonal=1e-12, const bool nonzero_starting=false)
#define AssertDimension(dim1, dim2)
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
AdditionalData(const double omega=1, const double min_diagonal=0, const unsigned int overlap=0)
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
virtual void Tvmult(VectorBase &dst, const VectorBase &src) const
double aggregation_threshold
std::vector< std::vector< bool > > constant_modes
unsigned int smoother_overlap
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
#define AssertThrow(cond, exc)
AdditionalData(const unsigned int ilu_fill=0, const double ilu_atol=0., const double ilu_rtol=1., const unsigned int overlap=0)
AdditionalData(const double omega=1, const double min_diagonal=0)
unsigned int smoother_sweeps
AdditionalData(const unsigned int overlap=0)
AdditionalData(const bool elliptic=true, const bool higher_order_elements=false, const unsigned int n_cycles=1, const bool w_cyle=false, const double aggregation_threshold=1e-4, const std::vector< std::vector< bool > > &constant_modes=std::vector< std::vector< bool > >(1), const unsigned int smoother_sweeps=2, const unsigned int smoother_overlap=0, const bool output_details=false)
#define Assert(cond, exc)
DeclException1(ExcNonMatchingMaps, std::string,<< "The sparse matrix the preconditioner is based on "<< "uses a map that is not compatible to the one in vector "<< arg1<< ". Check preconditioner and matrix setup.")
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void Tvmult(VectorBase &dst, const VectorBase &src) const
std_cxx1x::shared_ptr< Epetra_Map > vector_distributor
size_type local_size() const
AdditionalData(const double ilut_drop=0., const unsigned int ilut_fill=0, const double ilut_atol=0., const double ilut_rtol=1., const unsigned int overlap=0)
const Epetra_Map & vector_partitioner() const
AdditionalData(const double omega=1, const double min_diagonal=0, const unsigned int overlap=0)
std_cxx1x::shared_ptr< SparseMatrix > trilinos_matrix
void vmult(VectorBase &dst, const VectorBase &src) const
std_cxx1x::shared_ptr< Epetra_Operator > preconditioner
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
AdditionalData(const unsigned int ic_fill=0, const double ic_atol=0., const double ic_rtol=1., const unsigned int overlap=0)
virtual void vmult(VectorBase &dst, const VectorBase &src) const
const Epetra_MultiVector & trilinos_vector() const
::types::global_dof_index size_type
Epetra_MpiComm communicator
bool higher_order_elements
size_type memory_consumption() const