16 #ifndef dealii_trilinos_precondition_h 17 # define dealii_trilinos_precondition_h 20 # include <deal.II/base/config.h> 22 # ifdef DEAL_II_WITH_TRILINOS 24 # include <deal.II/base/subscriptor.h> 26 # include <deal.II/lac/la_parallel_vector.h> 27 # include <deal.II/lac/trilinos_vector.h> 31 # ifdef DEAL_II_WITH_MPI 32 # include <Epetra_MpiComm.h> 34 # include <Epetra_SerialComm.h> 36 # include <Epetra_Map.h> 37 # include <Epetra_RowMatrix.h> 38 # include <Epetra_Vector.h> 39 # include <Teuchos_ParameterList.hpp> 42 class Ifpack_Preconditioner;
43 class Ifpack_Chebyshev;
46 class MultiLevelPreconditioner;
50 DEAL_II_NAMESPACE_OPEN
53 template <
typename number>
55 template <
typename number>
121 get_mpi_communicator()
const;
160 const ::Vector<double> &src)
const;
168 const ::LinearAlgebra::distributed::Vector<double> &src)
const;
176 const ::LinearAlgebra::distributed::Vector<double> &src)
const;
188 trilinos_operator()
const;
201 locally_owned_domain_indices()
const;
209 locally_owned_range_indices()
const;
223 <<
"The sparse matrix the preconditioner is based on " 224 <<
"uses a map that is not compatible to the one in vector " 225 << arg1 <<
". Check preconditioner and matrix setup.");
241 # ifdef DEAL_II_WITH_MPI 244 Epetra_SerialComm communicator;
292 const double min_diagonal = 0,
293 const unsigned int n_sweeps = 1);
378 const double min_diagonal = 0,
379 const unsigned int overlap = 0,
380 const unsigned int n_sweeps = 1);
473 const double min_diagonal = 0,
474 const unsigned int overlap = 0,
475 const unsigned int n_sweeps = 1);
562 const std::string &block_creation_type =
"linear",
563 const double omega = 1,
564 const double min_diagonal = 0,
565 const unsigned int n_sweeps = 1);
666 const std::string &block_creation_type =
"linear",
667 const double omega = 1,
668 const double min_diagonal = 0,
669 const unsigned int overlap = 0,
670 const unsigned int n_sweeps = 1);
778 const std::string &block_creation_type =
"linear",
779 const double omega = 1,
780 const double min_diagonal = 0,
781 const unsigned int overlap = 0,
782 const unsigned int n_sweeps = 1);
904 const double ic_atol = 0.,
905 const double ic_rtol = 1.,
906 const unsigned int overlap = 0);
1024 const double ilu_atol = 0.,
1025 const double ilu_rtol = 1.,
1026 const unsigned int overlap = 0);
1130 const unsigned int ilut_fill = 0,
1131 const double ilut_atol = 0.,
1132 const double ilut_rtol = 1.,
1133 const unsigned int overlap = 0);
1255 const double max_eigenvalue = 10.,
1256 const double eigenvalue_ratio = 30.,
1257 const double min_eigenvalue = 1.,
1258 const double min_diagonal = 1e-12,
1259 const bool nonzero_starting =
false);
1373 const bool higher_order_elements =
false,
1374 const unsigned int n_cycles = 1,
1375 const bool w_cyle =
false,
1376 const double aggregation_threshold = 1e-4,
1377 const std::vector<std::vector<bool>> &constant_modes =
1378 std::vector<std::vector<bool>>(0),
1379 const unsigned int smoother_sweeps = 2,
1380 const unsigned int smoother_overlap = 0,
1381 const bool output_details =
false,
1382 const char * smoother_type =
"Chebyshev",
1383 const char * coarse_type =
"Amesos-KLU");
1543 initialize(
const Epetra_RowMatrix &matrix,
1560 const Teuchos::ParameterList &ml_parameters);
1570 initialize(
const Epetra_RowMatrix & matrix,
1571 const Teuchos::ParameterList &ml_parameters);
1579 template <
typename number>
1581 initialize(const ::SparseMatrix<number> &deal_ii_sparse_matrix,
1583 const double drop_tolerance = 1e-13,
1584 const ::SparsityPattern *use_this_sparsity =
nullptr);
1612 memory_consumption()
const;
1623 # if defined(DOXYGEN) || DEAL_II_TRILINOS_VERSION_GTE(11, 14, 0) 1657 const unsigned int n_cycles = 1,
1658 const bool w_cyle =
false,
1659 const double aggregation_threshold = 1e-4,
1660 const std::vector<std::vector<bool>> &constant_modes =
1661 std::vector<std::vector<bool>>(0),
1662 const unsigned int smoother_sweeps = 2,
1663 const unsigned int smoother_overlap = 0,
1664 const bool output_details =
false,
1665 const char * smoother_type =
"Chebyshev",
1666 const char * coarse_type =
"Amesos-KLU");
1802 initialize(
const Epetra_CrsMatrix &matrix,
1818 Teuchos::ParameterList &muelu_parameters);
1826 initialize(
const Epetra_CrsMatrix &matrix,
1827 Teuchos::ParameterList &muelu_parameters);
1835 template <
typename number>
1837 initialize(const ::SparseMatrix<number> &deal_ii_sparse_matrix,
1839 const double drop_tolerance = 1e-13,
1840 const ::SparsityPattern *use_this_sparsity =
nullptr);
1853 memory_consumption()
const;
1911 const ::Vector<double> &src)
const override;
1919 const ::Vector<double> &src)
const override;
1927 const ::LinearAlgebra::distributed::Vector<double> &src)
1937 const ::LinearAlgebra::distributed::Vector<double> &src)
1959 if (!preconditioner->UseTranspose())
1961 ierr = preconditioner->SetUseTranspose(
true);
1966 ierr = preconditioner->SetUseTranspose(
false);
1976 ExcNonMatchingMaps(
"dst"));
1978 ExcNonMatchingMaps(
"src"));
1989 ExcNonMatchingMaps(
"dst"));
1991 ExcNonMatchingMaps(
"src"));
1993 preconditioner->SetUseTranspose(
true);
1997 preconditioner->SetUseTranspose(
false);
2012 const ::Vector<double> &src)
const 2015 preconditioner->OperatorDomainMap().NumMyElements());
2016 AssertDimension(static_cast<TrilinosWrappers::types::int_type>(src.size()),
2017 preconditioner->OperatorRangeMap().NumMyElements());
2018 Epetra_Vector tril_dst(View,
2019 preconditioner->OperatorDomainMap(),
2021 Epetra_Vector tril_src(View,
2022 preconditioner->OperatorRangeMap(),
2023 const_cast<double *
>(src.begin()));
2025 const int ierr = preconditioner->ApplyInverse(tril_src, tril_dst);
2032 const ::Vector<double> &src)
const 2035 preconditioner->OperatorDomainMap().NumMyElements());
2036 AssertDimension(static_cast<TrilinosWrappers::types::int_type>(src.size()),
2037 preconditioner->OperatorRangeMap().NumMyElements());
2038 Epetra_Vector tril_dst(View,
2039 preconditioner->OperatorDomainMap(),
2041 Epetra_Vector tril_src(View,
2042 preconditioner->OperatorRangeMap(),
2043 const_cast<double *
>(src.begin()));
2045 preconditioner->SetUseTranspose(
true);
2046 const int ierr = preconditioner->ApplyInverse(tril_src, tril_dst);
2048 preconditioner->SetUseTranspose(
false);
2054 PreconditionBase::vmult(
2060 preconditioner->OperatorDomainMap().NumMyElements());
2063 preconditioner->OperatorRangeMap().NumMyElements());
2064 Epetra_Vector tril_dst(View,
2065 preconditioner->OperatorDomainMap(),
2067 Epetra_Vector tril_src(View,
2068 preconditioner->OperatorRangeMap(),
2069 const_cast<double *
>(src.
begin()));
2071 const int ierr = preconditioner->ApplyInverse(tril_src, tril_dst);
2076 PreconditionBase::Tvmult(
2082 preconditioner->OperatorDomainMap().NumMyElements());
2085 preconditioner->OperatorRangeMap().NumMyElements());
2086 Epetra_Vector tril_dst(View,
2087 preconditioner->OperatorDomainMap(),
2089 Epetra_Vector tril_src(View,
2090 preconditioner->OperatorRangeMap(),
2091 const_cast<double *
>(src.
begin()));
2093 preconditioner->SetUseTranspose(
true);
2094 const int ierr = preconditioner->ApplyInverse(tril_src, tril_dst);
2096 preconditioner->SetUseTranspose(
false);
2107 DEAL_II_NAMESPACE_CLOSE
2109 # endif // DEAL_II_WITH_TRILINOS 2111 #endif // trilinos_precondition_h
static::ExceptionBase & ExcTrilinosError(int arg1)
unsigned int smoother_sweeps
#define AssertDimension(dim1, dim2)
::types::global_dof_index size_type
double aggregation_threshold
std::vector< std::vector< bool > > constant_modes
unsigned int smoother_overlap
const char * smoother_type
std::shared_ptr< Epetra_Map > vector_distributor
#define AssertThrow(cond, exc)
unsigned int smoother_sweeps
size_type local_size() const
unsigned long long int global_dof_index
std::string block_creation_type
#define DeclException1(Exception1, type1, outsequence)
double aggregation_threshold
const char * smoother_type
#define Assert(cond, exc)
std::vector< std::vector< bool > > constant_modes
const Epetra_Map & vector_partitioner() const
SymmetricTensor< rank_, dim, Number > transpose(const SymmetricTensor< rank_, dim, Number > &t)
std::shared_ptr< SparseMatrix > trilinos_matrix
std::shared_ptr< SparseMatrix > trilinos_matrix
std::string block_creation_type
std::shared_ptr< Epetra_Operator > preconditioner
unsigned int smoother_overlap
Epetra_MpiComm communicator
const Epetra_MultiVector & trilinos_vector() const
bool higher_order_elements
std::string block_creation_type