www.mooseframework.org
Eigenvalue.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 // MOOSE includes
11 #include "Eigenvalue.h"
12 #include "EigenProblem.h"
13 #include "Factory.h"
14 #include "MooseApp.h"
15 #include "NonlinearEigenSystem.h"
16 #include "SlepcSupport.h"
17 #include "UserObject.h"
18 
19 registerMooseObject("MooseApp", Eigenvalue);
20 
23 {
25 
26  params.addClassDescription(
27  "Eigenvalue solves a standard/generalized linear or nonlinear eigenvalue problem");
28 
29  params += FEProblemSolve::validParams();
30  params.addParam<Real>("time", 0.0, "System time");
31 
32  // matrix_free will be an invalid for griffin once the integration is done.
33  // In this PR, we can not change it. It will still be a valid option when users
34  // use non-Newton algorithms
35  params.addParam<bool>(
36  "matrix_free",
37  false,
38  "Whether or not to use a matrix free fashion to form operators. "
39  "If true, shell matrices will be used and meanwhile a preconditioning matrix"
40  "may be formed as well.");
41 
42  params.addParam<bool>(
43  "precond_matrix_free",
44  false,
45  "Whether or not to use a matrix free fashion for forming the preconditioning matrix. "
46  "If true, a shell matrix will be used for preconditioner.");
47 
48  params.addParam<bool>("constant_matrices",
49  false,
50  "Whether or not to use constant matrices so that we can use them to form "
51  "residuals on both linear and "
52  "nonlinear iterations");
53 
54  params.addParam<bool>("precond_matrix_includes_eigen",
55  false,
56  "Whether or not to include eigen kernels in the preconditioning matrix. "
57  "If true, the preconditioning matrix will include eigen kernels.");
58 
59  params.addPrivateParam<bool>("_use_eigen_value", true);
60 
61  params.addParam<Real>("initial_eigenvalue", 1, "Initial eigenvalue");
62  params.addParam<PostprocessorName>(
63  "normalization", "Postprocessor evaluating norm of eigenvector for normalization");
64  params.addParam<Real>("normal_factor",
65  "Normalize eigenvector to make a defined norm equal to this factor");
66 
67  params.addParam<bool>("auto_initialization",
68  true,
69  "If true, we will set an initial eigen vector in moose, otherwise EPS "
70  "solver will initial eigen vector");
71 
72  params.addParamNamesToGroup("matrix_free precond_matrix_free constant_matrices "
73  "precond_matrix_includes_eigen",
74  "Matrix and Matrix-Free");
75  params.addParamNamesToGroup("initial_eigenvalue auto_initialization",
76  "Eigenvector and eigenvalue initialization");
77  params.addParamNamesToGroup("normalization normal_factor", "Solution normalization");
78 
79  // If Newton and Inverse Power is combined in SLEPc side
80  params.addPrivateParam<bool>("_newton_inverse_power", false);
81 
82 // Add slepc options and eigen problems
83 #ifdef LIBMESH_HAVE_SLEPC
85 
87 #endif
88  return params;
89 }
90 
92  : Executioner(parameters),
93  _eigen_problem(*getCheckedPointerParam<EigenProblem *>(
94  "_eigen_problem", "This might happen if you don't have a mesh")),
95  _feproblem_solve(*this),
96  _normalization(isParamValid("normalization") ? &getPostprocessorValue("normalization")
97  : nullptr),
98  _system_time(getParam<Real>("time")),
99  _time_step(_eigen_problem.timeStep()),
100  _time(_eigen_problem.time()),
101  _final_timer(registerTimedSection("final", 1))
102 {
103 // Extract and store SLEPc options
104 #ifdef LIBMESH_HAVE_SLEPC
106 
109 
110  // pass two control parameters to eigen problem
112  getParam<unsigned int>("free_power_iterations");
114  getParam<unsigned int>("extra_power_iterations");
115 
116  if (!isParamValid("normalization") && isParamValid("normal_factor"))
117  paramError("normal_factor",
118  "Cannot set scaling factor without defining normalization postprocessor.");
119 
120  _fixed_point_solve->setInnerSolve(_feproblem_solve);
122 
123  if (isParamValid("normalization"))
124  {
125  const auto & normpp = getParam<PostprocessorName>("normalization");
126  if (isParamValid("normal_factor"))
127  _eigen_problem.setNormalization(normpp, getParam<Real>("normal_factor"));
128  else
130  }
131 
132  _eigen_problem.setInitialEigenvalue(getParam<Real>("initial_eigenvalue"));
133 
134  // Set a flag to nonlinear eigen system
135  _eigen_problem.getNonlinearEigenSystem(/*nl_sys_num=*/0)
136  .precondMatrixIncludesEigenKernels(getParam<bool>("precond_matrix_includes_eigen"));
137 #else
138  mooseError("SLEPc is required to use Eigenvalue executioner, please use '--download-slepc in "
139  "PETSc configuration'");
140 #endif
141  // SLEPc older than 3.13.0 can not take initial guess from moose
142  // It may generate converge issues
143 #if PETSC_RELEASE_LESS_THAN(3, 13, 0)
145  "Please use SLEPc-3.13.0 or higher. Old versions of SLEPc likely produce bad convergence");
146 #endif
147 }
148 
149 #ifdef LIBMESH_HAVE_SLEPC
150 void
152 {
153  if (isParamValid("normalization"))
154  {
155  const auto & normpp = getParam<PostprocessorName>("normalization");
156  const auto & exec = _eigen_problem.getUserObject<UserObject>(normpp).getExecuteOnEnum();
157  if (!exec.contains(EXEC_LINEAR))
158  mooseError("Normalization postprocessor ", normpp, " requires execute_on = 'linear'");
159  }
160 
161  // Does not allow time kernels
162  checkIntegrity();
163 
164  // Provide vector of ones to solver
165  // "auto_initialization" is on by default and we init the vector values associated
166  // with eigen-variables as ones. If "auto_initialization" is turned off by users,
167  // it is up to users to provide an initial guess. If "auto_initialization" is off
168  // and users does not provide an initial guess, slepc will automatically generate
169  // a random vector as the initial guess. The motivation to offer this option is
170  // that we have to initialize ONLY eigen variables in multiphysics simulation.
171  // auto_initialization can be overriden by initial conditions.
172  if (getParam<bool>("auto_initialization"))
174 
175  // Some setup
178 
179  // Make sure all PETSc options are setup correctly
181 }
182 
183 void
185 {
186 #if PETSC_RELEASE_LESS_THAN(3, 12, 0)
187  // Make sure the SLEPc options are setup for this app
189 #else
190  // Options need to be setup once only
192  {
193  // Master app has the default data base
194  if (!_app.isUltimateMaster())
195  PetscOptionsPush(_eigen_problem.petscOptionsDatabase());
197  if (!_app.isUltimateMaster())
198  PetscOptionsPop();
200  }
201 #endif
202 }
203 
204 void
206 {
207  // check to make sure that we don't have any time kernels in eigenvalue simulation
209  mooseError("You have specified time kernels in your eigenvalue simulation");
210 }
211 
212 #endif
213 
214 void
216 {
217 #ifdef LIBMESH_HAVE_SLEPC
218  // Recovering makes sense for only transient simulations since the solution from
219  // the previous time steps is required.
220  if (_app.isRecovering())
221  {
222  _console << "\nCannot recover eigenvalue solves!\nExiting...\n" << std::endl;
223  return;
224  }
225 
226  // Outputs initial conditions set by users
227  // It is consistent with Steady
228  _time_step = 0;
229  _time = _time_step;
232 
233  preExecute();
234 
235  // The following code of this function is copied from "Steady"
236  // "Eigenvalue" implementation can be considered a one-time-step simulation to
237  // have the code compatible with the rest moose world.
239 
240  // First step in any eigenvalue state solve is always 1 (preserving backwards compatibility)
241  _time_step = 1;
242 
243 #ifdef LIBMESH_ENABLE_AMR
244 
245  // Define the refinement loop
246  auto steps = _eigen_problem.adaptivity().getSteps();
247  for (const auto r_step : make_range(steps + 1))
248  {
249 #endif // LIBMESH_ENABLE_AMR
251 
253  if (!lastSolveConverged())
254  {
255  _console << "Aborting as solve did not converge" << std::endl;
256  break;
257  }
258 
259  // Compute markers and indicators only when we do have at least one adaptivity step
260  if (steps)
261  {
264  }
265  // need to keep _time in sync with _time_step to get correct output
266  _time = _time_step;
269 
270 #ifdef LIBMESH_ENABLE_AMR
271  if (r_step < steps)
272  {
274  }
275 
276  _time_step++;
277  }
278 #endif
279 
280  {
281  TIME_SECTION(_final_timer)
286  _time = _time_step;
289  }
290 
291  postExecute();
292 
293 #else
294  mooseError("SLEPc is required for eigenvalue executioner, please use --download-slepc when "
295  "configuring PETSc ");
296 #endif
297 }
void finalizeMultiApps()
virtual void preExecute()
Override this for actions that should take place before execution.
Definition: Executioner.h:73
void timestepSetup() override
virtual bool lastSolveConverged() const override
Whether or not the last solve converged.
Definition: Eigenvalue.h:45
bool isUltimateMaster() const
Whether or not this app is the ultimate master app.
Definition: MooseApp.h:823
static InputParameters validParams()
Constructor.
Definition: Eigenvalue.C:22
SolverParams & solverParams()
Get the solver parameters.
T & getUserObject(const std::string &name, unsigned int tid=0) const
Get the user object by its name.
void precondMatrixIncludesEigenKernels(bool precond_matrix_includes_eigen)
If the preconditioning matrix includes eigen kernels.
void mooseDeprecated(Args &&... args) const
virtual void init() override
Initialize the executioner.
Definition: Eigenvalue.C:151
void addPrivateParam(const std::string &name, const T &value)
These method add a parameter to the InputParameters object which can be retrieved like any other para...
Eigenvalue executioner is used to drive the eigenvalue calculations.
Definition: Eigenvalue.h:30
virtual void postExecute()
Method called at the end of the simulation.
void initEigenvector(const Real initial_value)
For nonlinear eigen solver, a good initial value can help convergence.
Definition: EigenProblem.C:405
virtual void computeMarkers()
virtual void execute(const ExecFlagType &exec_type) override
Convenience function for performing execution of MOOSE systems.
Definition: EigenProblem.C:154
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
Real & _time
Definition: Eigenvalue.h:78
const ExecFlagType EXEC_TIMESTEP_END
Definition: Moose.C:32
bool _last_solve_converged
Definition: Eigenvalue.h:83
bool & petscOptionsInserted()
If PETSc options are already inserted.
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
virtual void advanceState()
Advance all of the state holding vectors / datastructures so that we can move to the next timestep...
InputParameters getSlepcValidParams(InputParameters &params)
Definition: SlepcSupport.C:37
void setNormalization(const PostprocessorName &pp, const Real value=std::numeric_limits< Real >::max())
Set postprocessor and normalization factor &#39;Postprocessor&#39; is often used to compute an integral of ph...
Definition: EigenProblem.C:604
Moose::EigenProblemType _eigen_problem_type
Definition: SolverParams.h:25
static InputParameters validParams()
Definition: Executioner.C:26
virtual void execute() override
Pure virtual execute function MUST be overridden by children classes.
Definition: Eigenvalue.C:215
virtual void checkIntegrity()
Eigenvalue executioner does not allow time kernels.
Definition: Eigenvalue.C:205
Eigenvalue(const InputParameters &parameters)
Definition: Eigenvalue.C:91
const ExecFlagType EXEC_PRE_MULTIAPP_SETUP
Definition: Moose.C:46
virtual void computeIndicators()
virtual void postExecute()
Override this for actions that should take place after execution.
Definition: Executioner.h:78
void initialSetup() override
Executioners are objects that do the actual work of solving your problem.
Definition: Executioner.h:30
NonlinearSystemBase & getNonlinearSystemBase(const unsigned int sys_num)
MooseApp & _app
The MOOSE application this is associated with.
Definition: MooseBase.h:69
void paramError(const std::string &param, Args... args) const
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
const ExecFlagType EXEC_LINEAR
Definition: Moose.C:29
FEProblemSolve _feproblem_solve
inner-most solve object to perform Newton solve with SLEPc
Definition: Eigenvalue.h:71
static InputParameters validParams()
registerMooseObject("MooseApp", Eigenvalue)
NonlinearEigenSystem & getNonlinearEigenSystem(const unsigned int nl_sys_num)
Definition: EigenProblem.h:302
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
int & _time_step
Definition: Eigenvalue.h:77
PetscOptions & petscOptionsDatabase()
void prepareSolverOptions()
Prepare right petsc options.
Definition: Eigenvalue.C:184
Real _system_time
Definition: Eigenvalue.h:76
unsigned int getSteps() const
Pull out the number of steps previously set by calling init()
Definition: Adaptivity.h:104
IntRange< T > make_range(T beg, T end)
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
unsigned int _free_power_iterations
Definition: SolverParams.h:31
const InputParameters & _pars
Parameters of this object, references the InputParameters stored in the InputParametersWarehouse.
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
const InputParameters & parameters() const
Get the parameters of the object.
PerfID _final_timer
Definition: Eigenvalue.h:80
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
void setInitialEigenvalue(const Real initial_eigenvalue)
Set an initial eigenvalue for initial normalization.
Definition: EigenProblem.h:72
void setEigenProblemSolverParams(EigenProblem &eigen_problem, const InputParameters &params)
Retrieve eigen problem params from &#39;params&#39;, and then set these params into SolverParams.
Definition: SlepcSupport.C:177
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
EigenProblem & _eigen_problem
Definition: Eigenvalue.h:68
bool execMultiApps(ExecFlagType type, bool auto_advance=true)
Execute the MultiApps associated with the ExecFlagType.
unsigned int _extra_power_iterations
Definition: SolverParams.h:32
Problem for solving eigenvalue problems.
Definition: EigenProblem.h:21
void storeSolveType(FEProblemBase &fe_problem, const InputParameters &params)
Set solve type into eigen problem (solverParams)
Definition: SlepcSupport.C:247
Adaptivity & adaptivity()
bool isRecovering() const
Whether or not this is a "recover" calculation.
Definition: MooseApp.C:1167
const ExecFlagType EXEC_FINAL
Definition: Moose.C:38
void setEigenproblemType(Moose::EigenProblemType eigen_problem_type)
Set eigen problem type.
Base class for user-specific data.
Definition: UserObject.h:39
std::unique_ptr< FixedPointSolve > _fixed_point_solve
Definition: Executioner.h:169
void slepcSetOptions(EigenProblem &eigen_problem, const InputParameters &params)
Push all SLEPc/PETSc options into SLEPc/PETSc side.
Definition: SlepcSupport.C:496
virtual void outputStep(ExecFlagType type)
Output the current step.
InputParameters getSlepcEigenProblemValidParams()
Retrieve valid params that allow users to specify eigen problem configuration.
Definition: SlepcSupport.C:66
virtual bool adaptMesh()
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)
This method takes a space delimited list of parameter names and adds them to the specified group name...
const ExecFlagType EXEC_INITIAL
Definition: Moose.C:28