DOLFIN
DOLFIN C++ interface
Loading...
Searching...
No Matches
SLEPcEigenSolver.h
1// Copyright (C) 2005-2017 Garth N. Wells
2//
3// This file is part of DOLFIN.
4//
5// DOLFIN is free software: you can redistribute it and/or modify
6// it under the terms of the GNU Lesser General Public License as published by
7// the Free Software Foundation, either version 3 of the License, or
8// (at your option) any later version.
9//
10// DOLFIN is distributed in the hope that it will be useful,
11// but WITHOUT ANY WARRANTY; without even the implied warranty of
12// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13// GNU Lesser General Public License for more details.
14//
15// You should have received a copy of the GNU Lesser General Public License
16// along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
17
18#ifndef __SLEPC_EIGEN_SOLVER_H
19#define __SLEPC_EIGEN_SOLVER_H
20
21#ifdef HAS_SLEPC
22
23#include <memory>
24#include <string>
25#include <slepceps.h>
26#include "dolfin/common/types.h"
27#include "dolfin/common/MPI.h"
28#include "PETScObject.h"
29
30namespace dolfin
31{
32
34 class GenericVector;
35 class PETScMatrix;
36 class PETScVector;
37 class VectorSpaceBasis;
38
117
118 class SLEPcEigenSolver : public Variable, public PETScObject
119 {
120 public:
121
123 explicit SLEPcEigenSolver(MPI_Comm comm);
124
126 explicit SLEPcEigenSolver(EPS eps);
127
129 explicit SLEPcEigenSolver(std::shared_ptr<const PETScMatrix> A);
130
132 SLEPcEigenSolver(MPI_Comm comm, std::shared_ptr<const PETScMatrix> A);
133
135 SLEPcEigenSolver(std::shared_ptr<const PETScMatrix> A,
136 std::shared_ptr<const PETScMatrix> B);
137
139 SLEPcEigenSolver(MPI_Comm comm, std::shared_ptr<const PETScMatrix> A,
140 std::shared_ptr<const PETScMatrix> B);
141
144
147 void set_operators(std::shared_ptr<const PETScMatrix> A,
148 std::shared_ptr<const PETScMatrix> B);
149
151 void solve();
152
154 void solve(std::size_t n);
155
157 void get_eigenvalue(double& lr, double& lc, std::size_t i) const;
158
160 void get_eigenpair(double& lr, double& lc,
161 GenericVector& r, GenericVector& c, std::size_t i) const;
162
164 void get_eigenpair(double& lr, double& lc,
165 PETScVector& r, PETScVector& c, std::size_t i) const;
166
168 std::size_t get_iteration_number() const;
169
171 std::size_t get_number_converged() const;
172
175 void set_deflation_space(const VectorSpaceBasis& deflation_space);
176
179 void set_initial_space(const VectorSpaceBasis& initial_space);
180
183 void set_options_prefix(std::string options_prefix);
184
187 std::string get_options_prefix() const;
188
190 void set_from_options() const;
191
193 EPS eps() const;
194
197 {
198 Parameters p("slepc_eigenvalue_solver");
199 p.add<std::string>("problem_type");
200 p.add<std::string>("spectrum");
201 p.add<std::string>("solver");
202 p.add<double>("tolerance");
203 p.add<int>("maximum_iterations");
204 p.add<std::string>("spectral_transform");
205 p.add<double>("spectral_shift");
206 p.add<bool>("verbose");
207
208 return p;
209 }
210
211 private:
212
214 void read_parameters();
215
216 // Set problem type (used for SLEPc internals)
217 void set_problem_type(std::string type);
218
219 // Set spectral transform
220 void set_spectral_transform(std::string transform, double shift);
221
222 // Set spectrum
223 void set_spectrum(std::string solver);
224
225 // Set solver
226 void set_solver(std::string spectrum);
227
228 // Set tolerance
229 void set_tolerance(double tolerance, int maxiter);
230
231 // SLEPc solver pointer
232 EPS _eps;
233
234 };
235
236}
237
238#endif
239
240#endif
This class defines a common interface for vectors.
Definition GenericVector.h:48
Definition PETScObject.h:34
Definition PETScVector.h:61
Definition Parameters.h:95
void add(std::string key)
Definition Parameters.h:128
Definition SLEPcEigenSolver.h:119
std::string get_options_prefix() const
Definition SLEPcEigenSolver.cpp:343
std::size_t get_number_converged() const
Get the number of converged eigenvalues.
Definition SLEPcEigenSolver.cpp:289
void solve()
Compute all eigenpairs of the matrix A (solve Ax = \lambda x)
Definition SLEPcEigenSolver.cpp:153
SLEPcEigenSolver(MPI_Comm comm)
Create eigenvalue solver.
Definition SLEPcEigenSolver.cpp:32
void set_from_options() const
Set options from PETSc options database.
Definition SLEPcEigenSolver.cpp:352
void set_deflation_space(const VectorSpaceBasis &deflation_space)
Definition SLEPcEigenSolver.cpp:297
EPS eps() const
Return SLEPc EPS pointer.
Definition SLEPcEigenSolver.cpp:540
void set_initial_space(const VectorSpaceBasis &initial_space)
Definition SLEPcEigenSolver.cpp:316
void set_options_prefix(std::string options_prefix)
Definition SLEPcEigenSolver.cpp:335
void set_operators(std::shared_ptr< const PETScMatrix > A, std::shared_ptr< const PETScMatrix > B)
Definition SLEPcEigenSolver.cpp:142
~SLEPcEigenSolver()
Destructor.
Definition SLEPcEigenSolver.cpp:135
static Parameters default_parameters()
Default parameter values.
Definition SLEPcEigenSolver.h:196
void get_eigenpair(double &lr, double &lc, GenericVector &r, GenericVector &c, std::size_t i) const
Get ith eigenpair.
Definition SLEPcEigenSolver.cpp:240
std::size_t get_iteration_number() const
Get the number of iterations used by the solver.
Definition SLEPcEigenSolver.cpp:532
void get_eigenvalue(double &lr, double &lc, std::size_t i) const
Get ith eigenvalue.
Definition SLEPcEigenSolver.cpp:220
Common base class for DOLFIN variables.
Definition Variable.h:36
Definition VectorSpaceBasis.h:34
Definition adapt.h:30