Reference documentation for deal.II version Git 7026f387cc 2019-10-15 14:19:01 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
petsc_precondition.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2018 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 
28 DEAL_II_NAMESPACE_OPEN
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 
83 
87  const PC &
88  get_pc() const;
89 
90  protected:
94  PC pc;
95 
99  Mat matrix;
100 
105  void
106  create_pc();
107 
113  operator Mat() const;
114 
115  // Make the solver class a friend, since it needs to call the conversion
116  // operator.
117  friend class SolverBase;
118  };
119 
120 
121 
134  {
135  public:
141  {};
142 
147  PreconditionJacobi() = default;
148 
149 
155  const MatrixBase & matrix,
156  const AdditionalData &additional_data = AdditionalData());
157 
163  const MPI_Comm communicator,
164  const AdditionalData &additional_data = AdditionalData());
165 
172  void
173  initialize(const MatrixBase & matrix,
174  const AdditionalData &additional_data = AdditionalData());
175 
176  protected:
181 
187  void
188  initialize();
189  };
190 
191 
192 
218  {
219  public:
225  {};
226 
231  PreconditionBlockJacobi() = default;
232 
238  const MatrixBase & matrix,
239  const AdditionalData &additional_data = AdditionalData());
240 
246  const MPI_Comm communicator,
247  const AdditionalData &additional_data = AdditionalData());
248 
249 
256  void
257  initialize(const MatrixBase & matrix,
258  const AdditionalData &additional_data = AdditionalData());
259 
260  protected:
265 
271  void
272  initialize();
273  };
274 
275 
276 
287  {
288  public:
294  {
298  AdditionalData(const double omega = 1);
299 
303  double omega;
304  };
305 
310  PreconditionSOR() = default;
311 
317  const AdditionalData &additional_data = AdditionalData());
318 
325  void
326  initialize(const MatrixBase & matrix,
327  const AdditionalData &additional_data = AdditionalData());
328 
329  protected:
334  };
335 
336 
337 
348  {
349  public:
355  {
359  AdditionalData(const double omega = 1);
360 
364  double omega;
365  };
366 
371  PreconditionSSOR() = default;
372 
377  PreconditionSSOR(const MatrixBase & matrix,
378  const AdditionalData &additional_data = AdditionalData());
379 
386  void
387  initialize(const MatrixBase & matrix,
388  const AdditionalData &additional_data = AdditionalData());
389 
390  protected:
395  };
396 
397 
398 
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 
474  {
475  public:
481  {
485  AdditionalData(const unsigned int levels = 0);
486 
490  unsigned int levels;
491  };
492 
497  PreconditionICC() = default;
498 
503  PreconditionICC(const MatrixBase & matrix,
504  const AdditionalData &additional_data = AdditionalData());
505 
512  void
513  initialize(const MatrixBase & matrix,
514  const AdditionalData &additional_data = AdditionalData());
515 
516  protected:
521  };
522 
523 
524 
535  {
536  public:
542  {
546  AdditionalData(const unsigned int levels = 0);
547 
551  unsigned int levels;
552  };
553 
558  PreconditionILU() = default;
559 
564  PreconditionILU(const MatrixBase & matrix,
565  const AdditionalData &additional_data = AdditionalData());
566 
573  void
574  initialize(const MatrixBase & matrix,
575  const AdditionalData &additional_data = AdditionalData());
576 
577  protected:
582  };
583 
584 
585 
601  {
602  public:
608  {
613  AdditionalData(const double pivoting = 1.e-6,
614  const double zero_pivot = 1.e-12,
615  const double damping = 0.0);
616 
622  double pivoting;
623 
628  double zero_pivot;
629 
634  double damping;
635  };
636 
641  PreconditionLU() = default;
642 
647  PreconditionLU(const MatrixBase & matrix,
648  const AdditionalData &additional_data = AdditionalData());
649 
656  void
657  initialize(const MatrixBase & matrix,
658  const AdditionalData &additional_data = AdditionalData());
659 
660  protected:
665  };
666 
667 
668 
681  {
682  public:
688  {
693  AdditionalData(const bool symmetric_operator = false,
694  const double strong_threshold = 0.25,
695  const double max_row_sum = 0.9,
696  const unsigned int aggressive_coarsening_num_levels = 0,
697  const bool output_details = false);
698 
706 
713 
725  double max_row_sum;
726 
733 
739  };
740 
745  PreconditionBoomerAMG() = default;
746 
752  const MatrixBase & matrix,
753  const AdditionalData &additional_data = AdditionalData());
754 
760  const MPI_Comm communicator,
761  const AdditionalData &additional_data = AdditionalData());
762 
763 
770  void
771  initialize(const MatrixBase & matrix,
772  const AdditionalData &additional_data = AdditionalData());
773 
774  protected:
779 
785  void
786  initialize();
787  };
788 
789 
790 
812  {
813  public:
819  {
823  AdditionalData(const unsigned int symmetric = 1,
824  const unsigned int n_levels = 1,
825  const double threshold = 0.1,
826  const double filter = 0.05,
827  const bool output_details = false);
828 
840  unsigned int symmetric;
841 
848  unsigned int n_levels;
849 
860  double threshold;
861 
872  double filter;
873 
879  };
880 
881 
882 
887  PreconditionParaSails() = default;
888 
894  const MatrixBase & matrix,
895  const AdditionalData &additional_data = AdditionalData());
896 
903  void
904  initialize(const MatrixBase & matrix,
905  const AdditionalData &additional_data = AdditionalData());
906 
907  private:
912  };
913 
914 
915 
923  {
924  public:
930  {};
931 
936  PreconditionNone() = default;
937 
943  PreconditionNone(const MatrixBase & matrix,
944  const AdditionalData &additional_data = AdditionalData());
945 
953  void
954  initialize(const MatrixBase & matrix,
955  const AdditionalData &additional_data = AdditionalData());
956 
957  private:
962  };
963 } // namespace PETScWrappers
964 
965 
966 
967 DEAL_II_NAMESPACE_CLOSE
968 
969 
970 # endif // DEAL_II_WITH_PETSC
971 
972 #endif
973 /*--------------------------- petsc_precondition.h --------------------------*/
void vmult(VectorBase &dst, const VectorBase &src) const