Reference documentation for deal.II version Git ac8d010384 2020-11-27 19:49:05 +0100
\(\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 - 2020 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 # ifdef DEAL_II_WITH_PETSC
23 
24 # include <deal.II/lac/exceptions.h>
25 
26 # include <petscpc.h>
27 
29 
30 
31 
32 namespace PETScWrappers
33 {
34  // forward declarations
35 # ifndef DOXYGEN
36  class MatrixBase;
37  class VectorBase;
38  class SolverBase;
39 # endif
40 
57  {
58  public:
63 
67  virtual ~PreconditionerBase();
68 
73  void
74  clear();
75 
79  void
80  vmult(VectorBase &dst, const VectorBase &src) const;
81 
85  void
86  Tvmult(VectorBase &dst, const VectorBase &src) const;
87 
88 
92  const PC &
93  get_pc() const;
94 
95  protected:
99  PC pc;
100 
104  Mat matrix;
105 
110  void
111  create_pc();
112 
118  operator Mat() const;
119 
120  // Make the solver class a friend, since it needs to call the conversion
121  // operator.
122  friend class SolverBase;
123  };
124 
125 
126 
138  {
139  public:
145  {};
146 
151  PreconditionJacobi() = default;
152 
153 
159  const MatrixBase & matrix,
160  const AdditionalData &additional_data = AdditionalData());
161 
167  const MPI_Comm & communicator,
168  const AdditionalData &additional_data = AdditionalData());
169 
176  void
177  initialize(const MatrixBase & matrix,
178  const AdditionalData &additional_data = AdditionalData());
179 
180  protected:
185 
191  void
192  initialize();
193  };
194 
195 
196 
221  {
222  public:
228  {};
229 
234  PreconditionBlockJacobi() = default;
235 
241  const MatrixBase & matrix,
242  const AdditionalData &additional_data = AdditionalData());
243 
249  const MPI_Comm & communicator,
250  const AdditionalData &additional_data = AdditionalData());
251 
252 
259  void
260  initialize(const MatrixBase & matrix,
261  const AdditionalData &additional_data = AdditionalData());
262 
263  protected:
268 
274  void
275  initialize();
276  };
277 
278 
279 
289  {
290  public:
296  {
300  AdditionalData(const double omega = 1);
301 
305  double omega;
306  };
307 
312  PreconditionSOR() = default;
313 
319  const AdditionalData &additional_data = AdditionalData());
320 
327  void
328  initialize(const MatrixBase & matrix,
329  const AdditionalData &additional_data = AdditionalData());
330 
331  protected:
336  };
337 
338 
339 
349  {
350  public:
356  {
360  AdditionalData(const double omega = 1);
361 
365  double omega;
366  };
367 
372  PreconditionSSOR() = default;
373 
378  PreconditionSSOR(const MatrixBase & matrix,
379  const AdditionalData &additional_data = AdditionalData());
380 
387  void
388  initialize(const MatrixBase & matrix,
389  const AdditionalData &additional_data = AdditionalData());
390 
391  protected:
396  };
397 
398 
399 
412  {
413  public:
419  {
423  AdditionalData(const double omega = 1);
424 
428  double omega;
429  };
430 
435  PreconditionEisenstat() = default;
436 
442  const MatrixBase & matrix,
443  const AdditionalData &additional_data = AdditionalData());
444 
451  void
452  initialize(const MatrixBase & matrix,
453  const AdditionalData &additional_data = AdditionalData());
454 
455  protected:
460  };
461 
462 
463 
473  {
474  public:
480  {
484  AdditionalData(const unsigned int levels = 0);
485 
489  unsigned int levels;
490  };
491 
496  PreconditionICC() = default;
497 
502  PreconditionICC(const MatrixBase & matrix,
503  const AdditionalData &additional_data = AdditionalData());
504 
511  void
512  initialize(const MatrixBase & matrix,
513  const AdditionalData &additional_data = AdditionalData());
514 
515  protected:
520  };
521 
522 
523 
533  {
534  public:
540  {
544  AdditionalData(const unsigned int levels = 0);
545 
549  unsigned int levels;
550  };
551 
556  PreconditionILU() = default;
557 
562  PreconditionILU(const MatrixBase & matrix,
563  const AdditionalData &additional_data = AdditionalData());
564 
571  void
572  initialize(const MatrixBase & matrix,
573  const AdditionalData &additional_data = AdditionalData());
574 
575  protected:
580  };
581 
582 
583 
598  {
599  public:
605  {
610  AdditionalData(const double pivoting = 1.e-6,
611  const double zero_pivot = 1.e-12,
612  const double damping = 0.0);
613 
619  double pivoting;
620 
625  double zero_pivot;
626 
631  double damping;
632  };
633 
638  PreconditionLU() = default;
639 
644  PreconditionLU(const MatrixBase & matrix,
645  const AdditionalData &additional_data = AdditionalData());
646 
653  void
654  initialize(const MatrixBase & matrix,
655  const AdditionalData &additional_data = AdditionalData());
656 
657  protected:
662  };
663 
664 
665 
677  {
678  public:
684  {
689  AdditionalData(const bool symmetric_operator = false,
690  const double strong_threshold = 0.25,
691  const double max_row_sum = 0.9,
692  const unsigned int aggressive_coarsening_num_levels = 0,
693  const bool output_details = false);
694 
702 
709 
721  double max_row_sum;
722 
729 
735  };
736 
741  PreconditionBoomerAMG() = default;
742 
748  const MatrixBase & matrix,
749  const AdditionalData &additional_data = AdditionalData());
750 
756  const MPI_Comm & communicator,
757  const AdditionalData &additional_data = AdditionalData());
758 
759 
766  void
767  initialize(const MatrixBase & matrix,
768  const AdditionalData &additional_data = AdditionalData());
769 
770  protected:
775 
781  void
782  initialize();
783  };
784 
785 
786 
807  {
808  public:
814  {
818  AdditionalData(const unsigned int symmetric = 1,
819  const unsigned int n_levels = 1,
820  const double threshold = 0.1,
821  const double filter = 0.05,
822  const bool output_details = false);
823 
835  unsigned int symmetric;
836 
843  unsigned int n_levels;
844 
855  double threshold;
856 
867  double filter;
868 
874  };
875 
876 
877 
882  PreconditionParaSails() = default;
883 
889  const MatrixBase & matrix,
890  const AdditionalData &additional_data = AdditionalData());
891 
898  void
899  initialize(const MatrixBase & matrix,
900  const AdditionalData &additional_data = AdditionalData());
901 
902  private:
907  };
908 
909 
910 
917  {
918  public:
924  {};
925 
930  PreconditionNone() = default;
931 
937  PreconditionNone(const MatrixBase & matrix,
938  const AdditionalData &additional_data = AdditionalData());
939 
947  void
948  initialize(const MatrixBase & matrix,
949  const AdditionalData &additional_data = AdditionalData());
950 
951  private:
956  };
957 } // namespace PETScWrappers
958 
959 
960 
962 
963 
964 # endif // DEAL_II_WITH_PETSC
965 
966 #endif
967 /*--------------------------- petsc_precondition.h --------------------------*/
Matrix is symmetric.
void vmult(VectorBase &dst, const VectorBase &src) const
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
PETScWrappers::PreconditionILU PreconditionILU
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:372
PETScWrappers::PreconditionJacobi PreconditionJacobi
void Tvmult(VectorBase &dst, const VectorBase &src) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:371
PreconditionSSOR< SparseMatrix > PreconditionSSOR