https://mooseframework.inl.gov
EigenProblem.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
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 #include "libmesh/libmesh_config.h"
11 
12 #include "EigenProblem.h"
13 #include "Assembly.h"
14 #include "AuxiliarySystem.h"
15 #include "DisplacedProblem.h"
16 #include "NonlinearEigenSystem.h"
17 #include "SlepcSupport.h"
18 #include "RandomData.h"
19 #include "OutputWarehouse.h"
20 #include "Function.h"
21 #include "MooseVariableScalar.h"
22 #include "UserObject.h"
23 
24 // libMesh includes
25 #include "libmesh/system.h"
26 #include "libmesh/eigen_solver.h"
27 #include "libmesh/enum_eigen_solver_type.h"
28 
29 // Needed for LIBMESH_CHECK_ERR
31 
33 
36 {
38  params.addClassDescription("Problem object for solving an eigenvalue problem.");
39  params.addParam<bool>("negative_sign_eigen_kernel",
40  true,
41  "Whether or not to use a negative sign for eigenvalue kernels. "
42  "Using a negative sign makes eigenvalue kernels consistent with "
43  "a nonlinear solver");
44 
45  params.addParam<unsigned int>(
46  "active_eigen_index",
47  0,
48  "Which eigenvector is used to compute residual and also associated to nonlinear variable");
49  params.addParam<PostprocessorName>("bx_norm", "A postprocessor describing the norm of Bx");
50 
51  params.addParamNamesToGroup("negative_sign_eigen_kernel active_eigen_index bx_norm",
52  "Eigenvalue solve");
53 
54  return params;
55 }
56 
58  : FEProblemBase(parameters)
59 #ifdef LIBMESH_HAVE_SLEPC
60  ,
61  // By default, we want to compute an eigenvalue only (smallest or largest)
62  _n_eigen_pairs_required(1),
63  _generalized_eigenvalue_problem(false),
64  _negative_sign_eigen_kernel(getParam<bool>("negative_sign_eigen_kernel")),
65  _active_eigen_index(getParam<unsigned int>("active_eigen_index")),
66  _do_free_power_iteration(false),
67  _output_inverse_eigenvalue(false),
68  _on_linear_solver(false),
69  _matrices_formed(false),
70  _constant_matrices(false),
71  _has_normalization(false),
72  _normal_factor(1.0),
73  _first_solve(declareRestartableData<bool>("first_solve", true)),
74  _bx_norm_name(isParamValid("bx_norm")
75  ? std::make_optional(getParam<PostprocessorName>("bx_norm"))
76  : std::nullopt)
77 #endif
78 {
79 #ifdef LIBMESH_HAVE_SLEPC
80  if (_nl_sys_names.size() > 1)
81  paramError("nl_sys_names",
82  "eigen problems do not currently support multiple nonlinear eigen systems");
83  if (_linear_sys_names.size())
84  paramError("linear_sys_names", "EigenProblem only works with a single nonlinear eigen system");
85 
86  for (const auto i : index_range(_nl_sys_names))
87  {
88  const auto & sys_name = _nl_sys_names[i];
89  auto & nl = _nl[i];
90  nl = std::make_shared<NonlinearEigenSystem>(*this, sys_name);
92  _current_nl_sys = nl.get();
95  }
96 
97  _aux = std::make_shared<AuxiliarySystem>(*this, "aux0");
98 
100 
102 
103  es().parameters.set<EigenProblem *>("_eigen_problem") = this;
104 #else
105  mooseError("Need to install SLEPc to solve eigenvalue problems, please reconfigure\n");
106 #endif /* LIBMESH_HAVE_SLEPC */
107 
108  // SLEPc older than 3.13.0 can not take initial guess from moose
109 #if PETSC_RELEASE_LESS_THAN(3, 13, 0)
111  "Please use SLEPc-3.13.0 or higher. Old versions of SLEPc likely produce bad convergence");
112 #endif
113  // Create extra vectors and matrices if any
115 
116  // Create extra solution vectors if any
118 }
119 
120 #ifdef LIBMESH_HAVE_SLEPC
121 void
123 {
124  switch (eigen_problem_type)
125  {
127  _nl_eigen->sys().set_eigenproblem_type(libMesh::HEP);
129  break;
130 
132  _nl_eigen->sys().set_eigenproblem_type(libMesh::NHEP);
134  break;
135 
137  _nl_eigen->sys().set_eigenproblem_type(libMesh::GHEP);
139  break;
140 
142  _nl_eigen->sys().set_eigenproblem_type(libMesh::GHIEP);
144  break;
145 
147  _nl_eigen->sys().set_eigenproblem_type(libMesh::GNHEP);
149  break;
150 
152  mooseError("libMesh does not support EPT_POS_GEN_NON_HERMITIAN currently \n");
153  break;
154 
157  break;
158 
159  default:
160  mooseError("Unknown eigen solver type \n");
161  }
162 }
163 
164 void
166 {
167  if (exec_type == EXEC_INITIAL && !_app.isRestarting())
168  // we need to scale the solution properly and we can do this only all initial setup of
169  // depending objects by the residual evaluations has been done to this point.
170  preScaleEigenVector(std::pair<Real, Real>(_initial_eigenvalue, 0));
171 
172  FEProblemBase::execute(exec_type);
173 }
174 
175 void
177  SparseMatrix<Number> & jacobian,
178  TagID tag)
179 {
180  TIME_SECTION("computeJacobianTag", 3);
181 
182  // Disassociate the default tags because we will associate vectors with only the
183  // specific system tags that we need for this instance
184  _nl_eigen->disassociateDefaultMatrixTags();
185 
186  // Clear FE tags and first add the specific tag associated with the Jacobian
187  _fe_matrix_tags.clear();
188  _fe_matrix_tags.insert(tag);
189 
190  // Add any other user-added matrix tags if they have associated matrices
191  const auto & matrix_tags = getMatrixTags();
192  for (const auto & matrix_tag : matrix_tags)
193  if (_nl_eigen->hasMatrix(matrix_tag.second))
194  _fe_matrix_tags.insert(matrix_tag.second);
195 
196  _nl_eigen->setSolution(soln);
197 
198  _nl_eigen->associateMatrixToTag(jacobian, tag);
199 
202 
203  _nl_eigen->disassociateMatrixFromTag(jacobian, tag);
204 }
205 
206 void
208  const std::vector<SparseMatrix<Number> *> & jacobians,
209  const std::set<TagID> & tags)
210 {
211  TIME_SECTION("computeMatricesTags", 3);
212 
213  if (jacobians.size() != tags.size())
214  mooseError("The number of matrices ",
215  jacobians.size(),
216  " does not equal the number of tags ",
217  tags.size());
218 
219  // Disassociate the default tags because we will associate vectors with only the
220  // specific system tags that we need for this instance
221  _nl_eigen->disassociateDefaultMatrixTags();
222 
223  _fe_matrix_tags.clear();
224 
225  _nl_eigen->setSolution(soln);
226 
227  unsigned int i = 0;
228  for (auto tag : tags)
229  _nl_eigen->associateMatrixToTag(*(jacobians[i++]), tag);
230 
232  computeJacobianTags(tags);
233 
234  i = 0;
235  for (auto tag : tags)
236  _nl_eigen->disassociateMatrixFromTag(*(jacobians[i++]), tag);
237 }
238 
239 void
240 EigenProblem::computeJacobianBlocks(std::vector<JacobianBlock *> & blocks,
241  const unsigned int nl_sys_num)
242 {
243  TIME_SECTION("computeJacobianBlocks", 3);
244  setCurrentNonlinearSystem(nl_sys_num);
245 
246  if (_displaced_problem)
248 
250 
252 
253  _current_nl_sys->computeJacobianBlocks(blocks, {_nl_eigen->precondMatrixTag()});
254 
256 }
257 
258 void
260  SparseMatrix<Number> & jacobianA,
261  SparseMatrix<Number> & jacobianB,
262  TagID tagA,
263  TagID tagB)
264 {
265  TIME_SECTION("computeJacobianAB", 3);
266 
267  // Disassociate the default tags because we will associate vectors with only the
268  // specific system tags that we need for this instance
269  _nl_eigen->disassociateDefaultMatrixTags();
270 
271  // Clear FE tags and first add the specific tags associated with the Jacobian
272  _fe_matrix_tags.clear();
273  _fe_matrix_tags.insert(tagA);
274  _fe_matrix_tags.insert(tagB);
275 
276  // Add any other user-added matrix tags if they have associated matrices
277  const auto & matrix_tags = getMatrixTags();
278  for (const auto & matrix_tag : matrix_tags)
279  if (_nl_eigen->hasMatrix(matrix_tag.second))
280  _fe_matrix_tags.insert(matrix_tag.second);
281 
282  _nl_eigen->setSolution(soln);
283 
284  _nl_eigen->associateMatrixToTag(jacobianA, tagA);
285  _nl_eigen->associateMatrixToTag(jacobianB, tagB);
286 
289 
290  _nl_eigen->disassociateMatrixFromTag(jacobianA, tagA);
291  _nl_eigen->disassociateMatrixFromTag(jacobianB, tagB);
292 }
293 
294 void
296  NumericVector<Number> & residual,
297  TagID tag)
298 {
299  TIME_SECTION("computeResidualTag", 3);
300 
301  // Disassociate the default tags because we will associate vectors with only the
302  // specific system tags that we need for this instance
303  _nl_eigen->disassociateDefaultVectorTags();
304 
305  // add the specific tag associated with the residual
306  mooseAssert(_fe_vector_tags.empty(), "This should be empty indicating a clean starting state");
307  _fe_vector_tags.insert(tag);
308 
309  // Add any other user-added vector residual tags if they have associated vectors
310  const auto & residual_vector_tags = getVectorTags(Moose::VECTOR_TAG_RESIDUAL);
311  for (const auto & vector_tag : residual_vector_tags)
312  if (_nl_eigen->hasVector(vector_tag._id))
313  _fe_vector_tags.insert(vector_tag._id);
314 
315  _nl_eigen->associateVectorToTag(residual, tag);
316 
317  _nl_eigen->setSolution(soln);
318 
321  _fe_vector_tags.clear();
322 
323  _nl_eigen->disassociateVectorFromTag(residual, tag);
324 }
325 
326 void
328  NumericVector<Number> & residualA,
329  NumericVector<Number> & residualB,
330  TagID tagA,
331  TagID tagB)
332 {
333  TIME_SECTION("computeResidualAB", 3);
334 
335  // Disassociate the default tags because we will associate vectors with only the
336  // specific system tags that we need for this instance
337  _nl_eigen->disassociateDefaultVectorTags();
338 
339  // add the specific tags associated with the residual
340  mooseAssert(_fe_vector_tags.empty(), "This should be empty indicating a clean starting state");
341  _fe_vector_tags.insert(tagA);
342  _fe_vector_tags.insert(tagB);
343 
344  // Add any other user-added vector residual tags if they have associated vectors
345  const auto & residual_vector_tags = getVectorTags(Moose::VECTOR_TAG_RESIDUAL);
346  for (const auto & vector_tag : residual_vector_tags)
347  if (_nl_eigen->hasVector(vector_tag._id))
348  _fe_vector_tags.insert(vector_tag._id);
349 
350  _nl_eigen->associateVectorToTag(residualA, tagA);
351  _nl_eigen->associateVectorToTag(residualB, tagB);
352 
353  _nl_eigen->setSolution(soln);
354 
356  _fe_vector_tags.clear();
357 
358  _nl_eigen->disassociateVectorFromTag(residualA, tagA);
359  _nl_eigen->disassociateVectorFromTag(residualB, tagB);
360 }
361 
362 Real
364 {
365  computeResidualAB(*_nl_eigen->currentSolution(),
366  _nl_eigen->residualVectorAX(),
367  _nl_eigen->residualVectorBX(),
368  _nl_eigen->nonEigenVectorTag(),
369  _nl_eigen->eigenVectorTag());
370 
371  Real eigenvalue = 1.0;
372 
373  if (_active_eigen_index < _nl_eigen->getNumConvergedEigenvalues())
374  eigenvalue = _nl_eigen->getConvergedEigenvalue(_active_eigen_index).first;
375 
376  // Scale BX with eigenvalue
377  _nl_eigen->residualVectorBX() *= eigenvalue;
378 
379  // Compute entire residual
381  _nl_eigen->residualVectorAX() += _nl_eigen->residualVectorBX();
382  else
383  _nl_eigen->residualVectorAX() -= _nl_eigen->residualVectorBX();
384 
385  return _nl_eigen->residualVectorAX().l2_norm();
386 }
387 
388 void
389 EigenProblem::adjustEigenVector(const Real value, bool scaling)
390 {
391  std::vector<VariableName> var_names = getVariableNames();
392  for (auto & vn : var_names)
393  {
394  MooseVariableBase * var = nullptr;
395  if (hasScalarVariable(vn))
396  var = &getScalarVariable(0, vn);
397  else
398  var = &getVariable(0, vn);
399  // Do operations for only eigen variable
400  if (var->eigen())
401  for (unsigned int vc = 0; vc < var->count(); ++vc)
402  {
403  std::set<dof_id_type> var_indices;
404  _nl_eigen->system().local_dof_indices(var->number() + vc, var_indices);
405  for (const auto & dof : var_indices)
406  _nl_eigen->solution().set(dof, scaling ? (_nl_eigen->solution()(dof) * value) : value);
407  }
408  }
409 
410  _nl_eigen->solution().close();
411  _nl_eigen->update();
412 }
413 
414 void
415 EigenProblem::scaleEigenvector(const Real scaling_factor)
416 {
417  adjustEigenVector(scaling_factor, true);
418 }
419 
420 void
421 EigenProblem::initEigenvector(const Real initial_value)
422 {
423  // Yaqi's note: the following code will set a flat solution for lagrange and
424  // constant monomial variables. For the first or higher order elemental variables,
425  // the solution is not flat. Fortunately, the initial guess does not affect
426  // the final solution as long as it is not perpendicular to the true solution.
427  // We, in general, do not need to worry about that.
428 
430 }
431 
432 void
433 EigenProblem::preScaleEigenVector(const std::pair<Real, Real> & eig)
434 {
435  // pre-scale the solution to make sure ||Bx||_2 is equal to inverse of eigenvalue
437  *_nl_eigen->currentSolution(), _nl_eigen->residualVectorBX(), _nl_eigen->eigenVectorTag());
438 
439  // Eigenvalue magnitude
440  Real v = std::sqrt(eig.first * eig.first + eig.second * eig.second);
441  // Scaling factor
442  Real factor = 1 / v / (bxNormProvided() ? formNorm() : _nl_eigen->residualVectorBX().l2_norm());
443  // Scale eigenvector
444  if (!MooseUtils::absoluteFuzzyEqual(factor, 1))
445  scaleEigenvector(factor);
446 }
447 
448 void
450 {
451  if (_has_normalization)
452  {
453  Real v;
455  {
456  if (_active_eigen_index >= _nl_eigen->getNumConvergedEigenvalues())
457  mooseError("Number of converged eigenvalues ",
458  _nl_eigen->getNumConvergedEigenvalues(),
459  " but you required eigenvalue ",
461 
462  // when normal factor is not provided, we use the inverse of the norm of
463  // the active eigenvalue for normalization
464  auto eig = _nl_eigen->getAllConvergedEigenvalues()[_active_eigen_index];
465  v = 1 / std::sqrt(eig.first * eig.first + eig.second * eig.second);
466  }
467  else
468  v = _normal_factor;
469 
471 
472  // We scale SLEPc eigen vector here, so we need to scale it back for optimal
473  // convergence if we call EPS solver again
474  mooseAssert(v != 0., "normal factor can not be zero");
475 
476  unsigned int itr = 0;
477 
478  while (!MooseUtils::relativeFuzzyEqual(v, c))
479  {
480  // If postprocessor is not defined on eigen variables, scaling might not work
481  if (itr > 10)
482  mooseError("Can not scale eigenvector to the required factor ",
483  v,
484  " please check if postprocessor is defined on only eigen variables");
485 
486  mooseAssert(c != 0., "postprocessor value used for scaling can not be zero");
487 
488  scaleEigenvector(v / c);
489 
490  // update all aux variables and user objects on linear
492 
494 
495  itr++;
496  }
497  }
498 }
499 
500 void
502 {
504  _nl_eigen->checkIntegrity();
505  if (_bx_norm_name)
506  {
508  paramWarning("bx_norm", "This parameter is only used for nonlinear solve types");
509  else if (auto & pp = getUserObjectBase(_bx_norm_name.value());
511  pp.paramError("execute_on",
512  "If providing the Bx norm, this postprocessor must execute on linear e.g. "
513  "during residual evaluations");
514  }
515 }
516 
517 void
518 EigenProblem::doFreeNonlinearPowerIterations(unsigned int free_power_iterations)
519 {
520  mooseAssert(_current_nl_sys, "This needs to be non-null");
521 
522  doFreePowerIteration(true);
523  // Set free power iterations
525 
526  // Call solver
529 
530  // Clear free power iterations
531  auto executioner = getMooseApp().getExecutioner();
532  if (executioner)
534  else
535  mooseError("There is no executioner for this moose app");
536 
537  doFreePowerIteration(false);
538 }
539 
540 void
541 EigenProblem::solve(const unsigned int nl_sys_num)
542 {
543 #if !PETSC_RELEASE_LESS_THAN(3, 12, 0)
544  // Master has the default database
545  if (!_app.isUltimateMaster())
546  LibmeshPetscCall(PetscOptionsPush(_petsc_option_data_base));
547 #endif
548 
549  setCurrentNonlinearSystem(nl_sys_num);
550 
551  if (_solve)
552  {
553  TIME_SECTION("solve", 1);
554 
555  // Set necessary slepc callbacks
556  // We delay this call as much as possible because libmesh
557  // could rebuild matrices due to mesh changes or something else.
558  _nl_eigen->attachSLEPcCallbacks();
559 
560  // If there is an eigenvalue, we scale 1/|Bx| to eigenvalue
561  if (_active_eigen_index < _nl_eigen->getNumConvergedEigenvalues())
562  {
563  std::pair<Real, Real> eig = _nl_eigen->getConvergedEigenvalue(_active_eigen_index);
564  preScaleEigenVector(eig);
565  }
566 
567  if (isNonlinearEigenvalueSolver(nl_sys_num) &&
568  solverParams(nl_sys_num)._eigen_solve_type != Moose::EST_NONLINEAR_POWER)
569  {
570  // Let do an initial solve if a nonlinear eigen solver but not power is used.
571  // The initial solver is a Inverse Power, and it is used to compute a good initial
572  // guess for Newton
573  if (solverParams(nl_sys_num)._free_power_iterations && _first_solve)
574  {
575  _console << std::endl << " -------------------------------" << std::endl;
576  _console << " Free power iteration starts ..." << std::endl;
577  _console << " -------------------------------" << std::endl << std::endl;
578  doFreeNonlinearPowerIterations(solverParams(nl_sys_num)._free_power_iterations);
579  _first_solve = false;
580  }
581 
582  // Let us do extra power iterations here if necessary
583  if (solverParams(nl_sys_num)._extra_power_iterations)
584  {
585  _console << std::endl << " --------------------------------------" << std::endl;
586  _console << " Extra Free power iteration starts ..." << std::endl;
587  _console << " --------------------------------------" << std::endl << std::endl;
588  doFreeNonlinearPowerIterations(solverParams(nl_sys_num)._extra_power_iterations);
589  }
590  }
591 
592  // We print this for only nonlinear solver
593  if (isNonlinearEigenvalueSolver(nl_sys_num))
594  {
595  _console << std::endl << " -------------------------------------" << std::endl;
596 
597  if (solverParams(nl_sys_num)._eigen_solve_type != Moose::EST_NONLINEAR_POWER)
598  _console << " Nonlinear Newton iteration starts ..." << std::endl;
599  else
600  _console << " Nonlinear power iteration starts ..." << std::endl;
601 
602  _console << " -------------------------------------" << std::endl << std::endl;
603  }
604 
607 
608  // with PJFNKMO solve type, we need to evaluate the linear objects to bring them up-to-date
609  if (solverParams(nl_sys_num)._eigen_solve_type == Moose::EST_PJFNKMO)
611 
612  // Scale eigen vector if users ask
614  }
615 
616 #if !PETSC_RELEASE_LESS_THAN(3, 12, 0)
617  if (!_app.isUltimateMaster())
618  LibmeshPetscCall(PetscOptionsPop());
619 #endif
620 
621  // sync solutions in displaced problem
622  if (_displaced_problem)
623  _displaced_problem->syncSolutions();
624 
625  // Reset the matrix flag, so that we reform matrix in next picard iteration
626  _matrices_formed = false;
627 }
628 
629 void
630 EigenProblem::setNormalization(const PostprocessorName & pp, const Real value)
631 {
632  _has_normalization = true;
633  _normalization = pp;
635 }
636 
637 void
639 {
640 #if PETSC_RELEASE_LESS_THAN(3, 13, 0)
641  // Prior to Slepc 3.13 we did not have a nonlinear eigenvalue solver so we must always assemble
642  // before the solve
643  _nl_eigen->sys().attach_assemble_function(Moose::assemble_matrix);
644 #else
645  mooseAssert(
646  numNonlinearSystems() == 1,
647  "We should have errored during construction if we had more than one nonlinear system");
648  mooseAssert(numLinearSystems() == 0,
649  "We should have errored during construction if we had any linear systems");
651  // We don't need to assemble before the solve
652  _nl_eigen->sys().assemble_before_solve = false;
653  else
654  _nl_eigen->sys().attach_assemble_function(Moose::assemble_matrix);
655 
656  // If matrix_free=true, this tells Libmesh to use shell matrices
657  _nl_eigen->sys().use_shell_matrices(solverParams(0)._eigen_matrix_free &&
658  !solverParams(0)._eigen_matrix_vector_mult);
659  // We need to tell libMesh if we are using a shell preconditioning matrix
660  _nl_eigen->sys().use_shell_precond_matrix(solverParams(0)._precond_matrix_free);
661 #endif
662 
664 }
665 
666 bool
668 {
669  if (_solve)
670  return _nl_eigen->converged();
671  else
672  return true;
673 }
674 
675 bool
676 EigenProblem::isNonlinearEigenvalueSolver(const unsigned int eigen_sys_num) const
677 {
678  const auto & solver_params = solverParams(eigen_sys_num);
679  return solver_params._eigen_solve_type == Moose::EST_NONLINEAR_POWER ||
680  solver_params._eigen_solve_type == Moose::EST_NEWTON ||
681  solver_params._eigen_solve_type == Moose::EST_PJFNK ||
682  solver_params._eigen_solve_type == Moose::EST_JFNK ||
683  solver_params._eigen_solve_type == Moose::EST_PJFNKMO;
684 }
685 
686 void
688 {
690 }
691 
692 Real
694 {
695  mooseAssert(_bx_norm_name,
696  "We should not get here unless a bx_norm postprocessor has been provided");
698 }
699 #endif
700 
701 std::string
702 EigenProblem::solverTypeString(const unsigned int solver_sys_num)
703 {
704  return Moose::stringify(solverParams(solver_sys_num)._eigen_solve_type);
705 }
Nonlinear eigenvalue system to be solved.
Generalized Non-Hermitian.
Definition: MooseTypes.h:876
virtual void computeResidualTag(const NumericVector< Number > &soln, NumericVector< Number > &residual, TagID tag) override
Form a vector for all kernels and BCs with a given tag.
Definition: EigenProblem.C:295
static InputParameters validParams()
Definition: EigenProblem.C:35
void computeJacobianBlocks(std::vector< JacobianBlock *> &blocks)
Computes several Jacobian blocks simultaneously, summing their contributions into smaller preconditio...
bool _matrices_formed
Whether or not matrices had been formed.
Definition: EigenProblem.h:260
Generalized Hermitian indefinite.
Definition: MooseTypes.h:875
Newton-based eigensolver with an assembled Jacobian matrix (fully coupled by default) ...
Definition: MooseTypes.h:861
bool isUltimateMaster() const
Whether or not this app is the ultimate master app.
Definition: MooseApp.h:837
virtual void initPetscOutputAndSomeSolverSettings() override
Hook up monitors for SNES and KSP.
Definition: EigenProblem.C:687
PostprocessorName _normalization
Postprocessor used to compute a factor from eigenvector.
Definition: EigenProblem.h:266
bool absoluteFuzzyEqual(const T &var1, const T2 &var2, const T3 &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Function to check whether two variables are equal within an absolute tolerance.
Definition: MooseUtils.h:380
void mooseDeprecated(Args &&... args) const
virtual void computeJacobianBlocks(std::vector< JacobianBlock *> &blocks, const unsigned int nl_sys_num) override
Computes several Jacobian blocks simultaneously, summing their contributions into smaller preconditio...
Definition: EigenProblem.C:240
unsigned int TagID
Definition: MooseTypes.h:210
virtual std::size_t numNonlinearSystems() const override
char ** blocks
unsigned int number() const
Get variable number coming from libMesh.
bool _currently_computing_jacobian
Flag to determine whether the problem is currently computing Jacobian.
Definition: SubProblem.h:1096
virtual void init() override
The same as PJFNK except that matrix-vector multiplication is employed to replace residual evaluation...
Definition: MooseTypes.h:863
void initEigenvector(const Real initial_value)
For nonlinear eigen solver, a good initial value can help convergence.
Definition: EigenProblem.C:421
void clearFreeNonlinearPowerIterations(const InputParameters &params)
Definition: SlepcSupport.C:419
virtual void newAssemblyArray(std::vector< std::shared_ptr< SolverSystem >> &solver_systems)
bool _has_normalization
Whether or not we normalize eigenvector.
Definition: EigenProblem.h:264
EigenProblem(const InputParameters &parameters)
Definition: EigenProblem.C:57
unsigned int count() const
Get the number of components Note: For standard and vector variables, the number is one...
virtual libMesh::System & system()=0
Get the reference to the libMesh system.
virtual void execute(const ExecFlagType &exec_type) override
Convenience function for performing execution of MOOSE systems.
Definition: EigenProblem.C:165
PetscOptions _petsc_option_data_base
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
void adjustEigenVector(const Real value, bool scaling)
Adjust eigen vector by either scaling the existing values or setting new values The operations are ap...
Definition: EigenProblem.C:389
void doFreeNonlinearPowerIterations(unsigned int free_power_iterations)
Do some free/extra power iterations.
Definition: EigenProblem.C:518
virtual bool hasScalarVariable(const std::string &var_name) const override
Returns a Boolean indicating whether any system contains a variable with the name provided...
std::vector< std::shared_ptr< SolverSystem > > _solver_systems
Combined container to base pointer of every solver system.
std::unique_ptr< T_DEST, T_DELETER > dynamic_pointer_cast(std::unique_ptr< T_SRC, T_DELETER > &src)
These are reworked from https://stackoverflow.com/a/11003103.
use whatever SLPEC has by default
Definition: MooseTypes.h:878
bool _negative_sign_eigen_kernel
Whether or not use negative sign for Bx.
Definition: EigenProblem.h:246
bool isRestarting() const
Whether or not this is a "restart" calculation.
Definition: MooseApp.C:1801
Generalized Non-Hermitian with positive (semi-)definite B.
Definition: MooseTypes.h:877
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
void preScaleEigenVector(const std::pair< Real, Real > &eig)
Eigenvector need to be scaled back if it was scaled in an earlier stage Scaling eigen vector does not...
Definition: EigenProblem.C:433
MooseApp & getMooseApp() const
Get the MooseApp this class is associated with.
Definition: MooseBase.h:45
bool isNonlinearEigenvalueSolver(unsigned int eigen_sys_num) const
Definition: EigenProblem.C:676
virtual const MooseVariableFieldBase & getVariable(const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type=Moose::VarKindType::VAR_ANY, Moose::VarFieldType expected_var_field_type=Moose::VarFieldType::VAR_FIELD_ANY) const override
Returns the variable reference for requested variable which must be of the expected_var_type (Nonline...
virtual void computeResidualTags(const std::set< TagID > &tags)
Form multiple residual vectors and each is associated with one tag.
auto max(const L &left, const R &right)
std::shared_ptr< NonlinearEigenSystem > _nl_eigen
Definition: EigenProblem.h:242
bool eigen() const
Whether or not this variable operates on an eigen kernel.
bool _generalized_eigenvalue_problem
Definition: EigenProblem.h:241
void update()
Update the system (doing libMesh magic)
Definition: SystemBase.C:1245
Real formNorm()
Form the Bx norm.
Definition: EigenProblem.C:693
bool & _first_solve
A flag to indicate if it is the first time calling the solve.
Definition: EigenProblem.h:271
std::vector< std::shared_ptr< NonlinearSystemBase > > _nl
The nonlinear systems.
void createTagSolutions()
Create extra tagged solution vectors.
bool contains(const std::string &value) const
Methods for seeing if a value is set in the MultiMooseEnum.
virtual void execute(const ExecFlagType &exec_type)
Convenience function for performing execution of MOOSE systems.
virtual std::vector< VariableName > getVariableNames()
Returns a list of all the variables in the problem (both from the NL and Aux systems.
virtual Real computeResidualL2Norm() override
Compute the residual of Ax - Bx.
Definition: EigenProblem.C:363
Preconditioned Jacobian-free Newton Krylov.
Definition: MooseTypes.h:862
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:630
const bool & _solve
Whether or not to actually solve the nonlinear system.
NonlinearSystemBase * _current_nl_sys
The current nonlinear system that we are solving.
void setCurrentNonlinearSystem(const unsigned int nl_sys_num)
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
void createTagVectors()
Create extra tagged vectors and matrices.
virtual void init() override
Definition: EigenProblem.C:638
const ExecFlagEnum & getExecuteOnEnum() const
Return the execute on MultiMooseEnum for this object.
virtual libMesh::EquationSystems & es() override
std::shared_ptr< AuxiliarySystem > _aux
The auxiliary system.
void computeResidualAB(const NumericVector< Number > &soln, NumericVector< Number > &residualA, NumericVector< Number > &residualB, TagID tagA, TagID tagB)
Form two vetors, where each is associated with one tag, through one element-loop. ...
Definition: EigenProblem.C:327
virtual MooseVariableScalar & getScalarVariable(const THREAD_ID tid, const std::string &var_name) override
Returns the scalar variable reference from whichever system contains it.
std::vector< VectorTag > getVectorTags(const std::set< TagID > &tag_ids) const
Definition: SubProblem.C:172
registerMooseObject("MooseApp", EigenProblem)
MooseApp & _app
The MOOSE application this is associated with.
Definition: MooseBase.h:84
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 ...
void assemble_matrix(EquationSystems &es, const std::string &system_name)
const ExecFlagType EXEC_LINEAR
Definition: Moose.C:29
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:64
virtual void checkProblemIntegrity()
Method called to perform a series of sanity checks before a simulation is run.
void computeMatricesTags(const NumericVector< Number > &soln, const std::vector< SparseMatrix< Number > *> &jacobians, const std::set< TagID > &tags)
Form several matrices simultaneously.
Definition: EigenProblem.C:207
Nonlinear inverse power.
Definition: MooseTypes.h:860
Real _initial_eigenvalue
A value used for initial normalization.
Definition: EigenProblem.h:273
unsigned int _active_eigen_index
Which eigenvalue is used to compute residual.
Definition: EigenProblem.h:249
const PostprocessorValue & getPostprocessorValueByName(const PostprocessorName &name, std::size_t t_index=0) const
Get a read-only reference to the value associated with a Postprocessor that exists.
virtual std::map< TagName, TagID > & getMatrixTags()
Return all matrix tags in the system, where a tag is represented by a map from name to ID...
Definition: SubProblem.h:253
virtual void initNullSpaceVectors(const InputParameters &parameters, std::vector< std::shared_ptr< NonlinearSystemBase >> &nl)
Executioner * getExecutioner() const
Retrieve the Executioner for this App.
Definition: MooseApp.C:2111
virtual void solve() override=0
Solve the system (using libMesh magic)
void solveSetup()
Calls the timestepSetup function for each of the output objects.
const ExecFlagType EXEC_PRE_DISPLACE
Definition: Moose.C:50
const std::vector< NonlinearSystemName > _nl_sys_names
The nonlinear system names.
const ExecFlagType EXEC_NONLINEAR
Definition: Moose.C:31
void prefer_hash_table_matrix_assembly(bool preference)
std::set< TagID > _fe_matrix_tags
void scaleEigenvector(const Real scaling_factor)
Scale eigenvector.
Definition: EigenProblem.C:415
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Class for containing MooseEnum item information.
Definition: MooseEnumItem.h:18
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template * sqrt(_arg)) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(tanh
void setFreeNonlinearPowerIterations(unsigned int free_power_iterations)
Set SLEPc/PETSc options to trigger free power iteration.
Definition: SlepcSupport.C:404
void postScaleEigenVector()
Normalize eigen vector.
Definition: EigenProblem.C:449
bool relativeFuzzyEqual(const T &var1, const T2 &var2, const T3 &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Function to check whether two variables are equal within a relative tolerance.
Definition: MooseUtils.h:492
T & set(const std::string &)
Jacobian-free Newton Krylov.
Definition: MooseTypes.h:864
const std::vector< LinearSystemName > _linear_sys_names
The linear system names.
Non-Hermitian.
Definition: MooseTypes.h:873
virtual std::string solverTypeString(unsigned int solver_sys_num=0) override
Return solver type as a human readable string.
Definition: EigenProblem.C:702
std::set< TagID > _fe_vector_tags
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
SolverParams & solverParams(unsigned int solver_sys_num=0)
Get the solver parameters.
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...
std::shared_ptr< DisplacedProblem > _displaced_problem
const InputParameters & parameters() const
Get the parameters of the object.
Real _normal_factor
Postprocessor target value.
Definition: EigenProblem.h:269
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an optional parameter and a documentation string to the InputParameters object...
static InputParameters validParams()
virtual std::size_t numLinearSystems() const override
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
const bool _use_hash_table_matrix_assembly
Whether to assemble matrices using hash tables instead of preallocating matrix memory.
std::optional< PostprocessorName > _bx_norm_name
The name of the Postprocessor providing the Bx norm.
Definition: EigenProblem.h:289
void paramWarning(const std::string &param, Args... args) const
Emits a warning prefixed with the file and line number of the given param (from the input file) along...
const UserObject & getUserObjectBase(const std::string &name, const THREAD_ID tid=0) const
Get the user object by its name.
Problem for solving eigenvalue problems.
Definition: EigenProblem.h:21
bool bxNormProvided() const
Whether a Bx norm postprocessor has been provided.
Definition: EigenProblem.h:232
virtual void checkProblemIntegrity() override
Method called to perform a series of sanity checks before a simulation is run.
Definition: EigenProblem.C:501
EigenProblemType
Type of the eigen problem.
Definition: MooseTypes.h:870
bool doFreePowerIteration() const
Whether or not we are doing free power iteration.
Definition: EigenProblem.h:79
void setEigenproblemType(Moose::EigenProblemType eigen_problem_type)
Set eigen problem type.
Number initial_value(const Point &, const Parameters &, const std::string &, const std::string &)
void ErrorVector unsigned int
virtual void solve(const unsigned int nl_sys_num) override
Definition: EigenProblem.C:541
auto index_range(const T &sizable)
Generalized Hermitian.
Definition: MooseTypes.h:874
OutputWarehouse & getOutputWarehouse()
Get the OutputWarehouse objects.
Definition: MooseApp.C:2407
void computeJacobianAB(const NumericVector< Number > &soln, SparseMatrix< Number > &jacobianA, SparseMatrix< Number > &jacobianB, TagID tagA, TagID tagB)
Form two Jacobian matrices, where each is associated with one tag, through one element-loop.
Definition: EigenProblem.C:259
void computeSystems(const ExecFlagType &type)
Do generic system computations.
virtual void computeJacobianTags(const std::set< TagID > &tags)
Form multiple matrices, and each is associated with a tag.
virtual void computeJacobianTag(const NumericVector< Number > &soln, SparseMatrix< Number > &jacobian, TagID tag) override
Form a Jacobian matrix for all kernels and BCs with a given tag.
Definition: EigenProblem.C:176
virtual bool solverSystemConverged(const unsigned int solver_sys_num) override
Definition: EigenProblem.C:667
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