Reference documentation for deal.II version Git 497f915867 2021-09-17 22:46:48 +0200
\(\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 
23 
24 # ifdef DEAL_II_WITH_PETSC
25 
26 # include <deal.II/lac/exceptions.h>
27 
28 # include <petscpc.h>
29 
31 
32 
33 
34 namespace PETScWrappers
35 {
36  // forward declarations
37 # ifndef DOXYGEN
38  class MatrixBase;
39  class VectorBase;
40  class SolverBase;
41 # endif
42 
59  {
60  public:
65 
69  virtual ~PreconditionBase();
70 
75  void
76  clear();
77 
81  void
82  vmult(VectorBase &dst, const VectorBase &src) const;
83 
87  void
88  Tvmult(VectorBase &dst, const VectorBase &src) const;
89 
90 
94  const PC &
95  get_pc() const;
96 
97  protected:
101  PC pc;
102 
106  Mat matrix;
107 
112  void
113  create_pc();
114 
120  operator Mat() const;
121 
122  // Make the solver class a friend, since it needs to call the conversion
123  // operator.
124  friend class SolverBase;
125  };
126 
127 
128 
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 
223  {
224  public:
230  {};
231 
236  PreconditionBlockJacobi() = default;
237 
243  const MatrixBase & matrix,
244  const AdditionalData &additional_data = AdditionalData());
245 
251  const MPI_Comm & communicator,
252  const AdditionalData &additional_data = AdditionalData());
253 
254 
261  void
262  initialize(const MatrixBase & matrix,
263  const AdditionalData &additional_data = AdditionalData());
264 
265  protected:
270 
276  void
277  initialize();
278  };
279 
280 
281 
291  {
292  public:
298  {
302  AdditionalData(const double omega = 1);
303 
307  double omega;
308  };
309 
314  PreconditionSOR() = default;
315 
321  const AdditionalData &additional_data = AdditionalData());
322 
329  void
330  initialize(const MatrixBase & matrix,
331  const AdditionalData &additional_data = AdditionalData());
332 
333  protected:
338  };
339 
340 
341 
351  {
352  public:
358  {
362  AdditionalData(const double omega = 1);
363 
367  double omega;
368  };
369 
374  PreconditionSSOR() = default;
375 
380  PreconditionSSOR(const MatrixBase & matrix,
381  const AdditionalData &additional_data = AdditionalData());
382 
389  void
390  initialize(const MatrixBase & matrix,
391  const AdditionalData &additional_data = AdditionalData());
392 
393  protected:
398  };
399 
409  {
410  public:
416  {
420  AdditionalData(const unsigned int levels = 0);
421 
425  unsigned int levels;
426  };
427 
432  PreconditionICC() = default;
433 
438  PreconditionICC(const MatrixBase & matrix,
439  const AdditionalData &additional_data = AdditionalData());
440 
447  void
448  initialize(const MatrixBase & matrix,
449  const AdditionalData &additional_data = AdditionalData());
450 
451  protected:
456  };
457 
458 
459 
469  {
470  public:
476  {
480  AdditionalData(const unsigned int levels = 0);
481 
485  unsigned int levels;
486  };
487 
492  PreconditionILU() = default;
493 
498  PreconditionILU(const MatrixBase & matrix,
499  const AdditionalData &additional_data = AdditionalData());
500 
507  void
508  initialize(const MatrixBase & matrix,
509  const AdditionalData &additional_data = AdditionalData());
510 
511  protected:
516  };
517 
518 
519 
534  {
535  public:
541  {
546  AdditionalData(const double pivoting = 1.e-6,
547  const double zero_pivot = 1.e-12,
548  const double damping = 0.0);
549 
555  double pivoting;
556 
561  double zero_pivot;
562 
567  double damping;
568  };
569 
574  PreconditionLU() = default;
575 
580  PreconditionLU(const MatrixBase & matrix,
581  const AdditionalData &additional_data = AdditionalData());
582 
589  void
590  initialize(const MatrixBase & matrix,
591  const AdditionalData &additional_data = AdditionalData());
592 
593  protected:
598  };
599 
600 
601 
613  {
614  public:
620  {
624  enum class RelaxationType
625  {
626  Jacobi,
627  sequentialGaussSeidel,
628  seqboundaryGaussSeidel,
629  SORJacobi,
630  backwardSORJacobi,
631  symmetricSORJacobi,
632  l1scaledSORJacobi,
633  GaussianElimination,
634  l1GaussSeidel,
635  backwardl1GaussSeidel,
636  CG,
637  Chebyshev,
638  FCFJacobi,
639  l1scaledJacobi,
640  None
641  };
642 
648  const bool symmetric_operator = false,
649  const double strong_threshold = 0.25,
650  const double max_row_sum = 0.9,
651  const unsigned int aggressive_coarsening_num_levels = 0,
652  const bool output_details = false,
653  const RelaxationType relaxation_type_up = RelaxationType::SORJacobi,
654  const RelaxationType relaxation_type_down = RelaxationType::SORJacobi,
655  const RelaxationType relaxation_type_coarse =
656  RelaxationType::GaussianElimination,
657  const unsigned int n_sweeps_coarse = 1,
658  const double tol = 0.0,
659  const unsigned int max_iter = 1,
660  const bool w_cycle = false);
661 
668 
675 
687  double max_row_sum;
688 
695 
701 
706 
711 
716 
720  unsigned int n_sweeps_coarse;
721 
725  double tol;
726 
730  unsigned int max_iter;
731 
736  bool w_cycle;
737  };
738 
743  PreconditionBoomerAMG() = default;
744 
750  const MatrixBase & matrix,
751  const AdditionalData &additional_data = AdditionalData());
752 
758  const MPI_Comm & communicator,
759  const AdditionalData &additional_data = AdditionalData());
760 
761 
768  void
769  initialize(const MatrixBase & matrix,
770  const AdditionalData &additional_data = AdditionalData());
771 
772  protected:
777 
783  void
784  initialize();
785  };
786 
787 
788 
809  {
810  public:
816  {
820  AdditionalData(const unsigned int symmetric = 1,
821  const unsigned int n_levels = 1,
822  const double threshold = 0.1,
823  const double filter = 0.05,
824  const bool output_details = false);
825 
837  unsigned int symmetric;
838 
845  unsigned int n_levels;
846 
857  double threshold;
858 
869  double filter;
870 
876  };
877 
878 
879 
884  PreconditionParaSails() = default;
885 
891  const MatrixBase & matrix,
892  const AdditionalData &additional_data = AdditionalData());
893 
900  void
901  initialize(const MatrixBase & matrix,
902  const AdditionalData &additional_data = AdditionalData());
903 
904  private:
909  };
910 
911 
912 
919  {
920  public:
926  {};
927 
932  PreconditionNone() = default;
933 
939  PreconditionNone(const MatrixBase & matrix,
940  const AdditionalData &additional_data = AdditionalData());
941 
949  void
950  initialize(const MatrixBase & matrix,
951  const AdditionalData &additional_data = AdditionalData());
952 
953  private:
958  };
959 
965 } // namespace PETScWrappers
966 
967 
968 
970 
971 
972 # endif // DEAL_II_WITH_PETSC
973 
974 #endif
975 /*--------------------------- petsc_precondition.h --------------------------*/
Matrix is symmetric.
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
PETScWrappers::PreconditionILU PreconditionILU
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
PETScWrappers::PreconditionJacobi PreconditionJacobi
void vmult(VectorBase &dst, const VectorBase &src) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
void Tvmult(VectorBase &dst, const VectorBase &src) const
PreconditionSSOR< SparseMatrix > PreconditionSSOR
#define DEAL_II_DEPRECATED
Definition: config.h:160