Reference documentation for deal.II version GIT 9042b9283b 2023-12-02 14:50:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
trilinos_precondition.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2008 - 2023 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_trilinos_precondition_h
17 #define dealii_trilinos_precondition_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #ifdef DEAL_II_WITH_TRILINOS
23 
25 
28 
29 # include <Epetra_Map.h>
30 # include <Epetra_MpiComm.h>
31 # include <Epetra_MultiVector.h>
32 # include <Epetra_RowMatrix.h>
33 # include <Epetra_Vector.h>
34 # include <Teuchos_ParameterList.hpp>
35 
36 # include <memory>
37 
38 // forward declarations
39 # ifndef DOXYGEN
40 class Ifpack_Preconditioner;
41 class Ifpack_Chebyshev;
42 namespace ML_Epetra
43 {
44  class MultiLevelPreconditioner;
45 }
46 # endif
47 
49 
50 // forward declarations
51 # ifndef DOXYGEN
52 template <typename number>
53 class SparseMatrix;
54 template <typename number>
55 class Vector;
56 class SparsityPattern;
57 # endif
58 
64 namespace TrilinosWrappers
65 {
66  // forward declarations
67  class SparseMatrix;
68  class BlockSparseMatrix;
69  class SolverBase;
70 
78  {
79  public:
84 
90  {};
91 
98 
103  void
104  clear();
105 
109  MPI_Comm
110  get_mpi_communicator() const;
111 
121  void
123 
127  virtual void
128  vmult(MPI::Vector &dst, const MPI::Vector &src) const;
129 
133  virtual void
134  Tvmult(MPI::Vector &dst, const MPI::Vector &src) const;
135 
140  virtual void
141  vmult(::Vector<double> &dst, const ::Vector<double> &src) const;
142 
147  virtual void
149  const ::Vector<double> &src) const;
150 
155  virtual void
157  const ::LinearAlgebra::distributed::Vector<double> &src) const;
158 
163  virtual void
165  const ::LinearAlgebra::distributed::Vector<double> &src) const;
166 
177  trilinos_operator() const;
189  IndexSet
191 
197  IndexSet
199 
210  std::string,
211  << "The sparse matrix the preconditioner is based on "
212  << "uses a map that is not compatible to the one in vector "
213  << arg1 << ". Check preconditioner and matrix setup.");
216  friend class SolverBase;
217 
218  protected:
223  Teuchos::RCP<Epetra_Operator> preconditioner;
224 
229  Epetra_MpiComm communicator;
230  };
231 
232 
249  {
250  public:
263  {
268  AdditionalData(const double omega = 1,
269  const double min_diagonal = 0,
270  const unsigned int n_sweeps = 1);
271 
275  double omega;
276 
284  double min_diagonal;
285 
290  unsigned int n_sweeps;
291  };
292 
297  void
299  const AdditionalData &additional_data = AdditionalData());
300  };
301 
302 
303 
330  {
331  public:
346  {
353  AdditionalData(const double omega = 1,
354  const double min_diagonal = 0,
355  const unsigned int overlap = 0,
356  const unsigned int n_sweeps = 1);
357 
362  double omega;
363 
371  double min_diagonal;
372 
377  unsigned int overlap;
378 
383  unsigned int n_sweeps;
384  };
385 
391  void
393  const AdditionalData &additional_data = AdditionalData());
394  };
395 
396 
397 
424  {
425  public:
440  {
447  AdditionalData(const double omega = 1,
448  const double min_diagonal = 0,
449  const unsigned int overlap = 0,
450  const unsigned int n_sweeps = 1);
451 
456  double omega;
457 
465  double min_diagonal;
466 
471  unsigned int overlap;
472 
477  unsigned int n_sweeps;
478  };
479 
485  void
487  const AdditionalData &additional_data = AdditionalData());
488  };
489 
490 
491 
509  {
510  public:
527  {
533  AdditionalData(const unsigned int block_size = 1,
534  const std::string &block_creation_type = "linear",
535  const double omega = 1,
536  const double min_diagonal = 0,
537  const unsigned int n_sweeps = 1);
538 
542  unsigned int block_size;
543 
551  std::string block_creation_type;
552 
557  double omega;
558 
566  double min_diagonal;
567 
572  unsigned int n_sweeps;
573  };
574 
579  void
581  const AdditionalData &additional_data = AdditionalData());
582  };
583 
584 
585 
608  {
609  public:
628  {
636  AdditionalData(const unsigned int block_size = 1,
637  const std::string &block_creation_type = "linear",
638  const double omega = 1,
639  const double min_diagonal = 0,
640  const unsigned int overlap = 0,
641  const unsigned int n_sweeps = 1);
642 
646  unsigned int block_size;
647 
655  std::string block_creation_type;
656 
661  double omega;
662 
670  double min_diagonal;
671 
676  unsigned int overlap;
677 
682  unsigned int n_sweeps;
683  };
684 
690  void
692  const AdditionalData &additional_data = AdditionalData());
693  };
694 
695 
696 
719  {
720  public:
739  {
747  AdditionalData(const unsigned int block_size = 1,
748  const std::string &block_creation_type = "linear",
749  const double omega = 1,
750  const double min_diagonal = 0,
751  const unsigned int overlap = 0,
752  const unsigned int n_sweeps = 1);
753 
757  unsigned int block_size;
758 
766  std::string block_creation_type;
767 
772  double omega;
773 
781  double min_diagonal;
782 
787  unsigned int overlap;
788 
793  unsigned int n_sweeps;
794  };
795 
801  void
803  const AdditionalData &additional_data = AdditionalData());
804  };
805 
806 
807 
843  {
844  public:
862  {
872  AdditionalData(const unsigned int ic_fill = 0,
873  const double ic_atol = 0.,
874  const double ic_rtol = 1.,
875  const unsigned int overlap = 0);
876 
885  unsigned int ic_fill;
886 
892  double ic_atol;
893 
898  double ic_rtol;
899 
904  unsigned int overlap;
905  };
906 
911  void
913  const AdditionalData &additional_data = AdditionalData());
914  };
915 
916 
917 
947  {
948  public:
987  {
991  AdditionalData(const unsigned int ilu_fill = 0,
992  const double ilu_atol = 0.,
993  const double ilu_rtol = 1.,
994  const unsigned int overlap = 0);
995 
999  unsigned int ilu_fill;
1000 
1005  double ilu_atol;
1006 
1011  double ilu_rtol;
1012 
1016  unsigned int overlap;
1017  };
1018 
1023  void
1025  const AdditionalData &additional_data = AdditionalData());
1026  };
1027 
1028 
1029 
1065  {
1066  public:
1085  {
1096  AdditionalData(const double ilut_drop = 0.,
1097  const unsigned int ilut_fill = 0,
1098  const double ilut_atol = 0.,
1099  const double ilut_rtol = 1.,
1100  const unsigned int overlap = 0);
1101 
1106  double ilut_drop;
1107 
1116  unsigned int ilut_fill;
1117 
1123  double ilut_atol;
1124 
1129  double ilut_rtol;
1130 
1135  unsigned int overlap;
1136  };
1137 
1142  void
1144  const AdditionalData &additional_data = AdditionalData());
1145  };
1146 
1147 
1148 
1167  {
1168  public:
1174  {
1178  AdditionalData(const unsigned int overlap = 0);
1179 
1180 
1185  unsigned int overlap;
1186  };
1187 
1192  void
1194  const AdditionalData &additional_data = AdditionalData());
1195  };
1196 
1197 
1198 
1208  {
1209  public:
1215  {
1219  AdditionalData(const unsigned int degree = 1,
1220  const double max_eigenvalue = 10.,
1221  const double eigenvalue_ratio = 30.,
1222  const double min_eigenvalue = 1.,
1223  const double min_diagonal = 1e-12,
1224  const bool nonzero_starting = false);
1225 
1231  unsigned int degree;
1232 
1238 
1243 
1249 
1255 
1267  };
1268 
1273  void
1275  const AdditionalData &additional_data = AdditionalData());
1276  };
1277 
1278 
1279 
1322  {
1323  public:
1331  {
1355  AdditionalData(const bool elliptic = true,
1356  const bool higher_order_elements = false,
1357  const unsigned int n_cycles = 1,
1358  const bool w_cycle = false,
1359  const double aggregation_threshold = 1e-4,
1360  const std::vector<std::vector<bool>> &constant_modes =
1361  std::vector<std::vector<bool>>(0),
1362  const unsigned int smoother_sweeps = 2,
1363  const unsigned int smoother_overlap = 0,
1364  const bool output_details = false,
1365  const char *smoother_type = "Chebyshev",
1366  const char *coarse_type = "Amesos-KLU");
1367 
1398  void
1400  Teuchos::ParameterList &parameter_list,
1401  std::unique_ptr<Epetra_MultiVector> &distributed_constant_modes,
1402  const Epetra_RowMatrix &matrix) const;
1403 
1411  void
1413  Teuchos::ParameterList &parameter_list,
1414  std::unique_ptr<Epetra_MultiVector> &distributed_constant_modes,
1415  const SparseMatrix &matrix) const;
1416 
1422  void
1424  Teuchos::ParameterList &parameter_list,
1425  std::unique_ptr<Epetra_MultiVector> &distributed_constant_modes,
1426  const Epetra_RowMatrix &matrix) const;
1427 
1433  void
1435  Teuchos::ParameterList &parameter_list,
1436  std::unique_ptr<Epetra_MultiVector> &distributed_constant_modes,
1437  const SparseMatrix &matrix) const;
1438 
1446  bool elliptic;
1447 
1453 
1458  unsigned int n_cycles;
1459 
1464  bool w_cycle;
1465 
1476 
1501  std::vector<std::vector<bool>> constant_modes;
1502 
1513  unsigned int smoother_sweeps;
1514 
1519  unsigned int smoother_overlap;
1520 
1527 
1562  const char *smoother_type;
1563 
1568  const char *coarse_type;
1569  };
1570 
1574  ~PreconditionAMG() override;
1575 
1576 
1582  void
1584  const AdditionalData &additional_data = AdditionalData());
1585 
1604  void
1605  initialize(const Epetra_RowMatrix &matrix,
1606  const AdditionalData &additional_data = AdditionalData());
1607 
1622  void
1624  const Teuchos::ParameterList &ml_parameters);
1625 
1633  void
1634  initialize(const Epetra_RowMatrix &matrix,
1635  const Teuchos::ParameterList &ml_parameters);
1636 
1643  template <typename number>
1644  void
1645  initialize(const ::SparseMatrix<number> &deal_ii_sparse_matrix,
1646  const AdditionalData &additional_data = AdditionalData(),
1647  const double drop_tolerance = 1e-13,
1648  const ::SparsityPattern *use_this_sparsity = nullptr);
1649 
1662  void
1663  reinit();
1664 
1669  void
1670  clear();
1671 
1675  size_type
1676  memory_consumption() const;
1677 
1678  private:
1682  std::shared_ptr<SparseMatrix> trilinos_matrix;
1683  };
1684 
1685 
1686 
1687 # if defined(DOXYGEN) || defined(DEAL_II_TRILINOS_WITH_MUELU)
1707  {
1708  public:
1716  {
1721  AdditionalData(const bool elliptic = true,
1722  const unsigned int n_cycles = 1,
1723  const bool w_cycle = false,
1724  const double aggregation_threshold = 1e-4,
1725  const std::vector<std::vector<bool>> &constant_modes =
1726  std::vector<std::vector<bool>>(0),
1727  const unsigned int smoother_sweeps = 2,
1728  const unsigned int smoother_overlap = 0,
1729  const bool output_details = false,
1730  const char *smoother_type = "Chebyshev",
1731  const char *coarse_type = "Amesos-KLU");
1732 
1740  bool elliptic;
1741 
1746  unsigned int n_cycles;
1747 
1752  bool w_cycle;
1753 
1764 
1771  std::vector<std::vector<bool>> constant_modes;
1772 
1783  unsigned int smoother_sweeps;
1784 
1789  unsigned int smoother_overlap;
1790 
1797 
1832  const char *smoother_type;
1833 
1838  const char *coarse_type;
1839  };
1840 
1845 
1849  virtual ~PreconditionAMGMueLu() override = default;
1850 
1856  void
1858  const AdditionalData &additional_data = AdditionalData());
1859 
1866  void
1867  initialize(const Epetra_CrsMatrix &matrix,
1868  const AdditionalData &additional_data = AdditionalData());
1869 
1883  void
1885  Teuchos::ParameterList &muelu_parameters);
1886 
1892  void
1893  initialize(const Epetra_CrsMatrix &matrix,
1894  Teuchos::ParameterList &muelu_parameters);
1895 
1902  template <typename number>
1903  void
1904  initialize(const ::SparseMatrix<number> &deal_ii_sparse_matrix,
1905  const AdditionalData &additional_data = AdditionalData(),
1906  const double drop_tolerance = 1e-13,
1907  const ::SparsityPattern *use_this_sparsity = nullptr);
1908 
1913  void
1914  clear();
1915 
1919  size_type
1920  memory_consumption() const;
1921 
1922  private:
1926  std::shared_ptr<SparseMatrix> trilinos_matrix;
1927  };
1928 # endif
1929 
1930 
1931 
1939  {
1940  public:
1946  {};
1947 
1954  void
1956  const AdditionalData &additional_data = AdditionalData());
1957 
1961  void
1962  vmult(MPI::Vector &dst, const MPI::Vector &src) const override;
1963 
1967  void
1968  Tvmult(MPI::Vector &dst, const MPI::Vector &src) const override;
1969 
1974  void
1975  vmult(::Vector<double> &dst,
1976  const ::Vector<double> &src) const override;
1977 
1982  void
1983  Tvmult(::Vector<double> &dst,
1984  const ::Vector<double> &src) const override;
1985 
1990  void
1992  const ::LinearAlgebra::distributed::Vector<double> &src)
1993  const override;
1994 
2000  void
2002  const ::LinearAlgebra::distributed::Vector<double> &src)
2003  const override;
2004  };
2005 
2006 
2007 
2008  // ----------------------- inline and template functions --------------------
2009 
2010 
2011 # ifndef DOXYGEN
2012 
2013 
2014  inline void
2016  {
2017  // This only flips a flag that tells
2018  // Trilinos that any vmult operation
2019  // should be done with the
2020  // transpose. However, the matrix
2021  // structure is not reset.
2022  int ierr;
2023 
2024  if (!preconditioner->UseTranspose())
2025  {
2026  ierr = preconditioner->SetUseTranspose(true);
2027  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
2028  }
2029  else
2030  {
2031  ierr = preconditioner->SetUseTranspose(false);
2032  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
2033  }
2034  }
2035 
2036 
2037  inline void
2038  PreconditionBase::vmult(MPI::Vector &dst, const MPI::Vector &src) const
2039  {
2040  Assert(dst.trilinos_partitioner().SameAs(
2041  preconditioner->OperatorRangeMap()),
2042  ExcNonMatchingMaps("dst"));
2043  Assert(src.trilinos_partitioner().SameAs(
2044  preconditioner->OperatorDomainMap()),
2045  ExcNonMatchingMaps("src"));
2046 
2047  const int ierr = preconditioner->ApplyInverse(src.trilinos_vector(),
2048  dst.trilinos_vector());
2049  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
2050  }
2051 
2052  inline void
2053  PreconditionBase::Tvmult(MPI::Vector &dst, const MPI::Vector &src) const
2054  {
2055  Assert(dst.trilinos_partitioner().SameAs(
2056  preconditioner->OperatorRangeMap()),
2057  ExcNonMatchingMaps("dst"));
2058  Assert(src.trilinos_partitioner().SameAs(
2059  preconditioner->OperatorDomainMap()),
2060  ExcNonMatchingMaps("src"));
2061 
2062  preconditioner->SetUseTranspose(true);
2063  const int ierr = preconditioner->ApplyInverse(src.trilinos_vector(),
2064  dst.trilinos_vector());
2065  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
2066  preconditioner->SetUseTranspose(false);
2067  }
2068 
2069 
2070  // For the implementation of the <code>vmult</code> function with deal.II
2071  // data structures we note that invoking a call of the Trilinos
2072  // preconditioner requires us to use Epetra vectors as well. We do this by
2073  // providing a view, i.e., feed Trilinos with a pointer to the data, so we
2074  // avoid copying the content of the vectors during the iteration (this
2075  // function is only useful when used in serial anyway). In the declaration
2076  // of the right hand side, we need to cast the source vector (that is
2077  // <code>const</code> in all deal.II calls) to non-constant value, as this
2078  // is the way Trilinos wants to have them.
2079  inline void
2081  const ::Vector<double> &src) const
2082  {
2083  AssertDimension(dst.size(),
2084  preconditioner->OperatorDomainMap().NumMyElements());
2085  AssertDimension(src.size(),
2086  preconditioner->OperatorRangeMap().NumMyElements());
2087  Epetra_Vector tril_dst(View,
2088  preconditioner->OperatorDomainMap(),
2089  dst.begin());
2090  Epetra_Vector tril_src(View,
2091  preconditioner->OperatorRangeMap(),
2092  const_cast<double *>(src.begin()));
2093 
2094  const int ierr = preconditioner->ApplyInverse(tril_src, tril_dst);
2095  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
2096  }
2097 
2098 
2099  inline void
2101  const ::Vector<double> &src) const
2102  {
2103  AssertDimension(dst.size(),
2104  preconditioner->OperatorDomainMap().NumMyElements());
2105  AssertDimension(src.size(),
2106  preconditioner->OperatorRangeMap().NumMyElements());
2107  Epetra_Vector tril_dst(View,
2108  preconditioner->OperatorDomainMap(),
2109  dst.begin());
2110  Epetra_Vector tril_src(View,
2111  preconditioner->OperatorRangeMap(),
2112  const_cast<double *>(src.begin()));
2113 
2114  preconditioner->SetUseTranspose(true);
2115  const int ierr = preconditioner->ApplyInverse(tril_src, tril_dst);
2116  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
2117  preconditioner->SetUseTranspose(false);
2118  }
2119 
2120 
2121 
2122  inline void
2126  {
2128  preconditioner->OperatorDomainMap().NumMyElements());
2130  preconditioner->OperatorRangeMap().NumMyElements());
2131  Epetra_Vector tril_dst(View,
2132  preconditioner->OperatorDomainMap(),
2133  dst.begin());
2134  Epetra_Vector tril_src(View,
2135  preconditioner->OperatorRangeMap(),
2136  const_cast<double *>(src.begin()));
2137 
2138  const int ierr = preconditioner->ApplyInverse(tril_src, tril_dst);
2139  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
2140  }
2141 
2142  inline void
2146  {
2148  preconditioner->OperatorDomainMap().NumMyElements());
2150  preconditioner->OperatorRangeMap().NumMyElements());
2151  Epetra_Vector tril_dst(View,
2152  preconditioner->OperatorDomainMap(),
2153  dst.begin());
2154  Epetra_Vector tril_src(View,
2155  preconditioner->OperatorRangeMap(),
2156  const_cast<double *>(src.begin()));
2157 
2158  preconditioner->SetUseTranspose(true);
2159  const int ierr = preconditioner->ApplyInverse(tril_src, tril_dst);
2160  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
2161  preconditioner->SetUseTranspose(false);
2162  }
2163 
2164 # endif
2165 
2166 } // namespace TrilinosWrappers
2167 
2168 
2173 
2174 #endif // DEAL_II_WITH_TRILINOS
2175 
2176 #endif
size_type locally_owned_size() const
std::shared_ptr< SparseMatrix > trilinos_matrix
virtual ~PreconditionAMGMueLu() override=default
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
std::shared_ptr< SparseMatrix > trilinos_matrix
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
virtual void Tvmult(::LinearAlgebra::distributed::Vector< double > &dst, const ::LinearAlgebra::distributed::Vector< double > &src) const
virtual void Tvmult(MPI::Vector &dst, const MPI::Vector &src) const
virtual void vmult(MPI::Vector &dst, const MPI::Vector &src) const
Epetra_Operator & trilinos_operator() const
Teuchos::RCP< Epetra_Operator > preconditioner
virtual void Tvmult(::Vector< double > &dst, const ::Vector< double > &src) const
virtual void vmult(::Vector< double > &dst, const ::Vector< double > &src) const
virtual void vmult(::LinearAlgebra::distributed::Vector< double > &dst, const ::LinearAlgebra::distributed::Vector< double > &src) const
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())
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())
void vmult(LinearAlgebra::distributed::Vector< double > &dst, const ::LinearAlgebra::distributed::Vector< double > &src) const override
void vmult(MPI::Vector &dst, const MPI::Vector &src) const override
void Tvmult(MPI::Vector &dst, const MPI::Vector &src) const override
void Tvmult(LinearAlgebra::distributed::Vector< double > &dst, const ::LinearAlgebra::distributed::Vector< double > &src) const override
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())
Definition: vector.h:110
virtual size_type size() const override
iterator begin()
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcTrilinosError(int arg1)
static ::ExceptionBase & ExcNonMatchingMaps(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1631
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1820
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:517
#define AssertThrow(cond, exc)
Definition: exceptions.h:1732
typename ::internal::FEInterfaceViews::ViewType< dim, spacedim, Extractor >::type View
@ matrix
Contents is actually a matrix.
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
unsigned int global_dof_index
Definition: types.h:82
AdditionalData(const bool elliptic=true, const unsigned int n_cycles=1, const bool w_cycle=false, const double aggregation_threshold=1e-4, const std::vector< std::vector< bool >> &constant_modes=std::vector< std::vector< bool >>(0), const unsigned int smoother_sweeps=2, const unsigned int smoother_overlap=0, const bool output_details=false, const char *smoother_type="Chebyshev", const char *coarse_type="Amesos-KLU")
AdditionalData(const bool elliptic=true, const bool higher_order_elements=false, const unsigned int n_cycles=1, const bool w_cycle=false, const double aggregation_threshold=1e-4, const std::vector< std::vector< bool >> &constant_modes=std::vector< std::vector< bool >>(0), const unsigned int smoother_sweeps=2, const unsigned int smoother_overlap=0, const bool output_details=false, const char *smoother_type="Chebyshev", const char *coarse_type="Amesos-KLU")
void set_operator_null_space(Teuchos::ParameterList &parameter_list, std::unique_ptr< Epetra_MultiVector > &distributed_constant_modes, const Epetra_RowMatrix &matrix) const
void set_parameters(Teuchos::ParameterList &parameter_list, std::unique_ptr< Epetra_MultiVector > &distributed_constant_modes, const Epetra_RowMatrix &matrix) const
AdditionalData(const unsigned int block_size=1, const std::string &block_creation_type="linear", const double omega=1, const double min_diagonal=0, const unsigned int n_sweeps=1)
AdditionalData(const unsigned int block_size=1, const std::string &block_creation_type="linear", const double omega=1, const double min_diagonal=0, const unsigned int overlap=0, const unsigned int n_sweeps=1)
AdditionalData(const unsigned int block_size=1, const std::string &block_creation_type="linear", const double omega=1, const double min_diagonal=0, const unsigned int overlap=0, const unsigned int n_sweeps=1)
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)
AdditionalData(const unsigned int ic_fill=0, const double ic_atol=0., const double ic_rtol=1., const unsigned int overlap=0)
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)
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, const unsigned int n_sweeps=1)
AdditionalData(const double omega=1, const double min_diagonal=0, const unsigned int overlap=0, const unsigned int n_sweeps=1)
AdditionalData(const double omega=1, const double min_diagonal=0, const unsigned int overlap=0, const unsigned int n_sweeps=1)