18 #ifdef DEAL_II_WITH_TRILINOS
24 # include <Epetra_MultiVector.h>
26 # include <Ifpack_Chebyshev.h>
27 # include <Teuchos_ParameterList.hpp>
28 # include <Teuchos_RCP.hpp>
35 : communicator(MPI_COMM_SELF)
58 ExcMessage(
"Trying to dereference a null pointer."));
80 const double min_diagonal,
81 const unsigned int n_sweeps)
83 , min_diagonal(min_diagonal)
96 Ifpack().Create(
"point relaxation",
97 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
100 Ifpack_Preconditioner *ifpack =
108 Teuchos::ParameterList parameter_list;
109 parameter_list.set(
"relaxation: sweeps",
110 static_cast<int>(additional_data.
n_sweeps));
111 parameter_list.set(
"relaxation: type",
"Jacobi");
112 parameter_list.set(
"relaxation: damping factor", additional_data.
omega);
113 parameter_list.set(
"relaxation: min diagonal value",
116 ierr = ifpack->SetParameters(parameter_list);
119 ierr = ifpack->Initialize();
122 ierr = ifpack->Compute();
131 const double min_diagonal,
132 const unsigned int overlap,
133 const unsigned int n_sweeps)
135 , min_diagonal(min_diagonal)
148 Ifpack().Create(
"point relaxation",
149 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
152 Ifpack_Preconditioner *ifpack =
160 Teuchos::ParameterList parameter_list;
161 parameter_list.set(
"relaxation: sweeps",
162 static_cast<int>(additional_data.
n_sweeps));
163 parameter_list.set(
"relaxation: type",
"symmetric Gauss-Seidel");
164 parameter_list.set(
"relaxation: damping factor", additional_data.
omega);
165 parameter_list.set(
"relaxation: min diagonal value",
167 parameter_list.set(
"schwarz: combine mode",
"Add");
169 ierr = ifpack->SetParameters(parameter_list);
172 ierr = ifpack->Initialize();
175 ierr = ifpack->Compute();
184 const double min_diagonal,
185 const unsigned int overlap,
186 const unsigned int n_sweeps)
188 , min_diagonal(min_diagonal)
201 Ifpack().Create(
"point relaxation",
202 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
205 Ifpack_Preconditioner *ifpack =
213 Teuchos::ParameterList parameter_list;
214 parameter_list.set(
"relaxation: sweeps",
215 static_cast<int>(additional_data.
n_sweeps));
216 parameter_list.set(
"relaxation: type",
"Gauss-Seidel");
217 parameter_list.set(
"relaxation: damping factor", additional_data.
omega);
218 parameter_list.set(
"relaxation: min diagonal value",
220 parameter_list.set(
"schwarz: combine mode",
"Add");
222 ierr = ifpack->SetParameters(parameter_list);
225 ierr = ifpack->Initialize();
228 ierr = ifpack->Compute();
238 const std::string &block_creation_type,
240 const double min_diagonal,
241 const unsigned int n_sweeps)
243 , block_creation_type(block_creation_type)
245 , min_diagonal(min_diagonal)
261 Ifpack().Create((
matrix.trilinos_matrix().NumMyRows() == 0) ?
264 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
267 Ifpack_Preconditioner *ifpack =
275 Teuchos::ParameterList parameter_list;
276 parameter_list.set(
"relaxation: sweeps",
277 static_cast<int>(additional_data.
n_sweeps));
278 parameter_list.set(
"relaxation: type",
"Jacobi");
279 parameter_list.set(
"relaxation: damping factor", additional_data.
omega);
280 parameter_list.set(
"relaxation: min diagonal value",
282 parameter_list.set(
"partitioner: type",
287 parameter_list.set(
"partitioner: local parts", n_local_parts);
289 ierr = ifpack->SetParameters(parameter_list);
292 ierr = ifpack->Initialize();
295 ierr = ifpack->Compute();
305 const std::string &block_creation_type,
307 const double min_diagonal,
308 const unsigned int overlap,
309 const unsigned int n_sweeps)
311 , block_creation_type(block_creation_type)
313 , min_diagonal(min_diagonal)
329 Ifpack().Create((
matrix.trilinos_matrix().NumMyRows() == 0) ?
332 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
335 Ifpack_Preconditioner *ifpack =
343 Teuchos::ParameterList parameter_list;
344 parameter_list.set(
"relaxation: sweeps",
345 static_cast<int>(additional_data.
n_sweeps));
346 parameter_list.set(
"relaxation: type",
"symmetric Gauss-Seidel");
347 parameter_list.set(
"relaxation: damping factor", additional_data.
omega);
348 parameter_list.set(
"relaxation: min diagonal value",
350 parameter_list.set(
"schwarz: combine mode",
"Add");
351 parameter_list.set(
"partitioner: type",
356 parameter_list.set(
"partitioner: local parts", n_local_parts);
358 ierr = ifpack->SetParameters(parameter_list);
361 ierr = ifpack->Initialize();
364 ierr = ifpack->Compute();
374 const std::string &block_creation_type,
376 const double min_diagonal,
377 const unsigned int overlap,
378 const unsigned int n_sweeps)
380 , block_creation_type(block_creation_type)
382 , min_diagonal(min_diagonal)
398 Ifpack().Create((
matrix.trilinos_matrix().NumMyRows() == 0) ?
401 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
404 Ifpack_Preconditioner *ifpack =
412 Teuchos::ParameterList parameter_list;
413 parameter_list.set(
"relaxation: sweeps",
414 static_cast<int>(additional_data.
n_sweeps));
415 parameter_list.set(
"relaxation: type",
"Gauss-Seidel");
416 parameter_list.set(
"relaxation: damping factor", additional_data.
omega);
417 parameter_list.set(
"relaxation: min diagonal value",
419 parameter_list.set(
"schwarz: combine mode",
"Add");
420 parameter_list.set(
"partitioner: type",
425 parameter_list.set(
"partitioner: local parts", n_local_parts);
427 ierr = ifpack->SetParameters(parameter_list);
430 ierr = ifpack->Initialize();
433 ierr = ifpack->Compute();
442 const double ic_atol,
443 const double ic_rtol,
444 const unsigned int overlap)
459 Ifpack().Create(
"IC",
460 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
463 Ifpack_Preconditioner *ifpack =
471 Teuchos::ParameterList parameter_list;
472 parameter_list.set(
"fact: level-of-fill", additional_data.
ic_fill);
473 parameter_list.set(
"fact: absolute threshold", additional_data.
ic_atol);
474 parameter_list.set(
"fact: relative threshold", additional_data.
ic_rtol);
475 parameter_list.set(
"schwarz: combine mode",
"Add");
477 ierr = ifpack->SetParameters(parameter_list);
480 ierr = ifpack->Initialize();
483 ierr = ifpack->Compute();
492 const double ilu_atol,
493 const double ilu_rtol,
494 const unsigned int overlap)
509 Ifpack().Create(
"ILU",
510 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
513 Ifpack_Preconditioner *ifpack =
521 Teuchos::ParameterList parameter_list;
522 parameter_list.set(
"fact: level-of-fill",
523 static_cast<int>(additional_data.
ilu_fill));
524 parameter_list.set(
"fact: absolute threshold", additional_data.
ilu_atol);
525 parameter_list.set(
"fact: relative threshold", additional_data.
ilu_rtol);
526 parameter_list.set(
"schwarz: combine mode",
"Add");
528 ierr = ifpack->SetParameters(parameter_list);
531 ierr = ifpack->Initialize();
534 ierr = ifpack->Compute();
543 const unsigned int ilut_fill,
544 const double ilut_atol,
545 const double ilut_rtol,
546 const unsigned int overlap)
547 : ilut_drop(ilut_drop)
548 , ilut_fill(ilut_fill)
549 , ilut_atol(ilut_atol)
550 , ilut_rtol(ilut_rtol)
562 Ifpack().Create(
"ILUT",
563 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
566 Ifpack_Preconditioner *ifpack =
574 Teuchos::ParameterList parameter_list;
575 parameter_list.set(
"fact: drop value", additional_data.
ilut_drop);
576 parameter_list.set(
"fact: level-of-fill",
577 static_cast<int>(additional_data.
ilut_fill));
578 parameter_list.set(
"fact: absolute threshold", additional_data.
ilut_atol);
579 parameter_list.set(
"fact: relative threshold", additional_data.
ilut_rtol);
580 parameter_list.set(
"schwarz: combine mode",
"Add");
582 ierr = ifpack->SetParameters(parameter_list);
585 ierr = ifpack->Initialize();
588 ierr = ifpack->Compute();
597 const unsigned int overlap)
609 Ifpack().Create(
"Amesos",
610 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
613 Ifpack_Preconditioner *ifpack =
621 Teuchos::ParameterList parameter_list;
622 parameter_list.set(
"schwarz: combine mode",
"Add");
624 ierr = ifpack->SetParameters(parameter_list);
627 ierr = ifpack->Initialize();
630 ierr = ifpack->Compute();
639 const unsigned int degree,
640 const double max_eigenvalue,
641 const double eigenvalue_ratio,
642 const double min_eigenvalue,
643 const double min_diagonal,
644 const bool nonzero_starting)
646 , max_eigenvalue(max_eigenvalue)
647 , eigenvalue_ratio(eigenvalue_ratio)
648 , min_eigenvalue(min_eigenvalue)
649 , min_diagonal(min_diagonal)
650 , nonzero_starting(nonzero_starting)
660 Teuchos::rcp(
new Ifpack_Chebyshev(&
matrix.trilinos_matrix()));
662 Ifpack_Chebyshev *ifpack =
670 Teuchos::ParameterList parameter_list;
671 parameter_list.set(
"chebyshev: ratio eigenvalue",
673 parameter_list.set(
"chebyshev: min eigenvalue",
675 parameter_list.set(
"chebyshev: max eigenvalue",
677 parameter_list.set(
"chebyshev: degree",
678 static_cast<int>(additional_data.
degree));
679 parameter_list.set(
"chebyshev: min diagonal value",
681 parameter_list.set(
"chebyshev: zero starting solution",
684 ierr = ifpack->SetParameters(parameter_list);
687 ierr = ifpack->Initialize();
690 ierr = ifpack->Compute();
714 Ifpack().Create(
"point relaxation",
715 const_cast<Epetra_CrsMatrix *
>(&
matrix.trilinos_matrix()),
718 Ifpack_Preconditioner *ifpack =
726 Teuchos::ParameterList parameter_list;
727 parameter_list.set(
"relaxation: sweeps", 1);
728 parameter_list.set(
"relaxation: type",
"Jacobi");
729 parameter_list.set(
"relaxation: damping factor", 1.0);
730 parameter_list.set(
"relaxation: min diagonal value", 0.0);
732 ierr = ifpack->SetParameters(parameter_list);
735 ierr = ifpack->Initialize();
738 ierr = ifpack->Compute();
756 const ::Vector<double> &src)
const
763 const ::Vector<double> &src)
const
IndexSet locally_owned_range_indices() const
Epetra_MpiComm communicator
Epetra_Operator & trilinos_operator() const
Teuchos::RCP< Epetra_Operator > preconditioner
MPI_Comm get_mpi_communicator() const
IndexSet locally_owned_domain_indices() 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(MPI::Vector &dst, const MPI::Vector &src) const override
void Tvmult(MPI::Vector &dst, const MPI::Vector &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())
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcTrilinosError(int arg1)
#define Assert(cond, exc)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
@ matrix
Contents is actually a matrix.
std::string block_creation_type
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)
std::string block_creation_type
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)
std::string block_creation_type
AdditionalData(const unsigned int overlap=0)
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)