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