Reference documentation for deal.II version GIT 5c75851574 2022-09-25 16:30:01+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\}}\)
petsc_precondition.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2021 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_petsc_precondition_h
17 #define dealii_petsc_precondition_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/point.h>
24 
25 #ifdef DEAL_II_WITH_PETSC
26 
27 # include <deal.II/lac/exceptions.h>
28 
29 # include <petscpc.h>
30 
32 
33 
34 
35 namespace PETScWrappers
36 {
37  // forward declarations
38 # ifndef DOXYGEN
39  class MatrixBase;
40  class VectorBase;
41  class SolverBase;
42 # endif
43 
60  {
61  public:
66 
70  virtual ~PreconditionBase();
71 
76  void
77  clear();
78 
82  void
83  vmult(VectorBase &dst, const VectorBase &src) const;
84 
88  void
89  Tvmult(VectorBase &dst, const VectorBase &src) const;
90 
91 
95  const PC &
96  get_pc() const;
97 
98  protected:
102  PC pc;
103 
107  Mat matrix;
108 
113  void
114  create_pc();
115 
121  operator Mat() const;
122 
123  // Make the solver class a friend, since it needs to call the conversion
124  // operator.
125  friend class SolverBase;
126  };
127 
128 
129 
141  {
142  public:
148  {};
149 
154  PreconditionJacobi() = default;
155 
156 
162  const MatrixBase & matrix,
164 
170  const MPI_Comm & communicator,
172 
179  void
180  initialize(const MatrixBase & matrix,
182 
183  protected:
188 
194  void
195  initialize();
196  };
197 
198 
199 
224  {
225  public:
231  {};
232 
238 
244  const MatrixBase & matrix,
246 
252  const MPI_Comm & communicator,
254 
255 
262  void
263  initialize(const MatrixBase & matrix,
265 
266  protected:
271 
277  void
278  initialize();
279  };
280 
281 
282 
292  {
293  public:
299  {
303  AdditionalData(const double omega = 1);
304 
308  double omega;
309  };
310 
315  PreconditionSOR() = default;
316 
323 
330  void
331  initialize(const MatrixBase & matrix,
333 
334  protected:
339  };
340 
341 
342 
352  {
353  public:
359  {
363  AdditionalData(const double omega = 1);
364 
368  double omega;
369  };
370 
375  PreconditionSSOR() = default;
376 
383 
390  void
391  initialize(const MatrixBase & matrix,
393 
394  protected:
399  };
400 
410  {
411  public:
417  {
421  AdditionalData(const unsigned int levels = 0);
422 
426  unsigned int levels;
427  };
428 
433  PreconditionICC() = default;
434 
441 
448  void
449  initialize(const MatrixBase & matrix,
451 
452  protected:
457  };
458 
459 
460 
470  {
471  public:
477  {
481  AdditionalData(const unsigned int levels = 0);
482 
486  unsigned int levels;
487  };
488 
493  PreconditionILU() = default;
494 
501 
508  void
509  initialize(const MatrixBase & matrix,
511 
512  protected:
517  };
518 
519 
520 
535  {
536  public:
542  {
547  AdditionalData(const double pivoting = 1.e-6,
548  const double zero_pivot = 1.e-12,
549  const double damping = 0.0);
550 
556  double pivoting;
557 
562  double zero_pivot;
563 
568  double damping;
569  };
570 
575  PreconditionLU() = default;
576 
583 
590  void
591  initialize(const MatrixBase & matrix,
593 
594  protected:
599  };
600 
601 
602 
614  {
615  public:
621  {
625  enum class RelaxationType
626  {
627  Jacobi,
630  SORJacobi,
637  CG,
638  Chebyshev,
639  FCFJacobi,
641  None
642  };
643 
649  const bool symmetric_operator = false,
650  const double strong_threshold = 0.25,
651  const double max_row_sum = 0.9,
652  const unsigned int aggressive_coarsening_num_levels = 0,
653  const bool output_details = false,
658  const unsigned int n_sweeps_coarse = 1,
659  const double tol = 0.0,
660  const unsigned int max_iter = 1,
661  const bool w_cycle = false);
662 
669 
676 
688  double max_row_sum;
689 
696 
702 
707 
712 
717 
721  unsigned int n_sweeps_coarse;
722 
726  double tol;
727 
731  unsigned int max_iter;
732 
737  bool w_cycle;
738  };
739 
745 
751  const MatrixBase & matrix,
753 
759  const MPI_Comm & communicator,
761 
762 
769  void
770  initialize(const MatrixBase & matrix,
772 
773  protected:
778 
784  void
785  initialize();
786  };
787 
788 
789 
810  {
811  public:
817  {
821  AdditionalData(const unsigned int symmetric = 1,
822  const unsigned int n_levels = 1,
823  const double threshold = 0.1,
824  const double filter = 0.05,
825  const bool output_details = false);
826 
838  unsigned int symmetric;
839 
846  unsigned int n_levels;
847 
858  double threshold;
859 
870  double filter;
871 
877  };
878 
879 
880 
886 
892  const MatrixBase & matrix,
894 
901  void
902  initialize(const MatrixBase & matrix,
904 
905  private:
910  };
911 
912 
913 
920  {
921  public:
927  {};
928 
933  PreconditionNone() = default;
934 
942 
950  void
951  initialize(const MatrixBase & matrix,
953 
954  private:
959  };
960 
978  template <int dim>
980  {
981  public:
987  {
992  AdditionalData(const bool use_vertices = true,
993  const bool use_edges = false,
994  const bool use_faces = false,
995  const bool symmetric = false,
996  const std::vector<Point<dim>> coords = {});
997 
1003 
1010 
1017 
1022 
1027  std::vector<Point<dim>> coords;
1028  };
1029 
1034  PreconditionBDDC() = default;
1035 
1042 
1047  PreconditionBDDC(const MPI_Comm communicator,
1049 
1056  void
1057  initialize(const MatrixBase & matrix,
1059 
1060  protected:
1065 
1071  void
1072  initialize();
1073  };
1074 
1080 } // namespace PETScWrappers
1081 
1082 
1084 
1085 
1086 #endif // DEAL_II_WITH_PETSC
1087 
1088 #endif
void Tvmult(VectorBase &dst, const VectorBase &src) const
void vmult(VectorBase &dst, const VectorBase &src) const
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
Definition: point.h:111
#define DEAL_II_DEPRECATED
Definition: config.h:164
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:457
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:458
AdditionalData(const bool use_vertices=true, const bool use_edges=false, const bool use_faces=false, const bool symmetric=false, const std::vector< Point< dim >> coords={})
AdditionalData(const bool symmetric_operator=false, const double strong_threshold=0.25, const double max_row_sum=0.9, const unsigned int aggressive_coarsening_num_levels=0, const bool output_details=false, const RelaxationType relaxation_type_up=RelaxationType::SORJacobi, const RelaxationType relaxation_type_down=RelaxationType::SORJacobi, const RelaxationType relaxation_type_coarse=RelaxationType::GaussianElimination, const unsigned int n_sweeps_coarse=1, const double tol=0.0, const unsigned int max_iter=1, const bool w_cycle=false)
AdditionalData(const double pivoting=1.e-6, const double zero_pivot=1.e-12, const double damping=0.0)
AdditionalData(const unsigned int symmetric=1, const unsigned int n_levels=1, const double threshold=0.1, const double filter=0.05, const bool output_details=false)