Reference documentation for deal.II version 9.6.0
\(\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\}}\)
Loading...
Searching...
No Matches
trilinos_tpetra_precondition.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_trilinos_tpetra_precondition_h
16#define dealii_trilinos_tpetra_precondition_h
17
18
19#include <deal.II/base/config.h>
20
23
24#include "deal.II/lac/vector.h"
25
26#include <Teuchos_BLAS_types.hpp>
27
28
29#ifdef DEAL_II_TRILINOS_WITH_TPETRA
30
32
35
36# include <Teuchos_ConfigDefs.hpp>
37# include <Teuchos_ParameterList.hpp>
38# include <Teuchos_RCPDecl.hpp>
39# include <Tpetra_Operator.hpp>
40
41
43
44namespace LinearAlgebra
45{
46 namespace TpetraWrappers
47 {
48
55 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
57 {
58 public:
65 {};
66
73 PreconditionBase() = default;
74
80 void
82
89 virtual void
91 const Vector<Number, MemorySpace> &src) const;
92
99 virtual void
101 const Vector<Number, MemorySpace> &src) const;
102
109 virtual void
111
115 virtual void
117
126
133 Teuchos::RCP<TpetraTypes::LinearOperator<Number, MemorySpace>>
135
149
159
172 std::string,
173 << "The sparse matrix the preconditioner is based on "
174 << "uses a map that is not compatible to the one in vector " << arg1
175 << ". Check preconditioner and matrix setup.");
176
182 "The chosen preconditioner does not support transposing the matrix.");
185 protected:
190 Teuchos::RCP<TpetraTypes::LinearOperator<Number, MemorySpace>>
192
204 Teuchos::ParameterList parameter_list;
205 };
206
207
208
209# ifdef DEAL_II_TRILINOS_WITH_IFPACK2
219 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
220 class PreconditionIdentity : public PreconditionBase<Number, MemorySpace>
221 {
222 public:
227 PreconditionIdentity() = default;
228
238 void
239 initialize(const SparseMatrix<Number, MemorySpace> &A);
240 };
241
242
243
250 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
251 class PreconditionIfpackBase : public PreconditionBase<Number, MemorySpace>
252 {
253 public:
263 PreconditionIfpackBase(const std::string &preconditioner_type);
264
270 void
271 initialize(const SparseMatrix<Number, MemorySpace> &A);
272
277 ExcTrilinosIpack2PreconditionerUnsupported,
278 std::string,
279 << "You tried to select the preconditioner type <" << arg1 << ">\n"
280 << "but this preconditioner is not supported by Trilinos/Ifpack22\n"
281 << "due to one of the following reasons:\n"
282 << "* This preconditioner does not exist\n"
283 << "* This preconditioner has a specialized constructor not supported by the Ifpack2 Factory.\n"
284 << "* This preconditioner is not (yet) supported by Trilinos/Ifpack2\n"
285 << "* Trilinos/Ifpack2 was not configured for its use.");
286
287 protected:
292 std::string preconditioner_type;
293 };
294
295
296
310 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
311 class PreconditionIfpack
312 : public PreconditionIfpackBase<Number, MemorySpace>
313 {
314 public:
331 PreconditionIfpack(const std::string &preconditioner_type);
332
342 void
343 set_parameter_list(Teuchos::ParameterList &parameter_list);
344 };
345
346
347
355 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
357 : public PreconditionIfpackBase<Number, MemorySpace>
358 {
359 public:
364 struct AdditionalData
365 {
377 AdditionalData(const double omega = 1.,
378 const bool fix_diagonal = false,
379 const double min_diagonal = 0.,
380 const int n_sweeps = 1);
381
386 double omega;
395 bool fix_diagonal;
403 double min_diagonal;
408 int n_sweeps;
409 };
415
422 void
423 initialize(const SparseMatrix<Number, MemorySpace> &A,
424 const AdditionalData &additional_data = AdditionalData());
425 };
426
427
428
439 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
440 class PreconditionL1Jacobi
441 : public PreconditionIfpackBase<Number, MemorySpace>
442 {
443 public:
448 struct AdditionalData
449 {
462 AdditionalData(const double omega = 1.,
463 const bool fix_diagonal = false,
464 const double min_diagonal = 0.,
465 const int n_sweeps = 1);
466
471 double omega;
480 bool fix_diagonal;
488 double min_diagonal;
493 int n_sweeps;
494 };
495
500 PreconditionL1Jacobi();
501
508 void
509 initialize(const SparseMatrix<Number, MemorySpace> &A,
510 const AdditionalData &additional_data = AdditionalData());
511 };
512
513
514
526 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
527 class PreconditionL1GaussSeidel
528 : public PreconditionIfpackBase<Number, MemorySpace>
529 {
530 public:
535 struct AdditionalData
536 {
550 AdditionalData(const double omega = 1,
551 const double eta = 1.5,
552 const bool fix_diagonal = false,
553 const double min_diagonal = 0,
554 const int n_sweeps = 1);
555
560 double omega;
561
569 double eta;
578 bool fix_diagonal;
586 double min_diagonal;
591 int n_sweeps;
592 };
593
598 PreconditionL1GaussSeidel();
599
606 void
607 initialize(const SparseMatrix<Number, MemorySpace> &A,
608 const AdditionalData &additional_data = AdditionalData());
609 };
610
611
612
622 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
623 class PreconditionSOR : public PreconditionIfpackBase<Number, MemorySpace>
624 {
625 public:
630 struct AdditionalData
631 {
644 AdditionalData(const double omega = 1.,
645 const int overlap = 0,
646 const bool fix_diagonal = false,
647 const double min_diagonal = 0.,
648 const int n_sweeps = 1);
649
654 double omega;
664 int overlap;
673 bool fix_diagonal;
681 double min_diagonal;
686 int n_sweeps;
687 };
688
694
701 void
702 initialize(const SparseMatrix<Number, MemorySpace> &A,
703 const AdditionalData &additional_data = AdditionalData());
704 };
705
706
707
714 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
715 class PreconditionSSOR : public PreconditionIfpackBase<Number, MemorySpace>
716 {
717 public:
722 struct AdditionalData
723 {
736 AdditionalData(const double omega = 1.,
737 const int overlap = 0,
738 const bool fix_diagonal = false,
739 const double min_diagonal = 0.,
740 const int n_sweeps = 1);
741
746 double omega;
756 int overlap;
765 bool fix_diagonal;
773 double min_diagonal;
778 int n_sweeps;
779 };
780
781
787
794 void
795 initialize(const SparseMatrix<Number, MemorySpace> &A,
796 const AdditionalData &additional_data = AdditionalData());
797 };
798
799
800
807 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
809 : public PreconditionIfpackBase<Number, MemorySpace>
810 {
811 public:
816 struct AdditionalData
817 {
828 AdditionalData(const int degree = 1,
829 const double max_eigenvalue = 10.,
830 const double min_eigenvalue = 1.,
831 const double eigenvalue_ratio = 30.,
832 const double min_diagonal = 1e-12,
833 const bool nonzero_starting = false);
842 int degree;
843 /*
844 * @brief Upper bound for the maximum eigenvalue of the matrix.
845 *
846 * This needs to be set properly for the appropriate performance of this
847 * preconditioner.
848 */
849 double max_eigenvalue;
850 /*
851 * @brief Lower bound for the minimum eigenvalue of the matrix.
852 *
853 */
854 double min_eigenvalue;
859 double eigenvalue_ratio;
867 double min_diagonal;
880 bool nonzero_starting;
881 };
882
888
895 void
896 initialize(const SparseMatrix<Number, MemorySpace> &A,
897 const AdditionalData &additional_data = AdditionalData());
898 };
899
900
901
919 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
920 class PreconditionILU : public PreconditionIfpackBase<Number, MemorySpace>
921 {
922 public:
927 struct AdditionalData
928 {
942 AdditionalData(const int ilu_fill = 0,
943 const double ilu_atol = 0.,
944 const double ilu_rtol = 1.,
945 const int overlap = 0);
946
954 int ilu_fill;
959 double ilu_atol;
964 double ilu_rtol;
974 int overlap;
975 };
976
982
989 void
990 initialize(const SparseMatrix<Number, MemorySpace> &A,
991 const AdditionalData &additional_data = AdditionalData());
992 };
993
994
995
1003 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
1004 class PreconditionILUT : public PreconditionIfpackBase<Number, MemorySpace>
1005 {
1006 public:
1011 struct AdditionalData
1012 {
1027 AdditionalData(const double ilut_drop = 0.,
1028 const double ilut_fill = 0.,
1029 const double ilut_atol = 0.,
1030 const double ilut_rtol = 1.,
1031 const int overlap = 0);
1032
1039 double ilut_drop;
1049 double ilut_fill;
1054 double ilut_atol;
1059 double ilut_rtol;
1069 int overlap;
1070 };
1071
1072
1073
1078 PreconditionILUT();
1079
1086 void
1087 initialize(const SparseMatrix<Number, MemorySpace> &A,
1088 const AdditionalData &additional_data = AdditionalData());
1089 };
1090
1091
1092
1106 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
1108 : public PreconditionIfpackBase<Number, MemorySpace>
1109 {
1110 public:
1115 struct AdditionalData
1116 {
1128 AdditionalData(const int n_local_parts = 1,
1129 const double omega = 1.,
1130 const int block_overlap = 0,
1131 const int n_sweeps = 1);
1132
1139 int n_local_parts;
1144 double omega;
1149 int block_overlap;
1154 int n_sweeps;
1155 };
1156
1157
1163
1170 void
1171 initialize(const SparseMatrix<Number, MemorySpace> &A,
1172 const AdditionalData &additional_data = AdditionalData());
1173 };
1174
1175
1176
1186 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
1188 : public PreconditionIfpackBase<Number, MemorySpace>
1189 {
1190 public:
1195 struct AdditionalData
1196 {
1207 AdditionalData(const int n_local_parts = 1,
1208 const double omega = 1,
1209 const int overlap = 0,
1210 const int n_sweeps = 1);
1211
1218 int n_local_parts;
1223 double omega;
1233 int overlap;
1238 int n_sweeps;
1239 };
1240
1242
1243 void
1244 initialize(const SparseMatrix<Number, MemorySpace> &A,
1245 const AdditionalData &additional_data = AdditionalData());
1246 };
1247
1248
1249
1259 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
1261 : public PreconditionIfpackBase<Number, MemorySpace>
1262 {
1263 public:
1268 struct AdditionalData
1269 {
1280 AdditionalData(const int n_local_parts = 1,
1281 const double omega = 1,
1282 const int overlap = 1,
1283 const int n_sweeps = 1);
1290 int n_local_parts;
1295 double omega;
1305 int overlap;
1310 int n_sweeps;
1311 };
1312
1314
1315 void
1316 initialize(const SparseMatrix<Number, MemorySpace> &A,
1317 const AdditionalData &additional_data = AdditionalData());
1318 };
1319# endif // DEAL_II_TRILINOS_WITH_IFPACK2
1320
1321 } // namespace TpetraWrappers
1322} // namespace LinearAlgebra
1323
1324
1326
1327#endif // DEAL_II_TRILINOS_WITH_TPETRA
1328
1329#endif
virtual void Tvmult(Vector< Number, MemorySpace > &dst, const Vector< Number, MemorySpace > &src) const
Apply the transpose preconditioner.
Teuchos::ParameterList parameter_list
The list of preconditioner parameters.
void clear()
Desctructor. Destroys the preconditioner, leaving an object like just after having called the constru...
virtual void vmult(::Vector< Number > &dst, ::Vector< Number > &src) const
Apply the preconditioner.
Teuchos::RCP< TpetraTypes::LinearOperator< Number, MemorySpace > > preconditioner
Teuchos::RCP< TpetraTypes::LinearOperator< Number, MemorySpace > > trilinos_rcp() const
Access to underlying Trilinos data.
virtual void Tvmult(::Vector< Number > &dst, ::Vector< Number > &src) const
Apply the transpose preconditioner.
PreconditionBase()=default
Constructor. Does not do anything. The initialize function of the derived classes will have to create...
const TpetraTypes::LinearOperator< Number, MemorySpace > & trilinos_operator() const
Access to underlying Trilinos data.
virtual void vmult(Vector< Number, MemorySpace > &dst, const Vector< Number, MemorySpace > &src) const
Apply the preconditioner.
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:503
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:504
static ::ExceptionBase & ExcNonMatchingMaps(std::string arg1)
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:494
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:516
static ::ExceptionBase & ExcTransposeNotSupported()
PETScWrappers::PreconditionILU PreconditionILU
Tpetra::Operator< Number, LO, GO, NodeType< MemorySpace > > LinearOperator
Standardized data struct to pipe additional flags to the preconditioner.