Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
petsc_solver.h
1 // ---------------------------------------------------------------------
2 // @f$Id: petsc_solver.h 30040 2013-07-18 17:06:48Z 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_solver_h
18 #define __deal2__petsc_solver_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 <deal.II/lac/solver_control.h>
27 # include <deal.II/base/std_cxx1x/shared_ptr.h>
28 
29 # include <petscksp.h>
30 
31 DEAL_II_NAMESPACE_OPEN
32 
33 
34 namespace PETScWrappers
35 {
36  // forward declarations
37  class MatrixBase;
38  class VectorBase;
39  class PreconditionerBase;
40 
41 
101  {
102  public:
120  const MPI_Comm &mpi_communicator);
121 
125  virtual ~SolverBase ();
126 
139  void
140  solve (const MatrixBase &A,
141  VectorBase &x,
142  const VectorBase &b,
143  const PreconditionerBase &preconditioner);
144 
145 
151  virtual void reset();
152 
153 
160  void set_prefix(const std::string &prefix);
161 
162 
167  SolverControl &control() const;
168 
172  DeclException1 (ExcPETScError,
173  int,
174  << "An error with error number " << arg1
175  << " occurred while calling a PETSc function");
176 
177  protected:
178 
190 
195  const MPI_Comm mpi_communicator;
196 
203  virtual void set_solver_type (KSP &ksp) const = 0;
204 
214  std::string prefix_name;
215 
216  private:
226  static
227  PetscErrorCode convergence_test (KSP ksp,
228  const PetscInt iteration,
229  const PetscReal residual_norm,
230  KSPConvergedReason *reason,
231  void *solver_control);
232 
262  struct SolverData
263  {
267  ~SolverData ();
268 
273  KSP ksp;
274  PC pc;
275  };
276 
283  std_cxx1x::shared_ptr<SolverData> solver_data;
284  };
285 
286 
287 
296  {
297  public:
304  {
310  AdditionalData (const double omega = 1);
311 
315  double omega;
316  };
317 
347  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
348  const AdditionalData &data = AdditionalData());
349 
350  protected:
356 
363  virtual void set_solver_type (KSP &ksp) const;
364  };
365 
366 
367 
376  {
377  public:
384  {};
385 
415  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
416  const AdditionalData &data = AdditionalData());
417 
418  protected:
424 
431  virtual void set_solver_type (KSP &ksp) const;
432  };
433 
434 
435 
443  class SolverCG : public SolverBase
444  {
445  public:
452  {};
453 
482  SolverCG (SolverControl &cn,
483  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
484  const AdditionalData &data = AdditionalData());
485 
486  protected:
492 
499  virtual void set_solver_type (KSP &ksp) const;
500  };
501 
502 
503 
511  class SolverBiCG : public SolverBase
512  {
513  public:
520  {};
521 
551  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
552  const AdditionalData &data = AdditionalData());
553 
554  protected:
560 
567  virtual void set_solver_type (KSP &ksp) const;
568  };
569 
570 
571 
578  class SolverGMRES : public SolverBase
579  {
580  public:
587  {
594  AdditionalData (const unsigned int restart_parameter = 30,
595  const bool right_preconditioning = false);
596 
601  unsigned int restart_parameter;
602 
608  };
609 
639  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
640  const AdditionalData &data = AdditionalData());
641 
642  protected:
648 
655  virtual void set_solver_type (KSP &ksp) const;
656  };
657 
658 
659 
667  class SolverBicgstab : public SolverBase
668  {
669  public:
676  {};
677 
707  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
708  const AdditionalData &data = AdditionalData());
709 
710  protected:
716 
723  virtual void set_solver_type (KSP &ksp) const;
724  };
725 
733  class SolverCGS : public SolverBase
734  {
735  public:
742  {};
743 
773  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
774  const AdditionalData &data = AdditionalData());
775 
776  protected:
782 
789  virtual void set_solver_type (KSP &ksp) const;
790  };
791 
792 
793 
801  class SolverTFQMR : public SolverBase
802  {
803  public:
810  {};
811 
841  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
842  const AdditionalData &data = AdditionalData());
843 
844  protected:
850 
857  virtual void set_solver_type (KSP &ksp) const;
858  };
859 
860 
861 
862 
874  class SolverTCQMR : public SolverBase
875  {
876  public:
883  {};
884 
914  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
915  const AdditionalData &data = AdditionalData());
916 
917  protected:
923 
930  virtual void set_solver_type (KSP &ksp) const;
931  };
932 
933 
934 
942  class SolverCR : public SolverBase
943  {
944  public:
951  {};
952 
981  SolverCR (SolverControl &cn,
982  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
983  const AdditionalData &data = AdditionalData());
984 
985  protected:
991 
998  virtual void set_solver_type (KSP &ksp) const;
999  };
1000 
1001 
1002 
1010  class SolverLSQR : public SolverBase
1011  {
1012  public:
1019  {};
1020 
1050  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
1051  const AdditionalData &data = AdditionalData());
1052 
1053  protected:
1059 
1066  virtual void set_solver_type (KSP &ksp) const;
1067  };
1068 
1069 
1082  class SolverPreOnly : public SolverBase
1083  {
1084  public:
1091  {};
1092 
1122  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
1123  const AdditionalData &data = AdditionalData());
1124 
1125  protected:
1131 
1138  virtual void set_solver_type (KSP &ksp) const;
1139  };
1140 
1166  {
1167  public:
1174  {};
1179  const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
1180  const AdditionalData &data = AdditionalData());
1181 
1186  void solve (const MatrixBase &A,
1187  VectorBase &x,
1188  const VectorBase &b);
1189 
1198  void set_symmetric_mode (const bool flag);
1199 
1200  protected:
1206 
1207  virtual void set_solver_type (KSP &ksp) const;
1208 
1209  private:
1218  static
1219  PetscErrorCode convergence_test (KSP ksp,
1220  const PetscInt iteration,
1221  const PetscReal residual_norm,
1222  KSPConvergedReason *reason,
1223  void *solver_control);
1224 
1235  {
1239  ~SolverDataMUMPS ();
1240 
1241  KSP ksp;
1242  PC pc;
1243  };
1244 
1245  std_cxx1x::shared_ptr<SolverDataMUMPS> solver_data;
1246 
1255  };
1256 }
1257 
1258 DEAL_II_NAMESPACE_CLOSE
1259 
1260 #endif // DEAL_II_WITH_PETSC
1261 
1262 /*---------------------------- petsc_solver.h ---------------------------*/
1263 
1264 #endif
1265 /*---------------------------- petsc_solver.h ---------------------------*/
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b, const PreconditionerBase &preconditioner)
virtual void set_solver_type(KSP &ksp) const =0
const AdditionalData additional_data
Definition: petsc_solver.h:491
DeclException1(ExcPETScError, int,<< "An error with error number "<< arg1<< " occurred while calling a PETSc function")
const AdditionalData additional_data
const AdditionalData additional_data
Definition: petsc_solver.h:355
SolverCG(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:423
const AdditionalData additional_data
Definition: petsc_solver.h:922
const AdditionalData additional_data
Definition: petsc_solver.h:647
virtual void set_solver_type(KSP &ksp) const
virtual void set_solver_type(KSP &ksp) const
virtual void set_solver_type(KSP &ksp) const
virtual void set_solver_type(KSP &ksp) const
SolverBicgstab(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const
SparseDirectMUMPS(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const
void set_prefix(const std::string &prefix)
const MPI_Comm mpi_communicator
Definition: petsc_solver.h:195
static PetscErrorCode convergence_test(KSP ksp, const PetscInt iteration, const PetscReal residual_norm, KSPConvergedReason *reason, void *solver_control)
SolverBase(SolverControl &cn, const MPI_Comm &mpi_communicator)
const AdditionalData additional_data
virtual void set_solver_type(KSP &ksp) const
virtual void set_solver_type(KSP &ksp) const
SolverCR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const
SolverControl & solver_control
Definition: petsc_solver.h:189
const AdditionalData additional_data
Definition: petsc_solver.h:559
SolverBiCG(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverLSQR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverChebychev(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverCGS(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverRichardson(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const
const AdditionalData additional_data
Definition: petsc_solver.h:715
virtual void set_solver_type(KSP &ksp) const
virtual void set_solver_type(KSP &ksp) const
const AdditionalData additional_data
AdditionalData(const unsigned int restart_parameter=30, const bool right_preconditioning=false)
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b)
SolverTCQMR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const
const AdditionalData additional_data
Definition: petsc_solver.h:781
SolverTFQMR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverGMRES(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverControl & control() const
void set_symmetric_mode(const bool flag)
std_cxx1x::shared_ptr< SolverData > solver_data
Definition: petsc_solver.h:283
SolverPreOnly(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
static PetscErrorCode convergence_test(KSP ksp, const PetscInt iteration, const PetscReal residual_norm, KSPConvergedReason *reason, void *solver_control)
const AdditionalData additional_data
Definition: petsc_solver.h:990
const AdditionalData additional_data
Definition: petsc_solver.h:849