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;
55 template <
typename number>
class Vector;
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 202 Epetra_SerialComm communicator;
269 const double min_diagonal = 0);
384 const double min_diagonal = 0,
385 const unsigned int overlap = 0);
512 const double min_diagonal = 0,
513 const unsigned int overlap = 0);
664 const double ic_atol = 0.,
665 const double ic_rtol = 1.,
666 const unsigned int overlap = 0);
832 const double ilu_atol = 0.,
833 const double ilu_rtol = 1.,
834 const unsigned int overlap = 0);
997 const unsigned int ilut_fill = 0,
998 const double ilut_atol = 0.,
999 const double ilut_rtol = 1.,
1000 const unsigned int overlap = 0);
1159 const double max_eigenvalue = 10.,
1160 const double eigenvalue_ratio = 30.,
1161 const double min_eigenvalue = 1.,
1162 const double min_diagonal = 1e-12,
1163 const bool nonzero_starting =
false);
1302 const bool higher_order_elements =
false,
1303 const unsigned int n_cycles = 1,
1304 const bool w_cyle =
false,
1305 const double aggregation_threshold = 1e-4,
1306 const std::vector<std::vector<bool> > &constant_modes = std::vector<std::vector<bool> > (1),
1307 const unsigned int smoother_sweeps = 2,
1308 const unsigned int smoother_overlap = 0,
1309 const bool output_details =
false);
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"));
1595 ExcNonMatchingMaps(
"dst"));
1597 ExcNonMatchingMaps(
"src"));
1599 preconditioner->SetUseTranspose(
true);
1603 preconditioner->SetUseTranspose(
false);
1633 const ::Vector<double> &src)
const 1636 preconditioner->OperatorDomainMap().NumMyElements());
1637 AssertDimension (static_cast<TrilinosWrappers::types::int_type>(src.size()),
1638 preconditioner->OperatorRangeMap().NumMyElements());
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 1654 preconditioner->OperatorDomainMap().NumMyElements());
1655 AssertDimension (static_cast<TrilinosWrappers::types::int_type>(src.size()),
1656 preconditioner->OperatorRangeMap().NumMyElements());
1657 Epetra_Vector tril_dst (View, preconditioner->OperatorDomainMap(),
1659 Epetra_Vector tril_src (View, preconditioner->OperatorRangeMap(),
1660 const_cast<double *
>(src.begin()));
1662 preconditioner->SetUseTranspose(
true);
1663 const int ierr = preconditioner->ApplyInverse (tril_src, tril_dst);
1665 preconditioner->SetUseTranspose(
false);
1676 preconditioner->OperatorDomainMap().NumMyElements());
1678 preconditioner->OperatorRangeMap().NumMyElements());
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);
1694 preconditioner->OperatorDomainMap().NumMyElements());
1696 preconditioner->OperatorRangeMap().NumMyElements());
1697 Epetra_Vector tril_dst (View, preconditioner->OperatorDomainMap(),
1699 Epetra_Vector tril_src (View, preconditioner->OperatorRangeMap(),
1700 const_cast<double *
>(src.
begin()));
1702 preconditioner->SetUseTranspose(
true);
1703 const int ierr = preconditioner->ApplyInverse (tril_src, tril_dst);
1705 preconditioner->SetUseTranspose(
false);
1716 DEAL_II_NAMESPACE_CLOSE
1718 #endif // DEAL_II_WITH_TRILINOS
#define AssertDimension(dim1, dim2)
double aggregation_threshold
std::vector< std::vector< bool > > constant_modes
unsigned int smoother_overlap
#define AssertThrow(cond, exc)
unsigned int smoother_sweeps
#define DeclException1(Exception1, type1, outsequence)
#define Assert(cond, exc)
std_cxx1x::shared_ptr< Epetra_Map > vector_distributor
const Epetra_Map & vector_partitioner() const
std_cxx1x::shared_ptr< SparseMatrix > trilinos_matrix
std_cxx1x::shared_ptr< Epetra_Operator > preconditioner
size_type local_size() const
const Epetra_MultiVector & trilinos_vector() const
::types::global_dof_index size_type
Epetra_MpiComm communicator
bool higher_order_elements