forked from idaholab/moose
/
NonlinearEigenSystem.h
182 lines (147 loc) · 4.3 KB
/
NonlinearEigenSystem.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
//* This file is part of the MOOSE framework
//* https://www.mooseframework.org
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html
#pragma once
#include "libmesh/libmesh_config.h"
#include "NonlinearSystemBase.h"
#include "libmesh/transient_system.h"
#include "libmesh/eigen_system.h"
// forward declarations
class EigenProblem;
class KernelBase;
#if LIBMESH_HAVE_SLEPC
/**
* Nonlinear eigenvalue system to be solved
*/
class NonlinearEigenSystem : public NonlinearSystemBase
{
public:
NonlinearEigenSystem(EigenProblem & problem, const std::string & name);
virtual void solve() override;
/**
* Quit the current solve as soon as possible.
*/
virtual void stopSolve() override;
/**
* Returns the current nonlinear iteration number. In libmesh, this is
* updated during the nonlinear solve, so it should be up-to-date.
*/
virtual unsigned int getCurrentNonlinearIterationNumber() override;
virtual void setupFiniteDifferencedPreconditioner() override;
/**
* Returns the convergence state
*
* @return true if converged, otherwise false
*/
virtual bool converged() override;
virtual NumericVector<Number> & RHS() override;
/*
* A residual vector at MOOSE side for AX
*/
NumericVector<Number> & residualVectorAX();
/*
* A residual vector at MOOSE side for BX
*/
NumericVector<Number> & residualVectorBX();
/**
* Add the eigen tag to the right kernels
*/
template <typename T>
void addEigenTagToMooseObjects(MooseObjectTagWarehouse<T> & warehouse);
virtual void initialSetup() override;
/**
* Get the number of converged eigenvalues
*
* @return The number of converged eigenvalues
*/
virtual unsigned int getNumConvergedEigenvalues() const
{
return _transient_sys.get_n_converged();
};
virtual NonlinearSolver<Number> * nonlinearSolver() override;
virtual TransientEigenSystem & sys() { return _transient_sys; }
/**
* For eigenvalue problems (including standard and generalized), inhomogeneous (Dirichlet or
* Neumann)
* boundary conditions are not allowed.
*/
void checkIntegrity();
/**
* Return the Nth converged eigenvalue.
*
* @return The Nth converged eigenvalue as a complex number, i.e. the first and the second number
* is the real and the imaginary part of
* the eigenvalue, respectively.
*/
virtual const std::pair<Real, Real> getNthConvergedEigenvalue(dof_id_type n);
/**
* Get the number of converged eigenvalues
*
* @return all converged eigenvalues as complex numbers
*/
virtual const std::vector<std::pair<Real, Real>> & getAllConvergedEigenvalues()
{
return _eigen_values;
}
/**
* Vector tag ID of right hand side
*/
TagID eigenVectorTag() { return _Bx_tag; }
/**
* Vector tag ID of left hand side
*/
TagID nonEigenVectorTag() { return _Ax_tag; }
/**
* Matrix tag ID of right hand side
*/
TagID eigenMatrixTag() { return _B_tag; }
/**
* Matrix tag ID of left hand side
*/
TagID nonEigenMatrixTag() { return _A_tag; }
protected:
NumericVector<Number> & solutionOldInternal() override
{
return *_transient_sys.old_local_solution;
}
const NumericVector<Number> & solutionOldInternal() const override
{
return *_transient_sys.old_local_solution;
}
NumericVector<Number> & solutionOlderInternal() override
{
return *_transient_sys.older_local_solution;
}
const NumericVector<Number> & solutionOlderInternal() const override
{
return *_transient_sys.older_local_solution;
}
TransientEigenSystem & _transient_sys;
EigenProblem & _eigen_problem;
std::vector<std::pair<Real, Real>> _eigen_values;
unsigned int _n_eigen_pairs_required;
NumericVector<Number> & _work_rhs_vector_AX;
NumericVector<Number> & _work_rhs_vector_BX;
TagID _Ax_tag;
TagID _Bx_tag;
TagID _A_tag;
TagID _B_tag;
};
#else
class NonlinearEigenSystem : public libMesh::ParallelObject
{
public:
NonlinearEigenSystem(EigenProblem & problem, const std::string & name);
/**
* Returns the convergence state
* @return true if converged, otherwise false
*/
bool converged() { return false; }
void checkIntegrity() {}
};
#endif