Reference documentation for deal.II version 8.1.0
petsc_precondition.h
1 // ---------------------------------------------------------------------
2 // @f$Id: petsc_precondition.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2004 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__petsc_precondition_h
18 #define __deal2__petsc_precondition_h
19 
20 
21 #include <deal.II/base/config.h>
22 
23 #ifdef DEAL_II_WITH_PETSC
24 
25 # include <deal.II/lac/exceptions.h>
26 # include <petscpc.h>
27 
29 
30 
31 
32 namespace PETScWrappers
33 {
34  // forward declarations
35  class MatrixBase;
36  class VectorBase;
37  class SolverBase;
38 
39 
57  {
58  public:
63 
67  virtual ~PreconditionerBase ();
68 
73  void vmult (VectorBase &dst,
74  const VectorBase &src) const;
75 
76 
81  const PC &get_pc () const;
82 
83  protected:
87  PC pc;
88 
93  Mat matrix;
94 
100  void create_pc ();
101 
110  operator Mat () const;
111 
117  friend class SolverBase;
118  };
119 
120 
121 
133  {
134  public:
141  {};
142 
149 
150 
158  const AdditionalData &additional_data = AdditionalData());
159 
171  void initialize (const MatrixBase &matrix,
172  const AdditionalData &additional_data = AdditionalData());
173 
174  protected:
180  };
181 
182 
183 
202  {
203  public:
210  {};
211 
218 
226  const AdditionalData &additional_data = AdditionalData());
227 
239  void initialize (const MatrixBase &matrix,
240  const AdditionalData &additional_data = AdditionalData());
241 
242  protected:
248  };
249 
250 
251 
263  {
264  public:
271  {
277  AdditionalData (const double omega = 1);
278 
282  double omega;
283  };
284 
290  PreconditionSOR ();
291 
299  const AdditionalData &additional_data = AdditionalData());
300 
312  void initialize (const MatrixBase &matrix,
313  const AdditionalData &additional_data = AdditionalData());
314 
315  protected:
321  };
322 
323 
324 
336  {
337  public:
344  {
350  AdditionalData (const double omega = 1);
351 
355  double omega;
356  };
357 
363  PreconditionSSOR ();
364 
371  PreconditionSSOR (const MatrixBase &matrix,
372  const AdditionalData &additional_data = AdditionalData());
373 
385  void initialize (const MatrixBase &matrix,
386  const AdditionalData &additional_data = AdditionalData());
387 
388  protected:
394  };
395 
396 
397 
409  {
410  public:
417  {
423  AdditionalData (const double omega = 1);
424 
428  double omega;
429  };
430 
437 
444  PreconditionEisenstat (const MatrixBase &matrix,
445  const AdditionalData &additional_data = AdditionalData());
446 
458  void initialize (const MatrixBase &matrix,
459  const AdditionalData &additional_data = AdditionalData());
460 
461  protected:
467  };
468 
469 
470 
482  {
483  public:
490  {
496  AdditionalData (const unsigned int levels = 0);
497 
501  unsigned int levels;
502  };
503 
509  PreconditionICC ();
510 
517  PreconditionICC (const MatrixBase &matrix,
518  const AdditionalData &additional_data = AdditionalData());
519 
531  void initialize (const MatrixBase &matrix,
532  const AdditionalData &additional_data = AdditionalData());
533 
534  protected:
540  };
541 
542 
543 
555  {
556  public:
563  {
569  AdditionalData (const unsigned int levels = 0);
570 
574  unsigned int levels;
575  };
576 
582  PreconditionILU ();
583 
590  PreconditionILU (const MatrixBase &matrix,
591  const AdditionalData &additional_data = AdditionalData());
592 
604  void initialize (const MatrixBase &matrix,
605  const AdditionalData &additional_data = AdditionalData());
606 
607  protected:
613  };
614 
615 
616 
630  {
631  public:
638  {
644  AdditionalData (const double pivoting = 1.e-6,
645  const double zero_pivot = 1.e-12,
646  const double damping = 0.0);
647 
656  double pivoting;
657 
664  double zero_pivot;
665 
671  double damping;
672  };
673 
679  PreconditionLU ();
680 
687  PreconditionLU (const MatrixBase &matrix,
688  const AdditionalData &additional_data = AdditionalData());
689 
701  void initialize (const MatrixBase &matrix,
702  const AdditionalData &additional_data = AdditionalData());
703 
704  protected:
710  };
711 
712 
713 
714 
727  {
728  public:
735  {
742  const bool symmetric_operator = false,
743  const double strong_threshold = 0.25,
744  const double max_row_sum = 0.9,
745  const unsigned int aggressive_coarsening_num_levels = 0,
746  const bool output_details = false
747  );
748 
760 
771 
793  double max_row_sum;
794 
802 
810  };
811 
818 
825  PreconditionBoomerAMG (const MatrixBase &matrix,
826  const AdditionalData &additional_data = AdditionalData());
827 
839  void initialize (const MatrixBase &matrix,
840  const AdditionalData &additional_data = AdditionalData());
841 
842  protected:
848  };
849 
850 
851 
874  {
875  public:
882  {
887  const unsigned int symmetric = 1,
888  const unsigned int n_levels = 1,
889  const double threshold = 0.1,
890  const double filter = 0.05,
891  const bool output_details = false
892  );
893 
904  unsigned int symmetric;
905 
915  unsigned int n_levels;
916 
934  double threshold;
935 
951  double filter;
952 
960  };
961 
962 
963 
970 
977  PreconditionParaSails (const MatrixBase &matrix,
978  const AdditionalData &additional_data = AdditionalData());
979 
991  void initialize (const MatrixBase &matrix,
992  const AdditionalData &additional_data = AdditionalData());
993 
994  private:
1000  };
1001 
1002 
1003 
1011  {
1012  public:
1019  {};
1020 
1026  PreconditionNone ();
1027 
1035  PreconditionNone (const MatrixBase &matrix,
1036  const AdditionalData &additional_data = AdditionalData());
1037 
1050  void initialize (const MatrixBase &matrix,
1051  const AdditionalData &additional_data = AdditionalData());
1052 
1053  private:
1059  };
1060 }
1061 
1062 
1063 
1064 DEAL_II_NAMESPACE_CLOSE
1065 
1066 
1067 #endif // DEAL_II_WITH_PETSC
1068 
1069 /*---------------------------- petsc_precondition.h ---------------------------*/
1070 
1071 #endif
1072 /*---------------------------- petsc_precondition.h ---------------------------*/
void vmult(VectorBase &dst, const VectorBase &src) const