Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
slepc_spectral_transformation.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2009 - 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 
17 #ifndef dealii_slepc_spectral_transformation_h
18 # define dealii_slepc_spectral_transformation_h
19 
20 
21 # include <deal.II/base/config.h>
22 
23 # ifdef DEAL_II_WITH_SLEPC
24 
25 # include <deal.II/lac/exceptions.h>
26 # include <deal.II/lac/petsc_solver.h>
27 
28 # include <petscksp.h>
29 
30 # include <slepceps.h>
31 
32 # include <memory>
33 
34 DEAL_II_NAMESPACE_OPEN
35 
36 namespace PETScWrappers
37 {
38  // forward declarations
39  class SolverBase;
40 } // namespace PETScWrappers
41 
42 namespace SLEPcWrappers
43 {
72  {
73  protected:
77  TransformationBase(const MPI_Comm &mpi_communicator);
78 
79  public:
83  virtual ~TransformationBase();
84 
93  void
94  set_matrix_mode(const STMatMode mode);
95 
100  void
101  set_solver(const PETScWrappers::SolverBase &solver);
102 
103  protected:
107  ST st;
108 
113  friend class SolverBase;
114  };
115 
123  {
124  public:
129  {
133  AdditionalData(const double shift_parameter = 0);
134 
138  const double shift_parameter;
139  };
140 
141 
145  TransformationShift(const MPI_Comm & mpi_communicator,
146  const AdditionalData &data = AdditionalData());
147 
148 
149  protected:
154  };
155 
164  {
165  public:
170  {
174  AdditionalData(const double shift_parameter = 0);
175 
179  const double shift_parameter;
180  };
181 
182 
186  TransformationShiftInvert(const MPI_Comm & mpi_communicator,
187  const AdditionalData &data = AdditionalData());
188 
189  protected:
194 
199  friend class SolverBase;
200  };
201 
211  {
212  public:
217  {
221  AdditionalData(const double shift_parameter = 0);
222 
226  const double shift_parameter;
227  };
228 
229 
234  const MPI_Comm & mpi_communicator,
235  const AdditionalData &data = AdditionalData());
236 
237  protected:
242  };
243 
251  {
252  public:
257  {
261  AdditionalData(const double shift_parameter = 0,
262  const double antishift_parameter = 0);
263 
267  const double shift_parameter;
268 
272  const double antishift_parameter;
273  };
274 
275 
279  TransformationCayley(const MPI_Comm & mpi_communicator,
280  const AdditionalData &data = AdditionalData());
281 
282  protected:
287  };
288 
289 } // namespace SLEPcWrappers
290 
291 DEAL_II_NAMESPACE_CLOSE
292 
293 # endif // DEAL_II_WITH_SLEPC
294 
295 /*-------------------- slepc_spectral_transformation.h ------------------*/
296 
297 #endif
298 
299 /*-------------------- slepc_spectral_transformation.h ------------------*/
void set_solver(const PETScWrappers::SolverBase &solver)
TransformationShiftInvert(const MPI_Comm &mpi_communicator, const AdditionalData &data=AdditionalData())
TransformationCayley(const MPI_Comm &mpi_communicator, const AdditionalData &data=AdditionalData())
TransformationBase(const MPI_Comm &mpi_communicator)
TransformationShift(const MPI_Comm &mpi_communicator, const AdditionalData &data=AdditionalData())
AdditionalData(const double shift_parameter=0, const double antishift_parameter=0)
TransformationSpectrumFolding(const MPI_Comm &mpi_communicator, const AdditionalData &data=AdditionalData())