Reference documentation for deal.II version GIT c14369f203 2023-10-01 07:40:02+00:00
\(\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_solver.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2023 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 
25 
26 # include <deal.II/lac/exceptions.h>
28 
29 # include <petscksp.h>
30 
31 # ifdef DEAL_II_WITH_SLEPC
33 # endif
34 
36 
37 // Forward declarations
38 # ifndef DOXYGEN
39 # ifdef DEAL_II_WITH_SLEPC
40 namespace SLEPcWrappers
41 {
42  // forward declarations
43  class TransformationBase;
44 } // namespace SLEPcWrappers
45 # endif
46 # endif
47 
48 namespace PETScWrappers
49 {
50  // forward declarations
51 # ifndef DOXYGEN
52  class MatrixBase;
53  class VectorBase;
54  class PreconditionBase;
55 # endif
56 
57 
92  class SolverBase
93  {
94  public:
99  SolverBase();
100 
105 
109  virtual ~SolverBase();
110 
118  void
119  solve(const MatrixBase &A,
120  VectorBase &x,
121  const VectorBase &b,
122  const PreconditionBase &preconditioner);
123 
128  virtual void
129  reset();
130 
135  void
136  set_prefix(const std::string &prefix);
137 
141  SolverControl &
142  control() const;
143 
148  void
149  initialize(const PreconditionBase &preconditioner);
150 
154  KSP
155  petsc_ksp();
156 
162  operator KSP() const;
163 
164  protected:
168  KSP ksp;
169 
177 
181  void
182  initialize_ksp_with_comm(const MPI_Comm comm);
183 
188  virtual void
189  set_solver_type(KSP &ksp) const;
190 
197  void
199 
206  std::string prefix_name;
207 
208  private:
215  static PetscErrorCode
216  convergence_test(KSP ksp,
217  const PetscInt iteration,
218  const PetscReal residual_norm,
219  KSPConvergedReason *reason,
220  void *solver_control);
221 
222 
223 # ifdef DEAL_II_WITH_SLEPC
224  // Make the transformation class a friend, since it needs to set the KSP
225  // solver.
227 # endif
228  };
229 
230 
231 
239  {
240  public:
245  {
249  explicit AdditionalData(const double omega = 1);
250 
254  double omega;
255  };
256 
265  const AdditionalData &data = AdditionalData());
266 
275  const MPI_Comm mpi_communicator,
276  const AdditionalData &data = AdditionalData());
277 
278  protected:
283 
288  virtual void
289  set_solver_type(KSP &ksp) const override;
290  };
291 
292 
293 
301  {
302  public:
307  {};
308 
317  const AdditionalData &data = AdditionalData());
318 
327  const MPI_Comm mpi_communicator,
328  const AdditionalData &data = AdditionalData());
329 
330  protected:
335 
340  virtual void
341  set_solver_type(KSP &ksp) const override;
342  };
343 
344 
345 
351  class SolverCG : public SolverBase
352  {
353  public:
358  {};
359 
367  SolverCG(SolverControl &cn, const AdditionalData &data = AdditionalData());
368 
377  const MPI_Comm mpi_communicator,
378  const AdditionalData &data = AdditionalData());
379 
380  protected:
385 
390  virtual void
391  set_solver_type(KSP &ksp) const override;
392  };
393 
394 
395 
401  class SolverBiCG : public SolverBase
402  {
403  public:
408  {};
409 
418  const AdditionalData &data = AdditionalData());
419 
428  const MPI_Comm mpi_communicator,
429  const AdditionalData &data = AdditionalData());
430 
431  protected:
436 
441  virtual void
442  set_solver_type(KSP &ksp) const override;
443  };
444 
445 
446 
452  class SolverGMRES : public SolverBase
453  {
454  public:
459  {
464  AdditionalData(const unsigned int restart_parameter = 30,
465  const bool right_preconditioning = false);
466 
470  unsigned int restart_parameter;
471 
476  };
477 
486  const AdditionalData &data = AdditionalData());
487 
496  const MPI_Comm mpi_communicator,
497  const AdditionalData &data = AdditionalData());
498 
499  protected:
504 
509  virtual void
510  set_solver_type(KSP &ksp) const override;
511  };
512 
513 
514 
521  class SolverBicgstab : public SolverBase
522  {
523  public:
528  {};
529 
538  const AdditionalData &data = AdditionalData());
539 
548  const MPI_Comm mpi_communicator,
549  const AdditionalData &data = AdditionalData());
550 
551  protected:
556 
561  virtual void
562  set_solver_type(KSP &ksp) const override;
563  };
564 
565 
566 
573  class SolverCGS : public SolverBase
574  {
575  public:
580  {};
581 
590 
599  const MPI_Comm mpi_communicator,
600  const AdditionalData &data = AdditionalData());
601 
602  protected:
607 
612  virtual void
613  set_solver_type(KSP &ksp) const override;
614  };
615 
616 
617 
623  class SolverTFQMR : public SolverBase
624  {
625  public:
630  {};
631 
640  const AdditionalData &data = AdditionalData());
641 
650  const MPI_Comm mpi_communicator,
651  const AdditionalData &data = AdditionalData());
652 
653  protected:
658 
663  virtual void
664  set_solver_type(KSP &ksp) const override;
665  };
666 
667 
668 
679  class SolverTCQMR : public SolverBase
680  {
681  public:
686  {};
687 
696  const AdditionalData &data = AdditionalData());
697 
706  const MPI_Comm mpi_communicator,
707  const AdditionalData &data = AdditionalData());
708 
709  protected:
714 
719  virtual void
720  set_solver_type(KSP &ksp) const override;
721  };
722 
723 
724 
730  class SolverCR : public SolverBase
731  {
732  public:
737  {};
738 
746  SolverCR(SolverControl &cn, const AdditionalData &data = AdditionalData());
747 
756  const MPI_Comm mpi_communicator,
757  const AdditionalData &data = AdditionalData());
758 
759  protected:
764 
769  virtual void
770  set_solver_type(KSP &ksp) const override;
771  };
772 
773 
774 
781  class SolverLSQR : public SolverBase
782  {
783  public:
788  {};
789 
798  const AdditionalData &data = AdditionalData());
799 
808  const MPI_Comm mpi_communicator,
809  const AdditionalData &data = AdditionalData());
810 
811  protected:
816 
821  virtual void
822  set_solver_type(KSP &ksp) const override;
823  };
824 
825 
837  class SolverPreOnly : public SolverBase
838  {
839  public:
844  {};
845 
854  const AdditionalData &data = AdditionalData());
855 
864  const MPI_Comm mpi_communicator,
865  const AdditionalData &data = AdditionalData());
866 
867  protected:
872 
877  virtual void
878  set_solver_type(KSP &ksp) const override;
879  };
880 
905  {
906  public:
911  {};
912 
917  const AdditionalData &data = AdditionalData());
918 
927  const MPI_Comm mpi_communicator,
928  const AdditionalData &data = AdditionalData());
929 
933  void
934  solve(const MatrixBase &A, VectorBase &x, const VectorBase &b);
935 
941  void
942  set_symmetric_mode(const bool flag);
943 
944  protected:
949 
950  virtual void
951  set_solver_type(KSP &ksp) const override;
952 
958  };
959 } // namespace PETScWrappers
960 
962 
963 #endif // DEAL_II_WITH_PETSC
964 
965 #endif
SolverControl & control() const
void initialize(const PreconditionBase &preconditioner)
void perhaps_set_convergence_test() const
static PetscErrorCode convergence_test(KSP ksp, const PetscInt iteration, const PetscReal residual_norm, KSPConvergedReason *reason, void *solver_control)
SmartPointer< SolverControl, SolverBase > solver_control
Definition: petsc_solver.h:176
void set_prefix(const std::string &prefix)
void initialize_ksp_with_comm(const MPI_Comm comm)
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b, const PreconditionBase &preconditioner)
Definition: petsc_solver.cc:84
virtual void set_solver_type(KSP &ksp) const
Definition: petsc_solver.cc:56
SolverBiCG(SolverControl &cn, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:435
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:555
SolverBicgstab(SolverControl &cn, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:606
virtual void set_solver_type(KSP &ksp) const override
SolverCGS(SolverControl &cn, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:384
SolverCG(SolverControl &cn, const AdditionalData &data=AdditionalData())
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:763
SolverCR(SolverControl &cn, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:334
virtual void set_solver_type(KSP &ksp) const override
SolverChebychev(SolverControl &cn, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:503
virtual void set_solver_type(KSP &ksp) const override
SolverGMRES(SolverControl &cn, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverLSQR(SolverControl &cn, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:815
SolverPreOnly(SolverControl &cn, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:871
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:282
SolverRichardson(SolverControl &cn, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:713
virtual void set_solver_type(KSP &ksp) const override
SolverTCQMR(SolverControl &cn, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:657
virtual void set_solver_type(KSP &ksp) const override
SolverTFQMR(SolverControl &cn, const AdditionalData &data=AdditionalData())
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b)
virtual void set_solver_type(KSP &ksp) const override
const AdditionalData additional_data
Definition: petsc_solver.h:948
void set_symmetric_mode(const bool flag)
SparseDirectMUMPS(SolverControl &cn, const AdditionalData &data=AdditionalData())
#define DEAL_II_DEPRECATED
Definition: config.h:178
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static const char A
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
AdditionalData(const unsigned int restart_parameter=30, const bool right_preconditioning=false)
const MPI_Comm comm