Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
petsc_solver.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_solver_h
17 # define dealii_petsc_solver_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 # include <deal.II/lac/solver_control.h>
26 
27 # include <petscksp.h>
28 
29 # include <memory>
30 
31 # ifdef DEAL_II_WITH_SLEPC
32 # include <deal.II/lac/slepc_spectral_transformation.h>
33 # endif
34 
35 DEAL_II_NAMESPACE_OPEN
36 
37 # ifdef DEAL_II_WITH_SLEPC
38 namespace SLEPcWrappers
39 {
40  // forward declarations
41  class TransformationBase;
42 } // namespace SLEPcWrappers
43 # endif
44 
45 namespace PETScWrappers
46 {
47  // forward declarations
48  class MatrixBase;
49  class VectorBase;
50  class PreconditionerBase;
51 
52 
109  {
110  public:
120  SolverBase(SolverControl &cn, const MPI_Comm &mpi_communicator);
121 
125  virtual ~SolverBase() = default;
126 
134  void
135  solve(const MatrixBase & A,
136  VectorBase & x,
137  const VectorBase & b,
138  const PreconditionerBase &preconditioner);
139 
140 
145  virtual void
146  reset();
147 
148 
153  void
154  set_prefix(const std::string &prefix);
155 
156 
160  SolverControl &
161  control() const;
162 
167  void
168  initialize(const PreconditionerBase &preconditioner);
169 
170  protected:
178 
182  const MPI_Comm mpi_communicator;
183 
188  virtual void
189  set_solver_type(KSP &ksp) const = 0;
190 
197  std::string prefix_name;
198 
199  private:
206  static PetscErrorCode
207  convergence_test(KSP ksp,
208  const PetscInt iteration,
209  const PetscReal residual_norm,
210  KSPConvergedReason *reason,
211  void * solver_control);
212 
230  struct SolverData
231  {
235  ~SolverData();
236 
240  KSP ksp;
241  };
242 
247  std::unique_ptr<SolverData> solver_data;
248 
249 # ifdef DEAL_II_WITH_SLEPC
250 
255 # endif
256  };
257 
258 
259 
268  {
269  public:
274  {
278  explicit AdditionalData(const double omega = 1);
279 
283  double omega;
284  };
285 
303  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
304  const AdditionalData &data = AdditionalData());
305 
306  protected:
311 
316  virtual void
317  set_solver_type(KSP &ksp) const override;
318  };
319 
320 
321 
330  {
331  public:
336  {};
337 
355  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
356  const AdditionalData &data = AdditionalData());
357 
358  protected:
363 
368  virtual void
369  set_solver_type(KSP &ksp) const override;
370  };
371 
372 
373 
380  class SolverCG : public SolverBase
381  {
382  public:
387  {};
388 
405  SolverCG(SolverControl & cn,
406  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
407  const AdditionalData &data = AdditionalData());
408 
409  protected:
414 
419  virtual void
420  set_solver_type(KSP &ksp) const override;
421  };
422 
423 
424 
431  class SolverBiCG : public SolverBase
432  {
433  public:
438  {};
439 
457  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
458  const AdditionalData &data = AdditionalData());
459 
460  protected:
465 
470  virtual void
471  set_solver_type(KSP &ksp) const override;
472  };
473 
474 
475 
482  class SolverGMRES : public SolverBase
483  {
484  public:
489  {
494  AdditionalData(const unsigned int restart_parameter = 30,
495  const bool right_preconditioning = false);
496 
500  unsigned int restart_parameter;
501 
506  };
507 
525  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
526  const AdditionalData &data = AdditionalData());
527 
528  protected:
533 
538  virtual void
539  set_solver_type(KSP &ksp) const override;
540  };
541 
542 
543 
551  class SolverBicgstab : public SolverBase
552  {
553  public:
558  {};
559 
577  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
578  const AdditionalData &data = AdditionalData());
579 
580  protected:
585 
590  virtual void
591  set_solver_type(KSP &ksp) const override;
592  };
593 
601  class SolverCGS : public SolverBase
602  {
603  public:
608  {};
609 
627  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
628  const AdditionalData &data = AdditionalData());
629 
630  protected:
635 
640  virtual void
641  set_solver_type(KSP &ksp) const override;
642  };
643 
644 
645 
652  class SolverTFQMR : public SolverBase
653  {
654  public:
659  {};
660 
678  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
679  const AdditionalData &data = AdditionalData());
680 
681  protected:
686 
691  virtual void
692  set_solver_type(KSP &ksp) const override;
693  };
694 
695 
696 
708  class SolverTCQMR : public SolverBase
709  {
710  public:
715  {};
716 
734  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
735  const AdditionalData &data = AdditionalData());
736 
737  protected:
742 
747  virtual void
748  set_solver_type(KSP &ksp) const override;
749  };
750 
751 
752 
759  class SolverCR : public SolverBase
760  {
761  public:
766  {};
767 
784  SolverCR(SolverControl & cn,
785  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
786  const AdditionalData &data = AdditionalData());
787 
788  protected:
793 
798  virtual void
799  set_solver_type(KSP &ksp) const override;
800  };
801 
802 
803 
811  class SolverLSQR : public SolverBase
812  {
813  public:
818  {};
819 
837  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
838  const AdditionalData &data = AdditionalData());
839 
840  protected:
845 
850  virtual void
851  set_solver_type(KSP &ksp) const override;
852  };
853 
854 
867  class SolverPreOnly : public SolverBase
868  {
869  public:
874  {};
875 
893  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
894  const AdditionalData &data = AdditionalData());
895 
896  protected:
901 
906  virtual void
907  set_solver_type(KSP &ksp) const override;
908  };
909 
935  {
936  public:
941  {};
946  const MPI_Comm & mpi_communicator = PETSC_COMM_SELF,
947  const AdditionalData &data = AdditionalData());
948 
952  void
953  solve(const MatrixBase &A, VectorBase &x, const VectorBase &b);
954 
960  void
961  set_symmetric_mode(const bool flag);
962 
963  protected:
968 
969  virtual void
970  set_solver_type(KSP &ksp) const override;
971 
972  private:
979  static PetscErrorCode
980  convergence_test(KSP ksp,
981  const PetscInt iteration,
982  const PetscReal residual_norm,
983  KSPConvergedReason *reason,
984  void * solver_control);
985 
992  {
997 
998  KSP ksp;
999  PC pc;
1000  };
1001 
1002  std::unique_ptr<SolverDataMUMPS> solver_data;
1003 
1009  };
1010 } // namespace PETScWrappers
1011 
1012 DEAL_II_NAMESPACE_CLOSE
1013 
1014 # endif // DEAL_II_WITH_PETSC
1015 
1016 /*---------------------------- petsc_solver.h ---------------------------*/
1017 
1018 #endif
1019 /*---------------------------- petsc_solver.h ---------------------------*/
virtual void set_solver_type(KSP &ksp) const override
virtual void set_solver_type(KSP &ksp) const =0
const AdditionalData additional_data
Definition: petsc_solver.h:413
virtual void set_solver_type(KSP &ksp) const override
virtual void set_solver_type(KSP &ksp) const override
const AdditionalData additional_data
Definition: petsc_solver.h:844
const AdditionalData additional_data
Definition: petsc_solver.h:310
void set_prefix(const std::string &prefix)
const AdditionalData additional_data
Definition: petsc_solver.h:362
AdditionalData(const unsigned int restart_parameter=30, const bool right_preconditioning=false)
SolverRichardson(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:741
const AdditionalData additional_data
Definition: petsc_solver.h:532
void initialize(const PreconditionerBase &preconditioner)
virtual void set_solver_type(KSP &ksp) const override
SparseDirectMUMPS(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
const MPI_Comm mpi_communicator
Definition: petsc_solver.h:182
std::unique_ptr< SolverData > solver_data
Definition: petsc_solver.h:247
virtual void set_solver_type(KSP &ksp) const override
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b, const PreconditionerBase &preconditioner)
Definition: petsc_solver.cc:53
virtual void set_solver_type(KSP &ksp) const override
virtual void set_solver_type(KSP &ksp) const override
const AdditionalData additional_data
Definition: petsc_solver.h:900
SolverBase(SolverControl &cn, const MPI_Comm &mpi_communicator)
Definition: petsc_solver.cc:45
virtual void set_solver_type(KSP &ksp) const override
virtual void set_solver_type(KSP &ksp) const override
virtual void set_solver_type(KSP &ksp) const override
SolverCR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverControl & solver_control
Definition: petsc_solver.h:177
const AdditionalData additional_data
Definition: petsc_solver.h:464
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())
const AdditionalData additional_data
Definition: petsc_solver.h:584
virtual ~SolverBase()=default
SolverBicgstab(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:967
SolverControl & control() const
virtual void set_solver_type(KSP &ksp) const override
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())
const AdditionalData additional_data
Definition: petsc_solver.h:634
SolverTFQMR(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)
void set_symmetric_mode(const bool flag)
static PetscErrorCode convergence_test(KSP ksp, const PetscInt iteration, const PetscReal residual_norm, KSPConvergedReason *reason, void *solver_control)
SolverCG(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverPreOnly(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:792
SolverGMRES(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
const AdditionalData additional_data
Definition: petsc_solver.h:685
virtual void set_solver_type(KSP &ksp) const override