Reference documentation for deal.II version 9.2.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\}}\)
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 
58  {
59  public:
64 
68  virtual ~PreconditionerBase();
69 
74  void
75  clear();
76 
80  void
81  vmult(VectorBase &dst, const VectorBase &src) const;
82 
86  void
87  Tvmult(VectorBase &dst, const VectorBase &src) const;
88 
89 
93  const PC &
94  get_pc() const;
95 
96  protected:
100  PC pc;
101 
105  Mat matrix;
106 
111  void
112  create_pc();
113 
119  operator Mat() const;
120 
121  // Make the solver class a friend, since it needs to call the conversion
122  // operator.
123  friend class SolverBase;
124  };
125 
126 
127 
140  {
141  public:
147  {};
148 
153  PreconditionJacobi() = default;
154 
155 
161  const MatrixBase & matrix,
162  const AdditionalData &additional_data = AdditionalData());
163 
169  const MPI_Comm communicator,
170  const AdditionalData &additional_data = AdditionalData());
171 
178  void
179  initialize(const MatrixBase & matrix,
180  const AdditionalData &additional_data = AdditionalData());
181 
182  protected:
187 
193  void
194  initialize();
195  };
196 
197 
198 
224  {
225  public:
231  {};
232 
237  PreconditionBlockJacobi() = default;
238 
244  const MatrixBase & matrix,
245  const AdditionalData &additional_data = AdditionalData());
246 
252  const MPI_Comm communicator,
253  const AdditionalData &additional_data = AdditionalData());
254 
255 
262  void
263  initialize(const MatrixBase & matrix,
264  const AdditionalData &additional_data = AdditionalData());
265 
266  protected:
271 
277  void
278  initialize();
279  };
280 
281 
282 
293  {
294  public:
300  {
304  AdditionalData(const double omega = 1);
305 
309  double omega;
310  };
311 
316  PreconditionSOR() = default;
317 
323  const AdditionalData &additional_data = AdditionalData());
324 
331  void
332  initialize(const MatrixBase & matrix,
333  const AdditionalData &additional_data = AdditionalData());
334 
335  protected:
340  };
341 
342 
343 
354  {
355  public:
361  {
365  AdditionalData(const double omega = 1);
366 
370  double omega;
371  };
372 
377  PreconditionSSOR() = default;
378 
383  PreconditionSSOR(const MatrixBase & matrix,
384  const AdditionalData &additional_data = AdditionalData());
385 
392  void
393  initialize(const MatrixBase & matrix,
394  const AdditionalData &additional_data = AdditionalData());
395 
396  protected:
401  };
402 
403 
404 
418  {
419  public:
425  {
429  AdditionalData(const double omega = 1);
430 
434  double omega;
435  };
436 
441  PreconditionEisenstat() = default;
442 
448  const MatrixBase & matrix,
449  const AdditionalData &additional_data = AdditionalData());
450 
457  void
458  initialize(const MatrixBase & matrix,
459  const AdditionalData &additional_data = AdditionalData());
460 
461  protected:
466  };
467 
468 
469 
480  {
481  public:
487  {
491  AdditionalData(const unsigned int levels = 0);
492 
496  unsigned int levels;
497  };
498 
503  PreconditionICC() = default;
504 
509  PreconditionICC(const MatrixBase & matrix,
510  const AdditionalData &additional_data = AdditionalData());
511 
518  void
519  initialize(const MatrixBase & matrix,
520  const AdditionalData &additional_data = AdditionalData());
521 
522  protected:
527  };
528 
529 
530 
541  {
542  public:
548  {
552  AdditionalData(const unsigned int levels = 0);
553 
557  unsigned int levels;
558  };
559 
564  PreconditionILU() = default;
565 
570  PreconditionILU(const MatrixBase & matrix,
571  const AdditionalData &additional_data = AdditionalData());
572 
579  void
580  initialize(const MatrixBase & matrix,
581  const AdditionalData &additional_data = AdditionalData());
582 
583  protected:
588  };
589 
590 
591 
607  {
608  public:
614  {
619  AdditionalData(const double pivoting = 1.e-6,
620  const double zero_pivot = 1.e-12,
621  const double damping = 0.0);
622 
628  double pivoting;
629 
634  double zero_pivot;
635 
640  double damping;
641  };
642 
647  PreconditionLU() = default;
648 
653  PreconditionLU(const MatrixBase & matrix,
654  const AdditionalData &additional_data = AdditionalData());
655 
662  void
663  initialize(const MatrixBase & matrix,
664  const AdditionalData &additional_data = AdditionalData());
665 
666  protected:
671  };
672 
673 
674 
687  {
688  public:
694  {
699  AdditionalData(const bool symmetric_operator = false,
700  const double strong_threshold = 0.25,
701  const double max_row_sum = 0.9,
702  const unsigned int aggressive_coarsening_num_levels = 0,
703  const bool output_details = false);
704 
712 
719 
731  double max_row_sum;
732 
739 
745  };
746 
751  PreconditionBoomerAMG() = default;
752 
758  const MatrixBase & matrix,
759  const AdditionalData &additional_data = AdditionalData());
760 
766  const MPI_Comm communicator,
767  const AdditionalData &additional_data = AdditionalData());
768 
769 
776  void
777  initialize(const MatrixBase & matrix,
778  const AdditionalData &additional_data = AdditionalData());
779 
780  protected:
785 
791  void
792  initialize();
793  };
794 
795 
796 
818  {
819  public:
825  {
829  AdditionalData(const unsigned int symmetric = 1,
830  const unsigned int n_levels = 1,
831  const double threshold = 0.1,
832  const double filter = 0.05,
833  const bool output_details = false);
834 
846  unsigned int symmetric;
847 
854  unsigned int n_levels;
855 
866  double threshold;
867 
878  double filter;
879 
885  };
886 
887 
888 
893  PreconditionParaSails() = default;
894 
900  const MatrixBase & matrix,
901  const AdditionalData &additional_data = AdditionalData());
902 
909  void
910  initialize(const MatrixBase & matrix,
911  const AdditionalData &additional_data = AdditionalData());
912 
913  private:
918  };
919 
920 
921 
929  {
930  public:
936  {};
937 
942  PreconditionNone() = default;
943 
949  PreconditionNone(const MatrixBase & matrix,
950  const AdditionalData &additional_data = AdditionalData());
951 
959  void
960  initialize(const MatrixBase & matrix,
961  const AdditionalData &additional_data = AdditionalData());
962 
963  private:
968  };
969 } // namespace PETScWrappers
970 
971 
972 
974 
975 
976 # endif // DEAL_II_WITH_PETSC
977 
978 #endif
979 /*--------------------------- 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:359
PETScWrappers::PreconditionJacobi PreconditionJacobi
void Tvmult(VectorBase &dst, const VectorBase &src) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
PreconditionSSOR< SparseMatrix > PreconditionSSOR