Line data Source code
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 "NonlinearEigenSystem.h"
11 :
12 : // MOOSE includes
13 : #include "DirichletBC.h"
14 : #include "EigenDirichletBC.h"
15 : #include "ArrayDirichletBC.h"
16 : #include "EigenArrayDirichletBC.h"
17 : #include "EigenProblem.h"
18 : #include "IntegratedBC.h"
19 : #include "KernelBase.h"
20 : #include "NodalBC.h"
21 : #include "TimeIntegrator.h"
22 : #include "SlepcSupport.h"
23 : #include "DGKernelBase.h"
24 : #include "ScalarKernelBase.h"
25 : #include "MooseVariableScalar.h"
26 : #include "ResidualObject.h"
27 :
28 : #include "libmesh/libmesh_config.h"
29 : #include "libmesh/petsc_matrix.h"
30 : #include "libmesh/sparse_matrix.h"
31 : #include "libmesh/diagonal_matrix.h"
32 : #include "libmesh/petsc_shell_matrix.h"
33 : #include "libmesh/petsc_solver_exception.h"
34 : #include "libmesh/slepc_eigen_solver.h"
35 :
36 : using namespace libMesh;
37 :
38 : #ifdef LIBMESH_HAVE_SLEPC
39 :
40 : namespace Moose
41 : {
42 :
43 : void
44 55 : assemble_matrix(EquationSystems & es, const std::string & system_name)
45 : {
46 55 : EigenProblem * p = es.parameters.get<EigenProblem *>("_eigen_problem");
47 55 : CondensedEigenSystem & eigen_system = es.get_system<CondensedEigenSystem>(system_name);
48 : NonlinearEigenSystem & eigen_nl =
49 55 : p->getNonlinearEigenSystem(/*nl_sys_num=*/eigen_system.number());
50 :
51 : // If this is a nonlinear eigenvalue problem,
52 : // we do not need to assemble anything
53 55 : if (p->isNonlinearEigenvalueSolver(eigen_nl.number()))
54 : {
55 : // If you want an efficient eigensolver,
56 : // please use PETSc 3.13 or newer.
57 : // We need to do an unnecessary assembly,
58 : // if you use PETSc that is older than 3.13.
59 : #if PETSC_RELEASE_LESS_THAN(3, 13, 0)
60 : if (eigen_system.has_matrix_B())
61 : p->computeJacobianTag(*eigen_system.current_local_solution,
62 : eigen_system.get_matrix_B(),
63 : eigen_nl.eigenMatrixTag());
64 : #endif
65 0 : return;
66 : }
67 :
68 : #if !PETSC_RELEASE_LESS_THAN(3, 13, 0)
69 : // If we use shell matrices and do not use a shell preconditioning matrix,
70 : // we only need to form a preconditioning matrix
71 55 : if (eigen_system.use_shell_matrices() && !eigen_system.use_shell_precond_matrix())
72 : {
73 0 : p->computeJacobianTag(*eigen_system.current_local_solution,
74 : eigen_system.get_precond_matrix(),
75 : eigen_nl.precondMatrixTag());
76 0 : return;
77 : }
78 : #endif
79 : // If it is a linear generalized eigenvalue problem,
80 : // we assemble A and B together
81 55 : if (eigen_system.generalized())
82 : {
83 22 : p->computeJacobianAB(*eigen_system.current_local_solution,
84 : eigen_system.get_matrix_A(),
85 : eigen_system.get_matrix_B(),
86 : eigen_nl.nonEigenMatrixTag(),
87 : eigen_nl.eigenMatrixTag());
88 : #if LIBMESH_HAVE_SLEPC
89 22 : if (p->negativeSignEigenKernel())
90 22 : LibmeshPetscCallA(
91 : p->comm().get(),
92 : MatScale(static_cast<PetscMatrix<Number> &>(eigen_system.get_matrix_B()).mat(), -1.0));
93 : #endif
94 22 : return;
95 : }
96 :
97 : // If it is a linear eigenvalue problem, we assemble matrix A
98 : {
99 33 : p->computeJacobianTag(*eigen_system.current_local_solution,
100 : eigen_system.get_matrix_A(),
101 : eigen_nl.nonEigenMatrixTag());
102 :
103 33 : return;
104 : }
105 : }
106 : }
107 :
108 565 : NonlinearEigenSystem::NonlinearEigenSystem(EigenProblem & eigen_problem, const std::string & name)
109 : : NonlinearSystemBase(
110 565 : eigen_problem, eigen_problem.es().add_system<CondensedEigenSystem>(name), name),
111 565 : _eigen_sys(eigen_problem.es().get_system<CondensedEigenSystem>(name)),
112 565 : _eigen_problem(eigen_problem),
113 565 : _solver_configuration(nullptr),
114 565 : _n_eigen_pairs_required(eigen_problem.getNEigenPairsRequired()),
115 1130 : _work_rhs_vector_AX(addVector("work_rhs_vector_Ax", false, PARALLEL)),
116 1130 : _work_rhs_vector_BX(addVector("work_rhs_vector_Bx", false, PARALLEL)),
117 565 : _precond_matrix_includes_eigen(false),
118 565 : _preconditioner(nullptr),
119 1695 : _num_constrained_dofs(0)
120 : {
121 : SlepcEigenSolver<Number> * solver =
122 565 : cast_ptr<SlepcEigenSolver<Number> *>(_eigen_sys.eigen_solver.get());
123 :
124 565 : if (!solver)
125 0 : mooseError("A slepc eigen solver is required");
126 :
127 : // setup of our class @SlepcSolverConfiguration
128 : _solver_configuration =
129 565 : std::make_unique<SlepcEigenSolverConfiguration>(eigen_problem, *solver, *this);
130 :
131 565 : solver->set_solver_configuration(*_solver_configuration);
132 :
133 565 : _Ax_tag = eigen_problem.addVectorTag("Ax_tag");
134 :
135 565 : _Bx_tag = eigen_problem.addVectorTag("Eigen");
136 :
137 565 : _A_tag = eigen_problem.addMatrixTag("A_tag");
138 :
139 565 : _B_tag = eigen_problem.addMatrixTag("Eigen");
140 :
141 : // By default, _precond_tag and _A_tag will share the same
142 : // objects. If we want to include eigen contributions to
143 : // the preconditioning matrix, and then _precond_tag will
144 : // point to part of "B" objects
145 565 : _precond_tag = eigen_problem.addMatrixTag("Eigen_precond");
146 :
147 : // We do not rely on creating submatrices in the solve routine
148 565 : _eigen_sys.dont_create_submatrices_in_solve();
149 565 : }
150 :
151 : void
152 3177 : NonlinearEigenSystem::postAddResidualObject(ResidualObject & object)
153 : {
154 : // If it is an eigen dirichlet boundary condition, we should skip it because their
155 : // contributions should be zero. If we do not skip it, preconditioning matrix will
156 : // be singular because boundary elements are zero.
157 3293 : if (_precond_matrix_includes_eigen && !dynamic_cast<EigenDirichletBC *>(&object) &&
158 116 : !dynamic_cast<EigenArrayDirichletBC *>(&object))
159 116 : object.useMatrixTag(_precond_tag, {});
160 :
161 3177 : auto & vtags = object.getVectorTags({});
162 3177 : auto & mtags = object.getMatrixTags({});
163 :
164 3177 : const bool eigen = (vtags.find(_Bx_tag) != vtags.end()) || (mtags.find(_B_tag) != mtags.end());
165 :
166 3177 : if (eigen && !_eigen_sys.generalized())
167 3 : object.mooseError("This object has been marked as contributing to B or Bx but the eigen "
168 : "problem type is not a generalized one");
169 :
170 : // If it is an eigen kernel, mark its variable as eigen
171 3174 : if (eigen)
172 : {
173 : // Note: the object may be on the displaced system
174 1344 : auto sys = object.parameters().get<SystemBase *>("_sys");
175 1344 : auto vname = object.variable().name();
176 1344 : if (hasScalarVariable(vname))
177 0 : sys->getScalarVariable(0, vname).eigen(true);
178 : else
179 1344 : sys->getVariable(0, vname).eigen(true);
180 :
181 : // Associate the eigen matrix tag and the vector tag
182 : // if this is a eigen kernel
183 1344 : object.useMatrixTag(_B_tag, {});
184 1344 : object.useVectorTag(_Bx_tag, {});
185 1344 : }
186 : else
187 : {
188 : // Noneigen Vector tag
189 1830 : object.useVectorTag(_Ax_tag, {});
190 : // Noneigen Matrix tag
191 1830 : object.useMatrixTag(_A_tag, {});
192 : // Noneigen Kernels
193 1830 : object.useMatrixTag(_precond_tag, {});
194 : }
195 3174 : }
196 :
197 : void
198 559 : NonlinearEigenSystem::initializeCondensedMatrices()
199 : {
200 559 : if (!(_num_constrained_dofs = dofMap().n_constrained_dofs()))
201 535 : return;
202 :
203 24 : _eigen_sys.initialize_condensed_dofs();
204 24 : const auto m = cast_int<numeric_index_type>(_eigen_sys.local_non_condensed_dofs_vector.size());
205 24 : auto M = m;
206 24 : _communicator.sum(M);
207 24 : if (_eigen_sys.has_condensed_matrix_A())
208 : {
209 12 : _eigen_sys.get_condensed_matrix_A().init(M, M, m, m);
210 : // A bit ludicrously MatCopy requires the matrix being copied to to be assembled
211 12 : _eigen_sys.get_condensed_matrix_A().close();
212 : }
213 24 : if (_eigen_sys.has_condensed_matrix_B())
214 : {
215 12 : _eigen_sys.get_condensed_matrix_B().init(M, M, m, m);
216 12 : _eigen_sys.get_condensed_matrix_B().close();
217 : }
218 24 : if (_eigen_sys.has_condensed_precond_matrix())
219 : {
220 12 : _eigen_sys.get_condensed_precond_matrix().init(M, M, m, m);
221 12 : _eigen_sys.get_condensed_precond_matrix().close();
222 : }
223 : }
224 :
225 : void
226 559 : NonlinearEigenSystem::postInit()
227 : {
228 559 : NonlinearSystemBase::postInit();
229 559 : initializeCondensedMatrices();
230 559 : }
231 :
232 : void
233 0 : NonlinearEigenSystem::reinit()
234 : {
235 0 : NonlinearSystemBase::reinit();
236 0 : initializeCondensedMatrices();
237 0 : }
238 :
239 : void
240 1244 : NonlinearEigenSystem::solve()
241 : {
242 1244 : const bool presolve_succeeded = preSolve();
243 1244 : if (!presolve_succeeded)
244 0 : return;
245 :
246 1244 : std::unique_ptr<NumericVector<Number>> subvec;
247 :
248 : // We apply initial guess for only nonlinear solver
249 1244 : if (_eigen_problem.isNonlinearEigenvalueSolver(number()))
250 : {
251 1189 : if (_num_constrained_dofs)
252 : {
253 44 : subvec = solution().get_subvector(_eigen_sys.local_non_condensed_dofs_vector);
254 44 : _eigen_sys.set_initial_space(*subvec);
255 : }
256 : else
257 1145 : _eigen_sys.set_initial_space(solution());
258 : }
259 :
260 1244 : const bool time_integrator_solve = std::any_of(_time_integrators.begin(),
261 : _time_integrators.end(),
262 0 : [](auto & ti) { return ti->overridesSolve(); });
263 : if (time_integrator_solve)
264 : mooseAssert(_time_integrators.size() == 1,
265 : "If solve is overridden, then there must be only one time integrator");
266 :
267 1244 : if (time_integrator_solve)
268 0 : _time_integrators.front()->solve();
269 : else
270 1244 : system().solve();
271 :
272 1244 : for (auto & ti : _time_integrators)
273 : {
274 0 : if (!ti->overridesSolve())
275 0 : ti->setNumIterationsLastSolve();
276 0 : ti->postSolve();
277 : }
278 :
279 : // store solve information
280 1244 : if (_eigen_problem.isNonlinearEigenvalueSolver(number()))
281 : {
282 1189 : auto snes = getSNES();
283 :
284 : // nonlinear iterations
285 : PetscInt nl_its;
286 1189 : LibmeshPetscCallA(_eigen_problem.comm().get(), SNESGetIterationNumber(snes, &nl_its));
287 1189 : _n_iters = nl_its;
288 :
289 : // linear iterations
290 : PetscInt l_its;
291 1189 : LibmeshPetscCallA(_eigen_problem.comm().get(), SNESGetLinearSolveIterations(snes, &l_its));
292 1189 : _n_linear_iters = l_its;
293 :
294 : // final residual
295 : PetscReal norm;
296 1189 : LibmeshPetscCall(SNESGetFunctionNorm(snes, &norm));
297 1189 : _final_residual = norm;
298 : }
299 :
300 : // store eigenvalues
301 1244 : unsigned int n_converged_eigenvalues = getNumConvergedEigenvalues();
302 :
303 1244 : _n_eigen_pairs_required = _eigen_problem.getNEigenPairsRequired();
304 :
305 1244 : if (_n_eigen_pairs_required < n_converged_eigenvalues)
306 11 : n_converged_eigenvalues = _n_eigen_pairs_required;
307 :
308 1244 : _eigen_values.resize(n_converged_eigenvalues);
309 2534 : for (unsigned int n = 0; n < n_converged_eigenvalues; n++)
310 1290 : _eigen_values[n] = getConvergedEigenvalue(n);
311 :
312 : // Update the solution vector to the active eigenvector
313 1244 : if (n_converged_eigenvalues)
314 1235 : getConvergedEigenpair(_eigen_problem.activeEigenvalueIndex());
315 :
316 1244 : if (_eigen_problem.isNonlinearEigenvalueSolver(number()) && _num_constrained_dofs)
317 44 : solution().restore_subvector(std::move(subvec), _eigen_sys.local_non_condensed_dofs_vector);
318 1244 : }
319 :
320 : unsigned int
321 24 : NonlinearEigenSystem::nNonlinearIterations() const
322 : {
323 24 : if (!_time_integrators.empty())
324 0 : mooseError("Not implemented for time integrators.");
325 24 : if (!_eigen_problem.isNonlinearEigenvalueSolver(number()))
326 0 : mooseError("Only implemented for nonlinear eigenvalue solvers.");
327 :
328 24 : return _n_iters;
329 : }
330 :
331 : unsigned int
332 6 : NonlinearEigenSystem::nLinearIterations() const
333 : {
334 6 : if (!_time_integrators.empty())
335 0 : mooseError("Not implemented for time integrators.");
336 6 : if (!_eigen_problem.isNonlinearEigenvalueSolver(number()))
337 0 : mooseError("Only implemented for nonlinear eigenvalue solvers.");
338 :
339 6 : return _n_linear_iters;
340 : }
341 :
342 : Real
343 6 : NonlinearEigenSystem::finalNonlinearResidual() const
344 : {
345 6 : if (!_eigen_problem.isNonlinearEigenvalueSolver(number()))
346 0 : mooseError("Only implemented for nonlinear eigenvalue solvers.");
347 :
348 6 : return _final_residual;
349 : }
350 :
351 : void
352 709 : NonlinearEigenSystem::attachSLEPcCallbacks()
353 : {
354 : // Tell libmesh not to close matrices before solve
355 709 : _eigen_sys.get_eigen_solver().set_close_matrix_before_solve(false);
356 :
357 709 : if (_num_constrained_dofs)
358 : {
359 : // Condensed Matrix A
360 22 : if (_eigen_sys.has_condensed_matrix_A())
361 : {
362 11 : Mat mat = static_cast<PetscMatrix<Number> &>(_eigen_sys.get_condensed_matrix_A()).mat();
363 :
364 11 : Moose::SlepcSupport::attachCallbacksToMat(_eigen_problem, mat, false);
365 : }
366 :
367 : // Condensed Matrix B
368 22 : if (_eigen_sys.has_condensed_matrix_B())
369 : {
370 11 : Mat mat = static_cast<PetscMatrix<Number> &>(_eigen_sys.get_condensed_matrix_B()).mat();
371 :
372 11 : Moose::SlepcSupport::attachCallbacksToMat(_eigen_problem, mat, true);
373 : }
374 :
375 : // Condensed Preconditioning matrix
376 22 : if (_eigen_sys.has_condensed_precond_matrix())
377 : {
378 11 : Mat mat = static_cast<PetscMatrix<Number> &>(_eigen_sys.get_condensed_precond_matrix()).mat();
379 :
380 11 : Moose::SlepcSupport::attachCallbacksToMat(_eigen_problem, mat, true);
381 : }
382 : }
383 : else
384 : {
385 : // Matrix A
386 687 : if (_eigen_sys.has_matrix_A())
387 : {
388 139 : Mat mat = static_cast<PetscMatrix<Number> &>(_eigen_sys.get_matrix_A()).mat();
389 :
390 139 : Moose::SlepcSupport::attachCallbacksToMat(_eigen_problem, mat, false);
391 : }
392 :
393 : // Matrix B
394 687 : if (_eigen_sys.has_matrix_B())
395 : {
396 106 : Mat mat = static_cast<PetscMatrix<Number> &>(_eigen_sys.get_matrix_B()).mat();
397 :
398 106 : Moose::SlepcSupport::attachCallbacksToMat(_eigen_problem, mat, true);
399 : }
400 :
401 : // Preconditioning matrix
402 687 : if (_eigen_sys.has_precond_matrix())
403 : {
404 512 : Mat mat = static_cast<PetscMatrix<Number> &>(_eigen_sys.get_precond_matrix()).mat();
405 :
406 512 : Moose::SlepcSupport::attachCallbacksToMat(_eigen_problem, mat, true);
407 : }
408 : }
409 :
410 : // Shell matrix A
411 709 : if (_eigen_sys.has_shell_matrix_A())
412 : {
413 559 : Mat mat = static_cast<PetscShellMatrix<Number> &>(_eigen_sys.get_shell_matrix_A()).mat();
414 :
415 : // Attach callbacks for nonlinear eigenvalue solver
416 559 : Moose::SlepcSupport::attachCallbacksToMat(_eigen_problem, mat, false);
417 :
418 : // Set MatMult operations for shell
419 559 : Moose::SlepcSupport::setOperationsForShellMat(_eigen_problem, mat, false);
420 : }
421 :
422 : // Shell matrix B
423 709 : if (_eigen_sys.has_shell_matrix_B())
424 : {
425 559 : Mat mat = static_cast<PetscShellMatrix<Number> &>(_eigen_sys.get_shell_matrix_B()).mat();
426 :
427 559 : Moose::SlepcSupport::attachCallbacksToMat(_eigen_problem, mat, true);
428 :
429 : // Set MatMult operations for shell
430 559 : Moose::SlepcSupport::setOperationsForShellMat(_eigen_problem, mat, true);
431 : }
432 :
433 : // Shell preconditioning matrix
434 709 : if (_eigen_sys.has_shell_precond_matrix())
435 : {
436 36 : Mat mat = static_cast<PetscShellMatrix<Number> &>(_eigen_sys.get_shell_precond_matrix()).mat();
437 :
438 36 : Moose::SlepcSupport::attachCallbacksToMat(_eigen_problem, mat, true);
439 : }
440 709 : }
441 :
442 : void
443 0 : NonlinearEigenSystem::stopSolve(const ExecFlagType &, const std::set<TagID> &)
444 : {
445 0 : mooseError("did not implement yet \n");
446 : }
447 :
448 : void
449 0 : NonlinearEigenSystem::setupFiniteDifferencedPreconditioner()
450 : {
451 0 : mooseError("did not implement yet \n");
452 : }
453 :
454 : bool
455 709 : NonlinearEigenSystem::converged()
456 : {
457 709 : return _eigen_sys.get_n_converged();
458 : }
459 :
460 : unsigned int
461 0 : NonlinearEigenSystem::getCurrentNonlinearIterationNumber()
462 : {
463 0 : mooseError("did not implement yet \n");
464 : return 0;
465 : }
466 :
467 : NumericVector<Number> &
468 27 : NonlinearEigenSystem::RHS()
469 : {
470 27 : return _work_rhs_vector_BX;
471 : }
472 :
473 : NumericVector<Number> &
474 396 : NonlinearEigenSystem::residualVectorAX()
475 : {
476 396 : return _work_rhs_vector_AX;
477 : }
478 :
479 : NumericVector<Number> &
480 1785 : NonlinearEigenSystem::residualVectorBX()
481 : {
482 1785 : return _work_rhs_vector_BX;
483 : }
484 :
485 : NonlinearSolver<Number> *
486 0 : NonlinearEigenSystem::nonlinearSolver()
487 : {
488 0 : mooseError("did not implement yet \n");
489 : return NULL;
490 : }
491 :
492 : SNES
493 7031 : NonlinearEigenSystem::getSNES()
494 : {
495 7031 : EPS eps = getEPS();
496 :
497 7031 : if (_eigen_problem.isNonlinearEigenvalueSolver(number()))
498 : {
499 7031 : SNES snes = nullptr;
500 7031 : LibmeshPetscCall(Moose::SlepcSupport::mooseSlepcEPSGetSNES(eps, &snes));
501 7031 : return snes;
502 : }
503 : else
504 0 : mooseError("There is no SNES in linear eigen solver");
505 : }
506 :
507 : EPS
508 12724 : NonlinearEigenSystem::getEPS()
509 : {
510 : SlepcEigenSolver<Number> * solver =
511 12724 : cast_ptr<SlepcEigenSolver<Number> *>(&(*_eigen_sys.eigen_solver));
512 :
513 12724 : if (!solver)
514 0 : mooseError("Unable to retrieve eigen solver");
515 :
516 12724 : return solver->eps();
517 : }
518 :
519 : void
520 559 : NonlinearEigenSystem::checkIntegrity()
521 : {
522 559 : if (_nodal_bcs.hasActiveObjects())
523 : {
524 491 : const auto & nodal_bcs = _nodal_bcs.getActiveObjects();
525 1800 : for (const auto & nodal_bc : nodal_bcs)
526 : {
527 : // If this is a dirichlet boundary condition
528 1315 : auto nbc = std::dynamic_pointer_cast<DirichletBC>(nodal_bc);
529 : // If this is a eigen Dirichlet boundary condition
530 1315 : auto eigen_nbc = std::dynamic_pointer_cast<EigenDirichletBC>(nodal_bc);
531 : // ArrayDirichletBC
532 1315 : auto anbc = std::dynamic_pointer_cast<ArrayDirichletBC>(nodal_bc);
533 : // EigenArrayDirichletBC
534 1315 : auto aeigen_nbc = std::dynamic_pointer_cast<EigenArrayDirichletBC>(nodal_bc);
535 : // If it is a Dirichlet boundary condition, then value has to be zero
536 2333 : if (nbc && nbc->variable().eigen() && nbc->getParam<Real>("value"))
537 3 : mooseError(
538 : "Can't set an inhomogeneous Dirichlet boundary condition for eigenvalue problems.");
539 : // If it is an array Dirichlet boundary condition, all values should be zero
540 1312 : else if (anbc)
541 : {
542 216 : auto & values = anbc->getParam<RealEigenVector>("values");
543 324 : for (MooseIndex(values) i = 0; i < values.size(); i++)
544 : {
545 216 : if (values(i))
546 0 : mooseError("Can't set an inhomogeneous array Dirichlet boundary condition for "
547 : "eigenvalue problems.");
548 : }
549 : }
550 1204 : else if (!nbc && !eigen_nbc && !anbc && !aeigen_nbc)
551 3 : mooseError(
552 : "Invalid NodalBC for eigenvalue problems, please use homogeneous (array) Dirichlet.");
553 1309 : }
554 : }
555 553 : }
556 :
557 : std::pair<Real, Real>
558 1598 : NonlinearEigenSystem::getConvergedEigenvalue(dof_id_type n) const
559 : {
560 1598 : unsigned int n_converged_eigenvalues = getNumConvergedEigenvalues();
561 1598 : if (n >= n_converged_eigenvalues)
562 0 : mooseError(n, " not in [0, ", n_converged_eigenvalues, ")");
563 3196 : return _eigen_sys.get_eigenvalue(n);
564 : }
565 :
566 : std::pair<Real, Real>
567 1235 : NonlinearEigenSystem::getConvergedEigenpair(dof_id_type n) const
568 : {
569 1235 : unsigned int n_converged_eigenvalues = getNumConvergedEigenvalues();
570 :
571 1235 : if (n >= n_converged_eigenvalues)
572 0 : mooseError(n, " not in [0, ", n_converged_eigenvalues, ")");
573 :
574 2470 : return _eigen_sys.get_eigenpair(n);
575 : }
576 :
577 : void
578 45 : NonlinearEigenSystem::attachPreconditioner(Preconditioner<Number> * preconditioner)
579 : {
580 45 : _preconditioner = preconditioner;
581 :
582 : // If we have a customized preconditioner,
583 : // We need to let PETSc know that
584 45 : if (_preconditioner)
585 : {
586 45 : LibmeshPetscCall(Moose::SlepcSupport::registerPCToPETSc());
587 : // Mark this, and then we can setup correct petsc options
588 45 : _eigen_problem.solverParams(number())._customized_pc_for_eigen = true;
589 45 : _eigen_problem.solverParams(number())._type = Moose::ST_JFNK;
590 : }
591 45 : }
592 :
593 : void
594 45 : NonlinearEigenSystem::turnOffJacobian()
595 : {
596 : // Let us do nothing at the current moment
597 45 : }
598 :
599 : void
600 0 : NonlinearEigenSystem::residualAndJacobianTogether()
601 : {
602 0 : mooseError(
603 : "NonlinearEigenSystem::residualAndJacobianTogether is not implemented. It might even be "
604 : "nonsensical. If it is sensical and you want this capability, please contact a MOOSE "
605 : "developer.");
606 : }
607 :
608 : void
609 9 : NonlinearEigenSystem::computeScalingJacobian()
610 : {
611 9 : _eigen_problem.computeJacobianTag(*_current_solution, *_scaling_matrix, precondMatrixTag());
612 9 : }
613 :
614 : void
615 9 : NonlinearEigenSystem::computeScalingResidual()
616 : {
617 9 : _eigen_problem.computeResidualTag(*_current_solution, RHS(), nonEigenVectorTag());
618 9 : }
619 :
620 : std::set<TagID>
621 30615 : NonlinearEigenSystem::defaultVectorTags() const
622 : {
623 30615 : auto tags = NonlinearSystemBase::defaultVectorTags();
624 30615 : tags.insert(eigenVectorTag());
625 30615 : tags.insert(nonEigenVectorTag());
626 30615 : return tags;
627 0 : }
628 :
629 : std::set<TagID>
630 3937 : NonlinearEigenSystem::defaultMatrixTags() const
631 : {
632 3937 : auto tags = NonlinearSystemBase::defaultMatrixTags();
633 3937 : tags.insert(eigenMatrixTag());
634 3937 : tags.insert(nonEigenMatrixTag());
635 3937 : return tags;
636 0 : }
637 :
638 : #else
639 :
640 : NonlinearEigenSystem::NonlinearEigenSystem(EigenProblem & eigen_problem,
641 : const std::string & /*name*/)
642 : : libMesh::ParallelObject(eigen_problem)
643 : {
644 : mooseError("Need to install SLEPc to solve eigenvalue problems, please reconfigure libMesh\n");
645 : }
646 :
647 : #endif /* LIBMESH_HAVE_SLEPC */
|
