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 "PetscSupport.h"
11 :
12 : // MOOSE includes
13 : #include "MooseApp.h"
14 : #include "FEProblem.h"
15 : #include "DisplacedProblem.h"
16 : #include "NonlinearSystem.h"
17 : #include "LinearSystem.h"
18 : #include "DisplacedProblem.h"
19 : #include "PenetrationLocator.h"
20 : #include "NearestNodeLocator.h"
21 : #include "MooseTypes.h"
22 : #include "MooseUtils.h"
23 : #include "CommandLine.h"
24 : #include "Console.h"
25 : #include "MultiMooseEnum.h"
26 : #include "Conversion.h"
27 : #include "Executioner.h"
28 : #include "MooseMesh.h"
29 : #include "ComputeLineSearchObjectWrapper.h"
30 : #include "Convergence.h"
31 : #include "ParallelParamObject.h"
32 :
33 : #include "libmesh/equation_systems.h"
34 : #include "libmesh/linear_implicit_system.h"
35 : #include "libmesh/nonlinear_implicit_system.h"
36 : #include "libmesh/petsc_linear_solver.h"
37 : #include "libmesh/petsc_matrix.h"
38 : #include "libmesh/petsc_nonlinear_solver.h"
39 : #include "libmesh/petsc_preconditioner.h"
40 : #include "libmesh/petsc_vector.h"
41 : #include "libmesh/sparse_matrix.h"
42 : #include "libmesh/petsc_solver_exception.h"
43 : #include "libmesh/simple_range.h"
44 :
45 : // PETSc includes
46 : #include <petsc.h>
47 : #include <petscsnes.h>
48 : #include <petscksp.h>
49 :
50 : // For graph coloring
51 : #include <petscmat.h>
52 : #include <petscis.h>
53 : #include <petscdm.h>
54 :
55 : // PetscDMMoose include
56 : #include "PetscDMMoose.h"
57 :
58 : // Standard includes
59 : #include <ostream>
60 : #include <fstream>
61 : #include <string>
62 :
63 : using namespace libMesh;
64 :
65 : void
66 0 : MooseVecView(NumericVector<Number> & vector)
67 : {
68 0 : PetscVector<Number> & petsc_vec = static_cast<PetscVector<Number> &>(vector);
69 0 : LibmeshPetscCallA(vector.comm().get(), VecView(petsc_vec.vec(), 0));
70 0 : }
71 :
72 : void
73 0 : MooseMatView(SparseMatrix<Number> & mat)
74 : {
75 0 : PetscMatrixBase<Number> & petsc_mat = static_cast<PetscMatrix<Number> &>(mat);
76 0 : LibmeshPetscCallA(mat.comm().get(), MatView(petsc_mat.mat(), 0));
77 0 : }
78 :
79 : void
80 0 : MooseVecView(const NumericVector<Number> & vector)
81 : {
82 0 : PetscVector<Number> & petsc_vec =
83 : static_cast<PetscVector<Number> &>(const_cast<NumericVector<Number> &>(vector));
84 0 : LibmeshPetscCallA(vector.comm().get(), VecView(petsc_vec.vec(), 0));
85 0 : }
86 :
87 : void
88 0 : MooseMatView(const SparseMatrix<Number> & mat)
89 : {
90 0 : PetscMatrixBase<Number> & petsc_mat =
91 : static_cast<PetscMatrix<Number> &>(const_cast<SparseMatrix<Number> &>(mat));
92 0 : LibmeshPetscCallA(mat.comm().get(), MatView(petsc_mat.mat(), 0));
93 0 : }
94 :
95 : namespace Moose
96 : {
97 : namespace PetscSupport
98 : {
99 :
100 : std::string
101 3061 : stringify(const LineSearchType & t)
102 : {
103 3061 : switch (t)
104 : {
105 2966 : case LS_BASIC:
106 5932 : return "basic";
107 0 : case LS_DEFAULT:
108 0 : return "default";
109 0 : case LS_NONE:
110 0 : return "none";
111 0 : case LS_SHELL:
112 0 : return "shell";
113 0 : case LS_L2:
114 0 : return "l2";
115 23 : case LS_BT:
116 46 : return "bt";
117 72 : case LS_CP:
118 144 : return "cp";
119 0 : case LS_CONTACT:
120 0 : return "contact";
121 0 : case LS_PROJECT:
122 0 : return "project";
123 0 : case LS_INVALID:
124 0 : mooseError("Invalid LineSearchType");
125 : }
126 0 : return "";
127 : }
128 :
129 : std::string
130 29878 : stringify(const MffdType & t)
131 : {
132 29878 : switch (t)
133 : {
134 29868 : case MFFD_WP:
135 59736 : return "wp";
136 10 : case MFFD_DS:
137 20 : return "ds";
138 0 : case MFFD_INVALID:
139 0 : mooseError("Invalid MffdType");
140 : }
141 0 : return "";
142 : }
143 :
144 : void
145 44134 : setSolverOptions(const SolverParams & solver_params, const MultiMooseEnum & dont_add_these_options)
146 : {
147 44134 : const auto prefix_with_dash = '-' + solver_params._prefix;
148 : // set PETSc options implied by a solve type
149 44134 : switch (solver_params._type)
150 : {
151 29743 : case Moose::ST_PJFNK:
152 59486 : setSinglePetscOptionIfAppropriate(dont_add_these_options,
153 59486 : prefix_with_dash + "snes_mf_operator");
154 29743 : setSinglePetscOptionIfAppropriate(dont_add_these_options,
155 59486 : prefix_with_dash + "mat_mffd_type",
156 59486 : stringify(solver_params._mffd_type));
157 29743 : break;
158 :
159 135 : case Moose::ST_JFNK:
160 135 : setSinglePetscOptionIfAppropriate(dont_add_these_options, prefix_with_dash + "snes_mf");
161 135 : setSinglePetscOptionIfAppropriate(dont_add_these_options,
162 270 : prefix_with_dash + "mat_mffd_type",
163 270 : stringify(solver_params._mffd_type));
164 135 : break;
165 :
166 11909 : case Moose::ST_NEWTON:
167 11909 : break;
168 :
169 0 : case Moose::ST_FD:
170 0 : setSinglePetscOptionIfAppropriate(dont_add_these_options, prefix_with_dash + "snes_fd");
171 0 : break;
172 :
173 2347 : case Moose::ST_LINEAR:
174 7041 : setSinglePetscOptionIfAppropriate(
175 4694 : dont_add_these_options, prefix_with_dash + "snes_type", "ksponly");
176 4694 : setSinglePetscOptionIfAppropriate(dont_add_these_options,
177 4694 : prefix_with_dash + "snes_monitor_cancel");
178 2347 : break;
179 : }
180 :
181 44134 : Moose::LineSearchType ls_type = solver_params._line_search;
182 44134 : if (ls_type == Moose::LS_NONE)
183 2812 : ls_type = Moose::LS_BASIC;
184 :
185 44134 : if (ls_type != Moose::LS_DEFAULT && ls_type != Moose::LS_CONTACT && ls_type != Moose::LS_PROJECT)
186 3061 : setSinglePetscOptionIfAppropriate(
187 6122 : dont_add_these_options, prefix_with_dash + "snes_linesearch_type", stringify(ls_type));
188 44134 : }
189 :
190 : void
191 44458 : addPetscOptionsFromCommandline()
192 : {
193 : // commandline options always win
194 : // the options from a user commandline will overwrite the existing ones if any conflicts
195 : { // Get any options specified on the command-line
196 : int argc;
197 : char ** args;
198 :
199 44458 : LibmeshPetscCallA(PETSC_COMM_WORLD, PetscGetArgs(&argc, &args));
200 : #if PETSC_VERSION_LESS_THAN(3, 7, 0)
201 : LibmeshPetscCallA(PETSC_COMM_WORLD, PetscOptionsInsert(&argc, &args, NULL));
202 : #else
203 44458 : LibmeshPetscCallA(PETSC_COMM_WORLD,
204 : PetscOptionsInsert(LIBMESH_PETSC_NULLPTR, &argc, &args, NULL));
205 : #endif
206 : }
207 44458 : }
208 :
209 : void
210 43900 : petscSetOptionsHelper(const PetscOptions & po, FEProblemBase * const problem)
211 : {
212 : // Add any additional options specified in the input file
213 46510 : for (const auto & flag : po.flags)
214 : // Need to use name method here to pass a str instead of an EnumItem because
215 : // we don't care if the id attributes match
216 2610 : if (!po.dont_add_these_options.contains(flag.name()) ||
217 0 : po.user_set_options.contains(flag.name()))
218 10440 : setSinglePetscOption(flag.rawName().c_str());
219 :
220 : // Add option pairs
221 185662 : for (auto & option : po.pairs)
222 143285 : if (!po.dont_add_these_options.contains(option.first) ||
223 1523 : po.user_set_options.contains(option.first))
224 140249 : setSinglePetscOption(option.first, option.second, problem);
225 :
226 43900 : addPetscOptionsFromCommandline();
227 43900 : }
228 :
229 : void
230 1815 : petscSetOptions(const PetscOptions & po,
231 : const SolverParams & solver_params,
232 : FEProblemBase * const problem)
233 : {
234 1815 : PetscCallAbort(PETSC_COMM_WORLD, PetscOptionsClear(LIBMESH_PETSC_NULLPTR));
235 1815 : setSolverOptions(solver_params, po.dont_add_these_options);
236 1815 : petscSetOptionsHelper(po, problem);
237 1815 : }
238 :
239 : void
240 42085 : petscSetOptions(const PetscOptions & po,
241 : const std::vector<SolverParams> & solver_params_vec,
242 : FEProblemBase * const problem)
243 : {
244 42085 : PetscCallAbort(PETSC_COMM_WORLD, PetscOptionsClear(LIBMESH_PETSC_NULLPTR));
245 84404 : for (const auto & solver_params : solver_params_vec)
246 42319 : setSolverOptions(solver_params, po.dont_add_these_options);
247 42085 : petscSetOptionsHelper(po, problem);
248 42085 : }
249 :
250 : PetscErrorCode
251 49856 : petscSetupOutput(CommandLine * cmd_line)
252 : {
253 : PetscFunctionBegin;
254 49856 : char code[10] = {45, 45, 109, 111, 111, 115, 101};
255 49856 : const std::vector<std::string> argv = cmd_line->getArguments();
256 580919 : for (const auto & arg : argv)
257 : {
258 531063 : if (arg.compare(code) == 0)
259 : {
260 0 : Console::petscSetupOutput();
261 0 : break;
262 : }
263 : }
264 49856 : PetscFunctionReturn(PETSC_SUCCESS);
265 49856 : }
266 :
267 : PetscErrorCode
268 827541 : petscNonlinearConverged(SNES /*snes*/,
269 : PetscInt it,
270 : PetscReal /*xnorm*/,
271 : PetscReal /*snorm*/,
272 : PetscReal /*fnorm*/,
273 : SNESConvergedReason * reason,
274 : void * ctx)
275 : {
276 : PetscFunctionBegin;
277 827541 : FEProblemBase & problem = *static_cast<FEProblemBase *>(ctx);
278 :
279 : // execute objects that may be used in convergence check
280 827541 : problem.execute(EXEC_NONLINEAR_CONVERGENCE);
281 :
282 : // perform the convergence check
283 : Convergence::MooseConvergenceStatus status;
284 827541 : if (problem.getFailNextNonlinearConvergenceCheck())
285 : {
286 212 : status = Convergence::MooseConvergenceStatus::DIVERGED;
287 212 : problem.resetFailNextNonlinearConvergenceCheck();
288 : }
289 : else
290 : {
291 827329 : auto & convergence = problem.getConvergence(
292 827329 : problem.getNonlinearConvergenceNames()[problem.currentNonlinearSystem().number()]);
293 827329 : status = convergence.checkConvergence(it);
294 : }
295 :
296 : // convert convergence status to PETSc converged reason
297 827537 : switch (status)
298 : {
299 515958 : case Convergence::MooseConvergenceStatus::ITERATING:
300 515958 : *reason = SNES_CONVERGED_ITERATING;
301 515958 : break;
302 :
303 311315 : case Convergence::MooseConvergenceStatus::CONVERGED:
304 311315 : *reason = SNES_CONVERGED_FNORM_ABS;
305 311315 : break;
306 :
307 264 : case Convergence::MooseConvergenceStatus::DIVERGED:
308 264 : *reason = SNES_DIVERGED_DTOL;
309 264 : break;
310 : }
311 :
312 827537 : PetscFunctionReturn(PETSC_SUCCESS);
313 : }
314 :
315 : PetscErrorCode
316 2634 : petscLinearConverged(
317 : KSP /*ksp*/, PetscInt it, PetscReal /*norm*/, KSPConvergedReason * reason, void * ctx)
318 : {
319 : PetscFunctionBegin;
320 2634 : FEProblemBase & problem = *static_cast<FEProblemBase *>(ctx);
321 :
322 : // execute objects that may be used in convergence check
323 : // Right now, setting objects to execute on this flag would be ignored except in the
324 : // linear-system-only use case.
325 2634 : problem.execute(EXEC_LINEAR_CONVERGENCE);
326 :
327 : // perform the convergence check
328 : Convergence::MooseConvergenceStatus status;
329 2634 : if (problem.getFailNextSystemConvergenceCheck())
330 : {
331 0 : status = Convergence::MooseConvergenceStatus::DIVERGED;
332 0 : problem.resetFailNextSystemConvergenceCheck();
333 : }
334 : else
335 : {
336 2634 : auto & convergence = problem.getConvergence(
337 2634 : problem.getLinearConvergenceNames()[problem.currentLinearSystem().number()]);
338 2634 : status = convergence.checkConvergence(it);
339 : }
340 :
341 : // convert convergence status to PETSc converged reason
342 2634 : switch (status)
343 : {
344 2201 : case Convergence::MooseConvergenceStatus::ITERATING:
345 2201 : *reason = KSP_CONVERGED_ITERATING;
346 2201 : break;
347 :
348 : // TODO: find a KSP code that works better for this case
349 433 : case Convergence::MooseConvergenceStatus::CONVERGED:
350 : #if PETSC_VERSION_LESS_THAN(3, 24, 0)
351 433 : *reason = KSP_CONVERGED_RTOL_NORMAL;
352 : #else
353 : *reason = KSP_CONVERGED_RTOL_NORMAL_EQUATIONS;
354 : #endif
355 433 : break;
356 :
357 0 : case Convergence::MooseConvergenceStatus::DIVERGED:
358 0 : *reason = KSP_DIVERGED_DTOL;
359 0 : break;
360 : }
361 :
362 2634 : PetscFunctionReturn(PETSC_SUCCESS);
363 : }
364 :
365 : PCSide
366 0 : getPetscPCSide(Moose::PCSideType pcs)
367 : {
368 0 : switch (pcs)
369 : {
370 0 : case Moose::PCS_LEFT:
371 0 : return PC_LEFT;
372 0 : case Moose::PCS_RIGHT:
373 0 : return PC_RIGHT;
374 0 : case Moose::PCS_SYMMETRIC:
375 0 : return PC_SYMMETRIC;
376 0 : default:
377 0 : mooseError("Unknown PC side requested.");
378 : break;
379 : }
380 : }
381 :
382 : KSPNormType
383 410221 : getPetscKSPNormType(Moose::MooseKSPNormType kspnorm)
384 : {
385 410221 : switch (kspnorm)
386 : {
387 0 : case Moose::KSPN_NONE:
388 0 : return KSP_NORM_NONE;
389 0 : case Moose::KSPN_PRECONDITIONED:
390 0 : return KSP_NORM_PRECONDITIONED;
391 410221 : case Moose::KSPN_UNPRECONDITIONED:
392 410221 : return KSP_NORM_UNPRECONDITIONED;
393 0 : case Moose::KSPN_NATURAL:
394 0 : return KSP_NORM_NATURAL;
395 0 : case Moose::KSPN_DEFAULT:
396 0 : return KSP_NORM_DEFAULT;
397 0 : default:
398 0 : mooseError("Unknown KSP norm type requested.");
399 : break;
400 : }
401 : }
402 :
403 : void
404 383563 : petscSetDefaultKSPNormType(FEProblemBase & problem, KSP ksp)
405 : {
406 793784 : for (const auto i : make_range(problem.numSolverSystems()))
407 : {
408 410221 : SolverSystem & sys = problem.getSolverSystem(i);
409 410221 : LibmeshPetscCallA(problem.comm().get(),
410 : KSPSetNormType(ksp, getPetscKSPNormType(sys.getMooseKSPNormType())));
411 : }
412 383563 : }
413 :
414 : void
415 383563 : petscSetDefaultPCSide(FEProblemBase & problem, KSP ksp)
416 : {
417 793784 : for (const auto i : make_range(problem.numSolverSystems()))
418 : {
419 410221 : SolverSystem & sys = problem.getSolverSystem(i);
420 :
421 : // PETSc 3.2.x+
422 410221 : if (sys.getPCSide() != Moose::PCS_DEFAULT)
423 0 : LibmeshPetscCallA(problem.comm().get(), KSPSetPCSide(ksp, getPetscPCSide(sys.getPCSide())));
424 : }
425 383563 : }
426 :
427 : void
428 382317 : petscSetKSPDefaults(FEProblemBase & problem, KSP ksp)
429 : {
430 382317 : auto & es = problem.es();
431 :
432 382317 : PetscReal rtol = es.parameters.get<Real>("linear solver tolerance");
433 382317 : PetscReal atol = es.parameters.get<Real>("linear solver absolute tolerance");
434 :
435 : // MOOSE defaults this to -1 for some dumb reason
436 382317 : if (atol < 0)
437 0 : atol = 1e-50;
438 :
439 382317 : PetscReal maxits = es.parameters.get<unsigned int>("linear solver maximum iterations");
440 :
441 : // 1e100 is because we don't use divtol currently
442 382317 : LibmeshPetscCallA(problem.comm().get(), KSPSetTolerances(ksp, rtol, atol, 1e100, maxits));
443 :
444 382317 : petscSetDefaultPCSide(problem, ksp);
445 :
446 382317 : petscSetDefaultKSPNormType(problem, ksp);
447 382317 : }
448 :
449 : void
450 374316 : petscSetDefaults(FEProblemBase & problem)
451 : {
452 : // We care about both nonlinear and linear systems when setting the SNES prefix because
453 : // SNESSetOptionsPrefix will also set its KSP prefix which could compete with linear system KSPs
454 756313 : for (const auto nl_index : make_range(problem.numNonlinearSystems()))
455 : {
456 381997 : NonlinearSystemBase & nl = problem.getNonlinearSystemBase(nl_index);
457 :
458 : //
459 : // prefix system matrices
460 : //
461 :
462 381997 : auto & lm_sys = nl.system();
463 : // This check is necessary because we still have at least the system matrix lying around even
464 : // when doing matrix-free, but critically even though the libmesh object(s) exist, the wrapped
465 : // PETSc Mat(s) do not
466 381997 : if (problem.solverParams(nl_index)._type != Moose::ST_JFNK)
467 764912 : for (auto & [mat_name, mat] : as_range(lm_sys.matrices_begin(), lm_sys.matrices_end()))
468 : {
469 383085 : libmesh_ignore(mat_name);
470 383085 : if (auto * const petsc_mat = dynamic_cast<PetscMatrixBase<Number> *>(mat.get()); petsc_mat)
471 : {
472 383085 : LibmeshPetscCallA(nl.comm().get(),
473 : MatSetOptionsPrefix(petsc_mat->mat(), (nl.name() + "_").c_str()));
474 : // We should call this here to ensure that options from the command line are properly
475 : // applied
476 383085 : LibmeshPetscCallA(nl.comm().get(), MatSetFromOptions(petsc_mat->mat()));
477 : }
478 : }
479 :
480 : //
481 : // prefix SNES/KSP
482 : //
483 :
484 : // dig out PETSc solver
485 381997 : auto * const petsc_solver = cast_ptr<PetscNonlinearSolver<Number> *>(nl.nonlinearSolver());
486 381997 : const char * snes_prefix = nullptr;
487 381997 : std::string snes_prefix_str;
488 381997 : if (nl.system().prefix_with_name())
489 : {
490 26658 : snes_prefix_str = nl.system().prefix();
491 26658 : snes_prefix = snes_prefix_str.c_str();
492 : }
493 381997 : SNES snes = petsc_solver->snes(snes_prefix);
494 : KSP ksp;
495 381997 : LibmeshPetscCallA(nl.comm().get(), SNESGetKSP(snes, &ksp));
496 :
497 381997 : LibmeshPetscCallA(nl.comm().get(), SNESSetMaxLinearSolveFailures(snes, 1000000));
498 :
499 381997 : LibmeshPetscCallA(nl.comm().get(), SNESSetCheckJacobianDomainError(snes, PETSC_TRUE));
500 :
501 : // In 3.0.0, the context pointer must actually be used, and the
502 : // final argument to KSPSetConvergenceTest() is a pointer to a
503 : // routine for destroying said private data context. In this case,
504 : // we use the default context provided by PETSc in addition to
505 : // a few other tests.
506 : {
507 381997 : LibmeshPetscCallA(
508 : nl.comm().get(),
509 : SNESSetConvergenceTest(snes, petscNonlinearConverged, &problem, LIBMESH_PETSC_NULLPTR));
510 : }
511 :
512 381997 : petscSetKSPDefaults(problem, ksp);
513 381997 : }
514 :
515 380413 : for (auto sys_index : make_range(problem.numLinearSystems()))
516 : {
517 : // dig out PETSc solver
518 6097 : LinearSystem & lin_sys = problem.getLinearSystem(sys_index);
519 6097 : PetscLinearSolver<Number> * petsc_solver = dynamic_cast<PetscLinearSolver<Number> *>(
520 6097 : lin_sys.linearImplicitSystem().get_linear_solver());
521 6097 : KSP ksp = petsc_solver->ksp();
522 :
523 6097 : if (problem.hasLinearConvergenceObjects())
524 447 : LibmeshPetscCallA(
525 : lin_sys.comm().get(),
526 : KSPSetConvergenceTest(ksp, petscLinearConverged, &problem, LIBMESH_PETSC_NULLPTR));
527 :
528 : // We dont set the KSP defaults here because they seem to clash with the linear solve parameters
529 : // set in FEProblemBase::solveLinearSystem
530 : }
531 374316 : }
532 :
533 : void
534 79102 : processSingletonMooseWrappedOptions(FEProblemBase & fe_problem, const InputParameters & params)
535 : {
536 79102 : setSolveTypeFromParams(fe_problem, params);
537 79102 : setLineSearchFromParams(fe_problem, params);
538 79102 : setMFFDTypeFromParams(fe_problem, params);
539 79102 : }
540 :
541 : #define checkPrefix(prefix) \
542 : mooseAssert(prefix[0] == '-', \
543 : "Leading prefix character must be a '-'. Current prefix is '" << prefix << "'"); \
544 : mooseAssert((prefix.size() == 1) || (prefix.back() == '_'), \
545 : "Terminating prefix character must be a '_'. Current prefix is '" << prefix << "'"); \
546 : mooseAssert(MooseUtils::isAllLowercase(prefix), "PETSc prefixes should be all lower-case")
547 :
548 : void
549 79102 : storePetscOptions(FEProblemBase & fe_problem,
550 : const std::string & prefix,
551 : const ParallelParamObject & param_object)
552 : {
553 79102 : const auto & params = param_object.parameters();
554 79102 : processSingletonMooseWrappedOptions(fe_problem, params);
555 :
556 : // The parameters contained in the Action
557 79102 : const auto & petsc_options = params.get<MultiMooseEnum>("petsc_options");
558 : const auto & petsc_pair_options =
559 237306 : params.get<MooseEnumItem, std::string>("petsc_options_iname", "petsc_options_value");
560 :
561 : // A reference to the PetscOptions object that contains the settings that will be used in the
562 : // solve
563 79102 : auto & po = fe_problem.getPetscOptions();
564 :
565 : // First process the single petsc options/flags
566 79102 : addPetscFlagsToPetscOptions(petsc_options, prefix, param_object, po);
567 :
568 : // Then process the option-value pairs
569 158204 : addPetscPairsToPetscOptions(
570 79102 : petsc_pair_options, fe_problem.mesh().dimension(), prefix, param_object, po);
571 79102 : }
572 :
573 : void
574 94119 : setSolveTypeFromParams(FEProblemBase & fe_problem, const InputParameters & params)
575 : {
576 : // Note: Options set in the Preconditioner block will override those set in the Executioner block
577 360529 : if (params.isParamValid("solve_type") && !params.isParamValid("_use_eigen_value"))
578 : {
579 : // Extract the solve type
580 38498 : const std::string & solve_type = params.get<MooseEnum>("solve_type");
581 77442 : for (const auto i : make_range(fe_problem.numNonlinearSystems()))
582 38944 : fe_problem.solverParams(i)._type = Moose::stringToEnum<Moose::SolveType>(solve_type);
583 38498 : }
584 94119 : }
585 :
586 : void
587 79102 : setLineSearchFromParams(FEProblemBase & fe_problem, const InputParameters & params)
588 : {
589 : // Note: Options set in the Preconditioner block will override those set in the Executioner block
590 158204 : if (params.isParamValid("line_search"))
591 : {
592 63959 : const auto & line_search = params.get<MooseEnum>("line_search");
593 127100 : for (const auto i : make_range(fe_problem.numNonlinearSystems()))
594 63141 : if (fe_problem.solverParams(i)._line_search == Moose::LS_INVALID || line_search != "default")
595 : {
596 : Moose::LineSearchType enum_line_search =
597 62755 : Moose::stringToEnum<Moose::LineSearchType>(line_search);
598 62755 : fe_problem.solverParams(i)._line_search = enum_line_search;
599 62755 : if (enum_line_search == LS_CONTACT || enum_line_search == LS_PROJECT)
600 : {
601 0 : NonlinearImplicitSystem * nl_system = dynamic_cast<NonlinearImplicitSystem *>(
602 0 : &fe_problem.getNonlinearSystemBase(i).system());
603 0 : if (!nl_system)
604 0 : mooseError("You've requested a line search but you must be solving an EigenProblem. "
605 : "These two things are not consistent.");
606 : PetscNonlinearSolver<Real> * petsc_nonlinear_solver =
607 0 : dynamic_cast<PetscNonlinearSolver<Real> *>(nl_system->nonlinear_solver.get());
608 0 : if (!petsc_nonlinear_solver)
609 0 : mooseError("Currently the MOOSE line searches all use Petsc, so you "
610 : "must use Petsc as your non-linear solver.");
611 : petsc_nonlinear_solver->linesearch_object =
612 0 : std::make_unique<ComputeLineSearchObjectWrapper>(fe_problem);
613 : }
614 : }
615 : }
616 79102 : }
617 :
618 : void
619 79102 : setMFFDTypeFromParams(FEProblemBase & fe_problem, const InputParameters & params)
620 : {
621 158204 : if (params.isParamValid("mffd_type"))
622 : {
623 78522 : const auto & mffd_type = params.get<MooseEnum>("mffd_type");
624 156494 : for (const auto i : make_range(fe_problem.numNonlinearSystems()))
625 77972 : fe_problem.solverParams(i)._mffd_type = Moose::stringToEnum<Moose::MffdType>(mffd_type);
626 : }
627 79102 : }
628 :
629 : template <typename T>
630 : void
631 60327 : checkUserProvidedPetscOption(const T & option, const ParallelParamObject & param_object)
632 : {
633 60327 : const auto & string_option = static_cast<const std::string &>(option);
634 60327 : if (string_option[0] != '-')
635 0 : param_object.mooseError("PETSc option '", string_option, "' does not begin with '-'");
636 60327 : }
637 :
638 : void
639 79102 : addPetscFlagsToPetscOptions(const MultiMooseEnum & petsc_flags,
640 : std::string prefix,
641 : const ParallelParamObject & param_object,
642 : PetscOptions & po)
643 : {
644 79102 : prefix.insert(prefix.begin(), '-');
645 : checkPrefix(prefix);
646 :
647 : // Update the PETSc single flags
648 81829 : for (const auto & option : petsc_flags)
649 : {
650 2727 : checkUserProvidedPetscOption(option, param_object);
651 :
652 2727 : const std::string & string_option = option.name();
653 :
654 : /**
655 : * "-log_summary" cannot be used in the input file. This option needs to be set when PETSc is
656 : * initialized
657 : * which happens before the parser is even created. We'll throw an error if somebody attempts
658 : * to add this option later.
659 : */
660 2727 : if (option == "-log_summary" || option == "-log_view")
661 0 : mooseError("The PETSc option \"-log_summary\" or \"-log_view\" can only be used on the "
662 : "command line. Please "
663 : "remove it from the input file");
664 :
665 : // Update the stored items, but do not create duplicates
666 2727 : const std::string prefixed_option = prefix + string_option.substr(1);
667 2727 : if (!po.flags.isValueSet(prefixed_option))
668 : {
669 2727 : po.flags.setAdditionalValue(prefixed_option);
670 2727 : po.user_set_options.setAdditionalValue(prefixed_option);
671 : }
672 2727 : }
673 79102 : }
674 :
675 : void
676 63959 : setConvergedReasonFlags(FEProblemBase & fe_problem, std::string prefix)
677 : {
678 63959 : prefix.insert(prefix.begin(), '-');
679 : checkPrefix(prefix);
680 63959 : libmesh_ignore(fe_problem); // avoid unused warnings for old PETSc
681 :
682 : #if !PETSC_VERSION_LESS_THAN(3, 14, 0)
683 : // the boolean in these pairs denote whether the user has specified any of the reason flags in the
684 : // input file
685 127918 : std::array<std::string, 2> reason_flags = {{"snes_converged_reason", "ksp_converged_reason"}};
686 :
687 63959 : auto & po = fe_problem.getPetscOptions();
688 :
689 191877 : for (const auto & reason_flag : reason_flags)
690 382897 : if (!po.flags.isValueSet(prefix + reason_flag) &&
691 127061 : (std::find_if(po.pairs.begin(),
692 : po.pairs.end(),
693 179326 : [&reason_flag, &prefix](auto & pair)
694 561366 : { return pair.first == (prefix + reason_flag); }) == po.pairs.end()))
695 127061 : po.pairs.emplace_back(prefix + reason_flag, "::failed");
696 : #endif
697 63959 : }
698 :
699 : void
700 79344 : addPetscPairsToPetscOptions(
701 : const std::vector<std::pair<MooseEnumItem, std::string>> & petsc_pair_options,
702 : const unsigned int mesh_dimension,
703 : std::string prefix,
704 : const ParallelParamObject & param_object,
705 : PetscOptions & po)
706 : {
707 79344 : prefix.insert(prefix.begin(), '-');
708 : checkPrefix(prefix);
709 :
710 : // Setup the name value pairs
711 79344 : bool boomeramg_found = false;
712 79344 : bool strong_threshold_found = false;
713 : #if !PETSC_VERSION_LESS_THAN(3, 7, 0)
714 79344 : bool superlu_dist_found = false;
715 79344 : bool fact_pattern_found = false;
716 79344 : bool tiny_pivot_found = false;
717 : #endif
718 79344 : std::string pc_description = "";
719 : #if !PETSC_VERSION_LESS_THAN(3, 12, 0)
720 : // If users use HMG, we would like to set
721 79344 : bool hmg_found = false;
722 79344 : bool matptap_found = false;
723 79344 : bool hmg_strong_threshold_found = false;
724 : #endif
725 79344 : std::vector<std::pair<std::string, std::string>> new_options;
726 :
727 136944 : for (const auto & [option_name, option_value] : petsc_pair_options)
728 : {
729 57600 : checkUserProvidedPetscOption(option_name, param_object);
730 :
731 57600 : new_options.clear();
732 : const std::string prefixed_option_name =
733 57600 : prefix + static_cast<const std::string &>(option_name).substr(1);
734 :
735 : // Do not add duplicate settings
736 57600 : if (auto it =
737 57600 : MooseUtils::findPair(po.pairs, po.pairs.begin(), prefixed_option_name, MooseUtils::Any);
738 57600 : it == po.pairs.end())
739 : {
740 : #if !PETSC_VERSION_LESS_THAN(3, 9, 0)
741 57480 : if (option_name == "-pc_factor_mat_solver_package")
742 156 : new_options.emplace_back(prefix + "pc_factor_mat_solver_type", option_value);
743 : #else
744 : if (option_name == "-pc_factor_mat_solver_type")
745 : new_options.push_back(prefix + "pc_factor_mat_solver_package", option_value);
746 : #endif
747 :
748 : // Look for a pc description
749 86531 : if (option_name == "-pc_type" || option_name == "-sub_pc_type" ||
750 29051 : option_name == "-pc_hypre_type")
751 50427 : pc_description += option_value + ' ';
752 :
753 : #if !PETSC_VERSION_LESS_THAN(3, 12, 0)
754 57480 : if (option_name == "-pc_type" && option_value == "hmg")
755 42 : hmg_found = true;
756 :
757 : // MPIAIJ for PETSc 3.12.0: -matptap_via
758 : // MAIJ for PETSc 3.12.0: -matmaijptap_via
759 : // MPIAIJ for PETSc 3.13 to 3.16: -matptap_via, -matproduct_ptap_via
760 : // MAIJ for PETSc 3.13 to 3.16: -matproduct_ptap_via
761 : // MPIAIJ for PETSc 3.17 and higher: -matptap_via, -mat_product_algorithm
762 : // MAIJ for PETSc 3.17 and higher: -mat_product_algorithm
763 : #if !PETSC_VERSION_LESS_THAN(3, 17, 0)
764 57606 : if (hmg_found && (option_name == "-matptap_via" || option_name == "-matmaijptap_via" ||
765 126 : option_name == "-matproduct_ptap_via"))
766 0 : new_options.emplace_back(prefix + "mat_product_algorithm", option_value);
767 : #elif !PETSC_VERSION_LESS_THAN(3, 13, 0)
768 : if (hmg_found && (option_name == "-matptap_via" || option_name == "-matmaijptap_via"))
769 : new_options.emplace_back(prefix + "matproduct_ptap_via", option_value);
770 : #else
771 : if (hmg_found && (option_name == "-matproduct_ptap_via"))
772 : {
773 : new_options.emplace_back(prefix + "matptap_via", option_value);
774 : new_options.emplace_back(prefix + "matmaijptap_via", option_value);
775 : }
776 : #endif
777 :
778 114960 : if (option_name == "-matptap_via" || option_name == "-matmaijptap_via" ||
779 114960 : option_name == "-matproduct_ptap_via" || option_name == "-mat_product_algorithm")
780 0 : matptap_found = true;
781 :
782 : // For 3D problems, we need to set this 0.7
783 57480 : if (option_name == "-hmg_inner_pc_hypre_boomeramg_strong_threshold")
784 0 : hmg_strong_threshold_found = true;
785 : #endif
786 : // This special case is common enough that we'd like to handle it for the user.
787 57480 : if (option_name == "-pc_hypre_type" && option_value == "boomeramg")
788 21998 : boomeramg_found = true;
789 57480 : if (option_name == "-pc_hypre_boomeramg_strong_threshold")
790 0 : strong_threshold_found = true;
791 : #if !PETSC_VERSION_LESS_THAN(3, 7, 0)
792 57480 : if ((option_name == "-pc_factor_mat_solver_package" ||
793 58442 : option_name == "-pc_factor_mat_solver_type") &&
794 962 : option_value == "superlu_dist")
795 311 : superlu_dist_found = true;
796 57480 : if (option_name == "-mat_superlu_dist_fact")
797 0 : fact_pattern_found = true;
798 57480 : if (option_name == "-mat_superlu_dist_replacetinypivot")
799 0 : tiny_pivot_found = true;
800 : #endif
801 :
802 57480 : if (!new_options.empty())
803 : {
804 156 : std::copy(new_options.begin(), new_options.end(), std::back_inserter(po.pairs));
805 312 : for (const auto & option : new_options)
806 156 : po.user_set_options.setAdditionalValue(option.first);
807 : }
808 : else
809 : {
810 57324 : po.pairs.push_back(std::make_pair(prefixed_option_name, option_value));
811 57324 : po.user_set_options.setAdditionalValue(prefixed_option_name);
812 : }
813 : }
814 : else
815 : {
816 : do
817 : {
818 120 : it->second = option_value;
819 120 : it = MooseUtils::findPair(po.pairs, std::next(it), prefixed_option_name, MooseUtils::Any);
820 120 : } while (it != po.pairs.end());
821 : }
822 57600 : }
823 :
824 : // When running a 3D mesh with boomeramg, it is almost always best to supply a strong threshold
825 : // value. We will provide that for the user here if they haven't supplied it themselves.
826 79344 : if (boomeramg_found && !strong_threshold_found && mesh_dimension == 3)
827 : {
828 1033 : po.pairs.emplace_back(prefix + "pc_hypre_boomeramg_strong_threshold", "0.7");
829 1033 : pc_description += "strong_threshold: 0.7 (auto)";
830 : }
831 :
832 : #if !PETSC_VERSION_LESS_THAN(3, 12, 0)
833 79344 : if (hmg_found && !hmg_strong_threshold_found && mesh_dimension == 3)
834 : {
835 28 : po.pairs.emplace_back(prefix + "hmg_inner_pc_hypre_boomeramg_strong_threshold", "0.7");
836 28 : pc_description += "strong_threshold: 0.7 (auto)";
837 : }
838 :
839 : // Default PETSc PtAP takes too much memory, and it is not quite useful
840 : // Let us switch to use new algorithm
841 79344 : if (hmg_found && !matptap_found)
842 : {
843 : #if !PETSC_VERSION_LESS_THAN(3, 17, 0)
844 42 : po.pairs.emplace_back(prefix + "mat_product_algorithm", "allatonce");
845 : #elif !PETSC_VERSION_LESS_THAN(3, 13, 0)
846 : po.pairs.emplace_back(prefix + "matproduct_ptap_via", "allatonce");
847 : #else
848 : po.pairs.emplace_back(prefix + "matptap_via", "allatonce");
849 : po.pairs.emplace_back(prefix + "matmaijptap_via", "allatonce");
850 : #endif
851 : }
852 : #endif
853 :
854 : #if !PETSC_VERSION_LESS_THAN(3, 7, 0)
855 : // In PETSc-3.7.{0--4}, there is a bug when using superlu_dist, and we have to use
856 : // SamePattern_SameRowPerm, otherwise we use whatever we have in PETSc
857 79344 : if (superlu_dist_found && !fact_pattern_found)
858 : {
859 311 : po.pairs.emplace_back(prefix + "mat_superlu_dist_fact",
860 : #if PETSC_VERSION_LESS_THAN(3, 7, 5)
861 : "SamePattern_SameRowPerm");
862 : pc_description += "mat_superlu_dist_fact: SamePattern_SameRowPerm ";
863 : #else
864 : "SamePattern");
865 311 : pc_description += "mat_superlu_dist_fact: SamePattern ";
866 : #endif
867 : }
868 :
869 : // restore this superlu option
870 79344 : if (superlu_dist_found && !tiny_pivot_found)
871 : {
872 311 : po.pairs.emplace_back(prefix + "mat_superlu_dist_replacetinypivot", "1");
873 311 : pc_description += " mat_superlu_dist_replacetinypivot: true ";
874 : }
875 : #endif
876 : // Set Preconditioner description
877 79344 : if (!pc_description.empty() && prefix.size() > 1)
878 810 : po.pc_description += "[" + prefix.substr(1, prefix.size() - 2) + "]: ";
879 79344 : po.pc_description += pc_description;
880 79344 : }
881 :
882 : std::set<std::string>
883 305175 : getPetscValidLineSearches()
884 : {
885 305175 : return {"default", "shell", "none", "basic", "l2", "bt", "cp"};
886 : }
887 :
888 : InputParameters
889 438474 : getPetscValidParams()
890 : {
891 438474 : InputParameters params = emptyInputParameters();
892 :
893 1753896 : MooseEnum solve_type("PJFNK JFNK NEWTON FD LINEAR");
894 1753896 : params.addParam<MooseEnum>("solve_type",
895 : solve_type,
896 : "PJFNK: Preconditioned Jacobian-Free Newton Krylov "
897 : "JFNK: Jacobian-Free Newton Krylov "
898 : "NEWTON: Full Newton Solve "
899 : "FD: Use finite differences to compute Jacobian "
900 : "LINEAR: Solving a linear problem");
901 :
902 1753896 : MooseEnum mffd_type("wp ds", "wp");
903 1753896 : params.addParam<MooseEnum>("mffd_type",
904 : mffd_type,
905 : "Specifies the finite differencing type for "
906 : "Jacobian-free solve types. Note that the "
907 : "default is wp (for Walker and Pernice).");
908 :
909 1315422 : params.addParam<MultiMooseEnum>(
910 876948 : "petsc_options", getCommonPetscFlags(), "Singleton PETSc options");
911 1315422 : params.addParam<MultiMooseEnum>(
912 876948 : "petsc_options_iname", getCommonPetscKeys(), "Names of PETSc name/value pairs");
913 1753896 : params.addParam<std::vector<std::string>>(
914 : "petsc_options_value",
915 : "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\"");
916 1315422 : params.addParamNamesToGroup("solve_type petsc_options petsc_options_iname petsc_options_value "
917 : "mffd_type",
918 : "PETSc");
919 :
920 876948 : return params;
921 438474 : }
922 :
923 : MultiMooseEnum
924 455543 : getCommonSNESFlags()
925 : {
926 : return MultiMooseEnum(
927 : "-ksp_monitor_snes_lg -snes_ksp_ew -ksp_snes_ew -snes_converged_reason "
928 : "-snes_ksp -snes_linesearch_monitor -snes_mf -snes_mf_operator -snes_monitor "
929 : "-snes_test_display -snes_view -snes_monitor_cancel",
930 : "",
931 911086 : true);
932 : }
933 :
934 : MultiMooseEnum
935 454879 : getCommonKSPFlags()
936 : {
937 : return MultiMooseEnum(
938 909758 : "-ksp_converged_reason -ksp_gmres_modifiedgramschmidt -ksp_monitor", "", true);
939 : }
940 :
941 : MultiMooseEnum
942 454879 : getCommonPetscFlags()
943 : {
944 454879 : auto options = MultiMooseEnum("-dm_moose_print_embedding -dm_view", "", true);
945 454879 : options.addValidName(getCommonKSPFlags());
946 454879 : options.addValidName(getCommonSNESFlags());
947 454879 : return options;
948 0 : }
949 :
950 : MultiMooseEnum
951 455543 : getCommonSNESKeys()
952 : {
953 : return MultiMooseEnum("-snes_atol -snes_linesearch_type -snes_ls -snes_max_it -snes_rtol "
954 : "-snes_divergence_tolerance -snes_type",
955 : "",
956 911086 : true);
957 : }
958 :
959 : MultiMooseEnum
960 454879 : getCommonKSPKeys()
961 : {
962 : return MultiMooseEnum("-ksp_atol -ksp_gmres_restart -ksp_max_it -ksp_pc_side -ksp_rtol "
963 : "-ksp_type -sub_ksp_type",
964 : "",
965 909758 : true);
966 : }
967 : MultiMooseEnum
968 454879 : getCommonPetscKeys()
969 : {
970 : auto options = MultiMooseEnum("-mat_fd_coloring_err -mat_fd_type -mat_mffd_type "
971 : "-pc_asm_overlap -pc_factor_levels "
972 : "-pc_factor_mat_ordering_type -pc_hypre_boomeramg_grid_sweeps_all "
973 : "-pc_hypre_boomeramg_max_iter "
974 : "-pc_hypre_boomeramg_strong_threshold -pc_hypre_type -pc_type "
975 : "-sub_pc_type",
976 : "",
977 454879 : true);
978 454879 : options.addValidName(getCommonKSPKeys());
979 454879 : options.addValidName(getCommonSNESKeys());
980 454879 : return options;
981 0 : }
982 :
983 : bool
984 604 : isSNESVI(FEProblemBase & fe_problem)
985 : {
986 604 : const PetscOptions & petsc = fe_problem.getPetscOptions();
987 :
988 : int argc;
989 : char ** args;
990 604 : LibmeshPetscCallA(fe_problem.comm().get(), PetscGetArgs(&argc, &args));
991 :
992 604 : std::vector<std::string> cml_arg;
993 7504 : for (int i = 0; i < argc; i++)
994 13800 : cml_arg.push_back(args[i]);
995 :
996 604 : if (MooseUtils::findPair(petsc.pairs, petsc.pairs.begin(), MooseUtils::Any, "vinewtonssls") ==
997 1208 : petsc.pairs.end() &&
998 604 : MooseUtils::findPair(petsc.pairs, petsc.pairs.begin(), MooseUtils::Any, "vinewtonrsls") ==
999 648 : petsc.pairs.end() &&
1000 1252 : std::find(cml_arg.begin(), cml_arg.end(), "vinewtonssls") == cml_arg.end() &&
1001 648 : std::find(cml_arg.begin(), cml_arg.end(), "vinewtonrsls") == cml_arg.end())
1002 44 : return false;
1003 :
1004 560 : return true;
1005 604 : }
1006 :
1007 : void
1008 239156 : setSinglePetscOption(const std::string & name,
1009 : const std::string & value /*=""*/,
1010 : FEProblemBase * const problem /*=nullptr*/)
1011 : {
1012 239156 : static const TIMPI::Communicator comm_world(PETSC_COMM_WORLD);
1013 239156 : const TIMPI::Communicator & comm = problem ? problem->comm() : comm_world;
1014 239156 : LibmeshPetscCallA(comm.get(),
1015 : PetscOptionsSetValue(LIBMESH_PETSC_NULLPTR,
1016 : name.c_str(),
1017 : value == "" ? LIBMESH_PETSC_NULLPTR : value.c_str()));
1018 239156 : const auto lower_case_name = MooseUtils::toLower(name);
1019 28 : auto check_problem = [problem, &lower_case_name]()
1020 : {
1021 28 : if (!problem)
1022 0 : mooseError(
1023 : "Setting the option '",
1024 : lower_case_name,
1025 : "' requires passing a 'problem' parameter. Contact a developer of your application "
1026 : "to have them update their code. If in doubt, reach out to the MOOSE team on Github "
1027 : "discussions");
1028 239184 : };
1029 :
1030 : // Select vector type from user-passed PETSc options
1031 239156 : if (lower_case_name.find("-vec_type") != std::string::npos)
1032 : {
1033 14 : check_problem();
1034 28 : for (auto & solver_system : problem->_solver_systems)
1035 : {
1036 14 : auto & lm_sys = solver_system->system();
1037 50 : for (auto & [vec_name, vec] : as_range(lm_sys.vectors_begin(), lm_sys.vectors_end()))
1038 : {
1039 36 : libmesh_ignore(vec_name);
1040 36 : auto * const petsc_vec = cast_ptr<PetscVector<Number> *>(vec.get());
1041 36 : LibmeshPetscCallA(comm.get(), VecSetFromOptions(petsc_vec->vec()));
1042 : }
1043 : // The solution vectors aren't included in the system vectors storage
1044 14 : auto * petsc_vec = cast_ptr<PetscVector<Number> *>(lm_sys.solution.get());
1045 14 : LibmeshPetscCallA(comm.get(), VecSetFromOptions(petsc_vec->vec()));
1046 14 : petsc_vec = cast_ptr<PetscVector<Number> *>(lm_sys.current_local_solution.get());
1047 14 : LibmeshPetscCallA(comm.get(), VecSetFromOptions(petsc_vec->vec()));
1048 : }
1049 : }
1050 : // Select matrix type from user-passed PETSc options
1051 239142 : else if (lower_case_name.find("mat_type") != std::string::npos)
1052 : {
1053 14 : check_problem();
1054 :
1055 14 : bool found_matching_prefix = false;
1056 14 : for (const auto i : index_range(problem->_solver_systems))
1057 : {
1058 14 : const auto & solver_sys_name = problem->_solver_sys_names[i];
1059 14 : if (lower_case_name.find("-" + MooseUtils::toLower(solver_sys_name) + "_mat_type") ==
1060 : std::string::npos)
1061 0 : continue;
1062 :
1063 14 : if (problem->solverParams(i)._type == Moose::ST_JFNK)
1064 0 : mooseError(
1065 : "Setting option '", lower_case_name, "' is incompatible with a JFNK 'solve_type'");
1066 :
1067 14 : auto & lm_sys = problem->_solver_systems[i]->system();
1068 28 : for (auto & [mat_name, mat] : as_range(lm_sys.matrices_begin(), lm_sys.matrices_end()))
1069 : {
1070 14 : libmesh_ignore(mat_name);
1071 14 : if (auto * const petsc_mat = dynamic_cast<PetscMatrixBase<Number> *>(mat.get()); petsc_mat)
1072 : {
1073 : #ifdef DEBUG
1074 : const char * mat_prefix;
1075 : LibmeshPetscCallA(comm.get(), MatGetOptionsPrefix(petsc_mat->mat(), &mat_prefix));
1076 : mooseAssert(strcmp(mat_prefix, (solver_sys_name + "_").c_str()) == 0,
1077 : "We should have prefixed these matrices previously");
1078 : #endif
1079 14 : LibmeshPetscCallA(comm.get(), MatSetFromOptions(petsc_mat->mat()));
1080 : }
1081 : }
1082 14 : found_matching_prefix = true;
1083 14 : break;
1084 : }
1085 :
1086 14 : if (!found_matching_prefix)
1087 0 : mooseError("We did not find a matching solver system name for the petsc option '",
1088 : lower_case_name,
1089 : "'");
1090 : }
1091 239156 : }
1092 :
1093 : void
1094 73422 : setSinglePetscOptionIfAppropriate(const MultiMooseEnum & dont_add_these_options,
1095 : const std::string & name,
1096 : const std::string & value /*=""*/,
1097 : FEProblemBase * const problem /*=nullptr*/)
1098 : {
1099 73422 : if (!dont_add_these_options.contains(name))
1100 71985 : setSinglePetscOption(name, value, problem);
1101 73422 : }
1102 :
1103 : void
1104 0 : colorAdjacencyMatrix(PetscScalar * adjacency_matrix,
1105 : unsigned int size,
1106 : unsigned int colors,
1107 : std::vector<unsigned int> & vertex_colors,
1108 : const char * coloring_algorithm)
1109 : {
1110 : // Mat A will be a dense matrix from the incoming data structure
1111 : Mat A;
1112 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatCreate(PETSC_COMM_SELF, &A));
1113 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatSetSizes(A, size, size, size, size));
1114 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatSetType(A, MATSEQDENSE));
1115 : // PETSc requires a non-const data array to populate the matrix
1116 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatSeqDenseSetPreallocation(A, adjacency_matrix));
1117 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
1118 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
1119 :
1120 : // Convert A to a sparse matrix
1121 : #if PETSC_VERSION_LESS_THAN(3, 7, 0)
1122 : LibmeshPetscCallA(PETSC_COMM_SELF, MatConvert(A, MATAIJ, MAT_REUSE_MATRIX, &A));
1123 : #else
1124 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatConvert(A, MATAIJ, MAT_INPLACE_MATRIX, &A));
1125 : #endif
1126 :
1127 : ISColoring iscoloring;
1128 : MatColoring mc;
1129 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringCreate(A, &mc));
1130 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringSetType(mc, coloring_algorithm));
1131 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringSetMaxColors(mc, static_cast<PetscInt>(colors)));
1132 :
1133 : // Petsc normally colors by distance two (neighbors of neighbors), we just want one
1134 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringSetDistance(mc, 1));
1135 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringSetFromOptions(mc));
1136 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringApply(mc, &iscoloring));
1137 :
1138 : PetscInt nn;
1139 : IS * is;
1140 : #if PETSC_RELEASE_LESS_THAN(3, 12, 0)
1141 : LibmeshPetscCallA(PETSC_COMM_SELF, ISColoringGetIS(iscoloring, &nn, &is));
1142 : #else
1143 0 : LibmeshPetscCallA(PETSC_COMM_SELF, ISColoringGetIS(iscoloring, PETSC_USE_POINTER, &nn, &is));
1144 : #endif
1145 :
1146 0 : if (nn > static_cast<PetscInt>(colors))
1147 0 : throw std::runtime_error("Not able to color with designated number of colors");
1148 :
1149 0 : for (int i = 0; i < nn; i++)
1150 : {
1151 : PetscInt isize;
1152 : const PetscInt * indices;
1153 0 : LibmeshPetscCallA(PETSC_COMM_SELF, ISGetLocalSize(is[i], &isize));
1154 0 : LibmeshPetscCallA(PETSC_COMM_SELF, ISGetIndices(is[i], &indices));
1155 0 : for (int j = 0; j < isize; j++)
1156 : {
1157 : mooseAssert(indices[j] < static_cast<PetscInt>(vertex_colors.size()), "Index out of bounds");
1158 0 : vertex_colors[indices[j]] = i;
1159 : }
1160 0 : LibmeshPetscCallA(PETSC_COMM_SELF, ISRestoreIndices(is[i], &indices));
1161 : }
1162 :
1163 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatDestroy(&A));
1164 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringDestroy(&mc));
1165 0 : LibmeshPetscCallA(PETSC_COMM_SELF, ISColoringDestroy(&iscoloring));
1166 0 : }
1167 :
1168 : void
1169 3750 : dontAddPetscFlag(const std::string & flag, PetscOptions & petsc_options)
1170 : {
1171 3750 : if (!petsc_options.dont_add_these_options.contains(flag))
1172 3750 : petsc_options.dont_add_these_options.setAdditionalValue(flag);
1173 3750 : }
1174 :
1175 : void
1176 657 : dontAddNonlinearConvergedReason(FEProblemBase & fe_problem)
1177 : {
1178 1314 : dontAddPetscFlag("-snes_converged_reason", fe_problem.getPetscOptions());
1179 657 : }
1180 :
1181 : void
1182 667 : dontAddLinearConvergedReason(FEProblemBase & fe_problem)
1183 : {
1184 1334 : dontAddPetscFlag("-ksp_converged_reason", fe_problem.getPetscOptions());
1185 667 : }
1186 :
1187 : void
1188 0 : dontAddCommonKSPOptions(FEProblemBase & fe_problem)
1189 : {
1190 0 : auto & petsc_options = fe_problem.getPetscOptions();
1191 0 : for (const auto & flag : getCommonKSPFlags().getNames())
1192 0 : dontAddPetscFlag(flag, petsc_options);
1193 0 : for (const auto & key : getCommonKSPKeys().getNames())
1194 0 : dontAddPetscFlag(key, petsc_options);
1195 0 : }
1196 :
1197 : void
1198 664 : dontAddCommonSNESOptions(FEProblemBase & fe_problem)
1199 : {
1200 664 : auto & petsc_options = fe_problem.getPetscOptions();
1201 8632 : for (const auto & flag : getCommonSNESFlags().getNames())
1202 7968 : if (!petsc_options.dont_add_these_options.contains(flag))
1203 4648 : petsc_options.dont_add_these_options.setAdditionalValue(flag);
1204 5312 : for (const auto & key : getCommonSNESKeys().getNames())
1205 4648 : if (!petsc_options.dont_add_these_options.contains(key))
1206 2988 : petsc_options.dont_add_these_options.setAdditionalValue(key);
1207 664 : }
1208 :
1209 : std::unique_ptr<PetscMatrix<Number>>
1210 48 : createMatrixFromFile(const libMesh::Parallel::Communicator & comm,
1211 : Mat & mat,
1212 : const std::string & binary_mat_file,
1213 : const unsigned int mat_number_to_load)
1214 : {
1215 48 : LibmeshPetscCallA(comm.get(), MatCreate(comm.get(), &mat));
1216 : PetscViewer matviewer;
1217 48 : LibmeshPetscCallA(
1218 : comm.get(),
1219 : PetscViewerBinaryOpen(comm.get(), binary_mat_file.c_str(), FILE_MODE_READ, &matviewer));
1220 96 : for (unsigned int i = 0; i < mat_number_to_load; ++i)
1221 48 : LibmeshPetscCallA(comm.get(), MatLoad(mat, matviewer));
1222 48 : LibmeshPetscCallA(comm.get(), PetscViewerDestroy(&matviewer));
1223 :
1224 96 : return std::make_unique<PetscMatrix<Number>>(mat, comm);
1225 : }
1226 :
1227 : } // Namespace PetscSupport
1228 : } // Namespace MOOSE
|