Reference documentation for deal.II version Git 7026f387cc 2019-10-15 14:19:01 -0400
\(\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 // Forward declaration
37 # ifndef DOXYGEN
38 namespace PETScWrappers
39 {
40  // forward declarations
41  class SolverBase;
42 } // namespace PETScWrappers
43 # endif
44 
45 namespace SLEPcWrappers
46 {
75  {
76  protected:
80  TransformationBase(const MPI_Comm &mpi_communicator);
81 
82  public:
86  virtual ~TransformationBase();
87 
96  void
97  set_matrix_mode(const STMatMode mode);
98 
103  void
104  set_solver(const PETScWrappers::SolverBase &solver);
105 
106  protected:
110  ST st;
111 
112  // Make the solver class a friend, since it needs to set spectral
113  // transformation object.
114  friend class SolverBase;
115  };
116 
124  {
125  public:
130  {
134  AdditionalData(const double shift_parameter = 0);
135 
139  const double shift_parameter;
140  };
141 
142 
146  TransformationShift(const MPI_Comm & mpi_communicator,
147  const AdditionalData &data = AdditionalData());
148 
149 
150  protected:
155  };
156 
165  {
166  public:
171  {
175  AdditionalData(const double shift_parameter = 0);
176 
180  const double shift_parameter;
181  };
182 
183 
187  TransformationShiftInvert(const MPI_Comm & mpi_communicator,
188  const AdditionalData &data = AdditionalData());
189 
190  protected:
195 
196  // Make the solver class a friend, since it may need to set target
197  // equal the provided shift value.
198  friend class SolverBase;
199  };
200 
210  {
211  public:
216  {
220  AdditionalData(const double shift_parameter = 0);
221 
225  const double shift_parameter;
226  };
227 
228 
233  const MPI_Comm & mpi_communicator,
234  const AdditionalData &data = AdditionalData());
235 
236  protected:
241  };
242 
250  {
251  public:
256  {
260  AdditionalData(const double shift_parameter = 0,
261  const double antishift_parameter = 0);
262 
266  const double shift_parameter;
267 
271  const double antishift_parameter;
272  };
273 
274 
278  TransformationCayley(const MPI_Comm & mpi_communicator,
279  const AdditionalData &data = AdditionalData());
280 
281  protected:
286  };
287 
288 } // namespace SLEPcWrappers
289 
290 DEAL_II_NAMESPACE_CLOSE
291 
292 # endif // DEAL_II_WITH_SLEPC
293 
294 /*-------------------- slepc_spectral_transformation.h ------------------*/
295 
296 #endif
297 
298 /*-------------------- slepc_spectral_transformation.h ------------------*/