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 29350 : stringify(const MffdType & t)
131 : {
132 29350 : switch (t)
133 : {
134 29340 : case MFFD_WP:
135 58680 : 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 43286 : setSolverOptions(const SolverParams & solver_params, const MultiMooseEnum & dont_add_these_options)
146 : {
147 43286 : const auto prefix_with_dash = '-' + solver_params._prefix;
148 : // set PETSc options implied by a solve type
149 43286 : switch (solver_params._type)
150 : {
151 29228 : case Moose::ST_PJFNK:
152 58456 : setSinglePetscOptionIfAppropriate(dont_add_these_options,
153 58456 : prefix_with_dash + "snes_mf_operator");
154 29228 : setSinglePetscOptionIfAppropriate(dont_add_these_options,
155 58456 : prefix_with_dash + "mat_mffd_type",
156 58456 : stringify(solver_params._mffd_type));
157 29228 : break;
158 :
159 122 : case Moose::ST_JFNK:
160 122 : setSinglePetscOptionIfAppropriate(dont_add_these_options, prefix_with_dash + "snes_mf");
161 122 : setSinglePetscOptionIfAppropriate(dont_add_these_options,
162 244 : prefix_with_dash + "mat_mffd_type",
163 244 : stringify(solver_params._mffd_type));
164 122 : break;
165 :
166 11619 : case Moose::ST_NEWTON:
167 11619 : 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 2317 : case Moose::ST_LINEAR:
174 6951 : setSinglePetscOptionIfAppropriate(
175 4634 : dont_add_these_options, prefix_with_dash + "snes_type", "ksponly");
176 4634 : setSinglePetscOptionIfAppropriate(dont_add_these_options,
177 4634 : prefix_with_dash + "snes_monitor_cancel");
178 2317 : break;
179 : }
180 :
181 43286 : Moose::LineSearchType ls_type = solver_params._line_search;
182 43286 : if (ls_type == Moose::LS_NONE)
183 2812 : ls_type = Moose::LS_BASIC;
184 :
185 43286 : 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 43286 : }
189 :
190 : void
191 43610 : 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 43610 : 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 43610 : LibmeshPetscCallA(PETSC_COMM_WORLD,
204 : PetscOptionsInsert(LIBMESH_PETSC_NULLPTR, &argc, &args, NULL));
205 : #endif
206 : }
207 43610 : }
208 :
209 : void
210 43052 : petscSetOptionsHelper(const PetscOptions & po, FEProblemBase * const problem)
211 : {
212 : // Add any additional options specified in the input file
213 45662 : 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 182279 : for (auto & option : po.pairs)
222 140750 : if (!po.dont_add_these_options.contains(option.first) ||
223 1523 : po.user_set_options.contains(option.first))
224 137714 : setSinglePetscOption(option.first, option.second, problem);
225 :
226 43052 : addPetscOptionsFromCommandline();
227 43052 : }
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 41237 : petscSetOptions(const PetscOptions & po,
241 : const std::vector<SolverParams> & solver_params_vec,
242 : FEProblemBase * const problem)
243 : {
244 41237 : PetscCallAbort(PETSC_COMM_WORLD, PetscOptionsClear(LIBMESH_PETSC_NULLPTR));
245 82708 : for (const auto & solver_params : solver_params_vec)
246 41471 : setSolverOptions(solver_params, po.dont_add_these_options);
247 41237 : petscSetOptionsHelper(po, problem);
248 41237 : }
249 :
250 : PetscErrorCode
251 48795 : petscSetupOutput(CommandLine * cmd_line)
252 : {
253 : PetscFunctionBegin;
254 48795 : char code[10] = {45, 45, 109, 111, 111, 115, 101};
255 48795 : const std::vector<std::string> argv = cmd_line->getArguments();
256 567234 : for (const auto & arg : argv)
257 : {
258 518439 : if (arg.compare(code) == 0)
259 : {
260 0 : Console::petscSetupOutput();
261 0 : break;
262 : }
263 : }
264 48795 : PetscFunctionReturn(PETSC_SUCCESS);
265 48795 : }
266 :
267 : PetscErrorCode
268 816586 : 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 816586 : FEProblemBase & problem = *static_cast<FEProblemBase *>(ctx);
278 :
279 : // execute objects that may be used in convergence check
280 816586 : problem.execute(EXEC_NONLINEAR_CONVERGENCE);
281 :
282 : // perform the convergence check
283 : Convergence::MooseConvergenceStatus status;
284 816586 : if (problem.getFailNextNonlinearConvergenceCheck())
285 : {
286 212 : status = Convergence::MooseConvergenceStatus::DIVERGED;
287 212 : problem.resetFailNextNonlinearConvergenceCheck();
288 : }
289 : else
290 : {
291 816374 : auto & convergence = problem.getConvergence(
292 816374 : problem.getNonlinearConvergenceNames()[problem.currentNonlinearSystem().number()]);
293 816374 : status = convergence.checkConvergence(it);
294 : }
295 :
296 : // convert convergence status to PETSc converged reason
297 816582 : switch (status)
298 : {
299 508626 : case Convergence::MooseConvergenceStatus::ITERATING:
300 508626 : *reason = SNES_CONVERGED_ITERATING;
301 508626 : break;
302 :
303 307692 : case Convergence::MooseConvergenceStatus::CONVERGED:
304 307692 : *reason = SNES_CONVERGED_FNORM_ABS;
305 307692 : break;
306 :
307 264 : case Convergence::MooseConvergenceStatus::DIVERGED:
308 264 : *reason = SNES_DIVERGED_DTOL;
309 264 : break;
310 : }
311 :
312 816582 : 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 433 : *reason = KSP_CONVERGED_RTOL_NORMAL;
351 433 : break;
352 :
353 0 : case Convergence::MooseConvergenceStatus::DIVERGED:
354 0 : *reason = KSP_DIVERGED_DTOL;
355 0 : break;
356 : }
357 :
358 2634 : PetscFunctionReturn(PETSC_SUCCESS);
359 : }
360 :
361 : PCSide
362 0 : getPetscPCSide(Moose::PCSideType pcs)
363 : {
364 0 : switch (pcs)
365 : {
366 0 : case Moose::PCS_LEFT:
367 0 : return PC_LEFT;
368 0 : case Moose::PCS_RIGHT:
369 0 : return PC_RIGHT;
370 0 : case Moose::PCS_SYMMETRIC:
371 0 : return PC_SYMMETRIC;
372 0 : default:
373 0 : mooseError("Unknown PC side requested.");
374 : break;
375 : }
376 : }
377 :
378 : KSPNormType
379 406284 : getPetscKSPNormType(Moose::MooseKSPNormType kspnorm)
380 : {
381 406284 : switch (kspnorm)
382 : {
383 0 : case Moose::KSPN_NONE:
384 0 : return KSP_NORM_NONE;
385 0 : case Moose::KSPN_PRECONDITIONED:
386 0 : return KSP_NORM_PRECONDITIONED;
387 406284 : case Moose::KSPN_UNPRECONDITIONED:
388 406284 : return KSP_NORM_UNPRECONDITIONED;
389 0 : case Moose::KSPN_NATURAL:
390 0 : return KSP_NORM_NATURAL;
391 0 : case Moose::KSPN_DEFAULT:
392 0 : return KSP_NORM_DEFAULT;
393 0 : default:
394 0 : mooseError("Unknown KSP norm type requested.");
395 : break;
396 : }
397 : }
398 :
399 : void
400 379626 : petscSetDefaultKSPNormType(FEProblemBase & problem, KSP ksp)
401 : {
402 785910 : for (const auto i : make_range(problem.numSolverSystems()))
403 : {
404 406284 : SolverSystem & sys = problem.getSolverSystem(i);
405 406284 : LibmeshPetscCallA(problem.comm().get(),
406 : KSPSetNormType(ksp, getPetscKSPNormType(sys.getMooseKSPNormType())));
407 : }
408 379626 : }
409 :
410 : void
411 379626 : petscSetDefaultPCSide(FEProblemBase & problem, KSP ksp)
412 : {
413 785910 : for (const auto i : make_range(problem.numSolverSystems()))
414 : {
415 406284 : SolverSystem & sys = problem.getSolverSystem(i);
416 :
417 : // PETSc 3.2.x+
418 406284 : if (sys.getPCSide() != Moose::PCS_DEFAULT)
419 0 : LibmeshPetscCallA(problem.comm().get(), KSPSetPCSide(ksp, getPetscPCSide(sys.getPCSide())));
420 : }
421 379626 : }
422 :
423 : void
424 378380 : petscSetKSPDefaults(FEProblemBase & problem, KSP ksp)
425 : {
426 378380 : auto & es = problem.es();
427 :
428 378380 : PetscReal rtol = es.parameters.get<Real>("linear solver tolerance");
429 378380 : PetscReal atol = es.parameters.get<Real>("linear solver absolute tolerance");
430 :
431 : // MOOSE defaults this to -1 for some dumb reason
432 378380 : if (atol < 0)
433 0 : atol = 1e-50;
434 :
435 378380 : PetscReal maxits = es.parameters.get<unsigned int>("linear solver maximum iterations");
436 :
437 : // 1e100 is because we don't use divtol currently
438 378380 : LibmeshPetscCallA(problem.comm().get(), KSPSetTolerances(ksp, rtol, atol, 1e100, maxits));
439 :
440 378380 : petscSetDefaultPCSide(problem, ksp);
441 :
442 378380 : petscSetDefaultKSPNormType(problem, ksp);
443 378380 : }
444 :
445 : void
446 370379 : petscSetDefaults(FEProblemBase & problem)
447 : {
448 : // We care about both nonlinear and linear systems when setting the SNES prefix because
449 : // SNESSetOptionsPrefix will also set its KSP prefix which could compete with linear system KSPs
450 748439 : for (const auto nl_index : make_range(problem.numNonlinearSystems()))
451 : {
452 378060 : NonlinearSystemBase & nl = problem.getNonlinearSystemBase(nl_index);
453 :
454 : //
455 : // prefix system matrices
456 : //
457 :
458 378060 : auto & lm_sys = nl.system();
459 : // This check is necessary because we still have at least the system matrix lying around even
460 : // when doing matrix-free, but critically even though the libmesh object(s) exist, the wrapped
461 : // PETSc Mat(s) do not
462 378060 : if (problem.solverParams(nl_index)._type != Moose::ST_JFNK)
463 756850 : for (auto & [mat_name, mat] : as_range(lm_sys.matrices_begin(), lm_sys.matrices_end()))
464 : {
465 378923 : libmesh_ignore(mat_name);
466 378923 : if (auto * const petsc_mat = dynamic_cast<PetscMatrixBase<Number> *>(mat.get()); petsc_mat)
467 : {
468 378923 : LibmeshPetscCallA(nl.comm().get(),
469 : MatSetOptionsPrefix(petsc_mat->mat(), (nl.name() + "_").c_str()));
470 : // We should call this here to ensure that options from the command line are properly
471 : // applied
472 378923 : LibmeshPetscCallA(nl.comm().get(), MatSetFromOptions(petsc_mat->mat()));
473 : }
474 : }
475 :
476 : //
477 : // prefix SNES/KSP
478 : //
479 :
480 : // dig out PETSc solver
481 378060 : auto * const petsc_solver = cast_ptr<PetscNonlinearSolver<Number> *>(nl.nonlinearSolver());
482 378060 : const char * snes_prefix = nullptr;
483 378060 : std::string snes_prefix_str;
484 378060 : if (nl.system().prefix_with_name())
485 : {
486 26658 : snes_prefix_str = nl.system().prefix();
487 26658 : snes_prefix = snes_prefix_str.c_str();
488 : }
489 378060 : SNES snes = petsc_solver->snes(snes_prefix);
490 : KSP ksp;
491 378060 : LibmeshPetscCallA(nl.comm().get(), SNESGetKSP(snes, &ksp));
492 :
493 378060 : LibmeshPetscCallA(nl.comm().get(), SNESSetMaxLinearSolveFailures(snes, 1000000));
494 :
495 378060 : LibmeshPetscCallA(nl.comm().get(), SNESSetCheckJacobianDomainError(snes, PETSC_TRUE));
496 :
497 : // In 3.0.0, the context pointer must actually be used, and the
498 : // final argument to KSPSetConvergenceTest() is a pointer to a
499 : // routine for destroying said private data context. In this case,
500 : // we use the default context provided by PETSc in addition to
501 : // a few other tests.
502 : {
503 378060 : LibmeshPetscCallA(
504 : nl.comm().get(),
505 : SNESSetConvergenceTest(snes, petscNonlinearConverged, &problem, LIBMESH_PETSC_NULLPTR));
506 : }
507 :
508 378060 : petscSetKSPDefaults(problem, ksp);
509 378060 : }
510 :
511 376476 : for (auto sys_index : make_range(problem.numLinearSystems()))
512 : {
513 : // dig out PETSc solver
514 6097 : LinearSystem & lin_sys = problem.getLinearSystem(sys_index);
515 6097 : PetscLinearSolver<Number> * petsc_solver = dynamic_cast<PetscLinearSolver<Number> *>(
516 6097 : lin_sys.linearImplicitSystem().get_linear_solver());
517 6097 : KSP ksp = petsc_solver->ksp();
518 :
519 6097 : if (problem.hasLinearConvergenceObjects())
520 447 : LibmeshPetscCallA(
521 : lin_sys.comm().get(),
522 : KSPSetConvergenceTest(ksp, petscLinearConverged, &problem, LIBMESH_PETSC_NULLPTR));
523 :
524 : // We dont set the KSP defaults here because they seem to clash with the linear solve parameters
525 : // set in FEProblemBase::solveLinearSystem
526 : }
527 370379 : }
528 :
529 : void
530 77633 : processSingletonMooseWrappedOptions(FEProblemBase & fe_problem, const InputParameters & params)
531 : {
532 77633 : setSolveTypeFromParams(fe_problem, params);
533 77633 : setLineSearchFromParams(fe_problem, params);
534 77633 : setMFFDTypeFromParams(fe_problem, params);
535 77633 : }
536 :
537 : #define checkPrefix(prefix) \
538 : mooseAssert(prefix[0] == '-', \
539 : "Leading prefix character must be a '-'. Current prefix is '" << prefix << "'"); \
540 : mooseAssert((prefix.size() == 1) || (prefix.back() == '_'), \
541 : "Terminating prefix character must be a '_'. Current prefix is '" << prefix << "'"); \
542 : mooseAssert(MooseUtils::isAllLowercase(prefix), "PETSc prefixes should be all lower-case")
543 :
544 : void
545 77633 : storePetscOptions(FEProblemBase & fe_problem,
546 : const std::string & prefix,
547 : const ParallelParamObject & param_object)
548 : {
549 77633 : const auto & params = param_object.parameters();
550 77633 : processSingletonMooseWrappedOptions(fe_problem, params);
551 :
552 : // The parameters contained in the Action
553 77633 : const auto & petsc_options = params.get<MultiMooseEnum>("petsc_options");
554 : const auto & petsc_pair_options =
555 232899 : params.get<MooseEnumItem, std::string>("petsc_options_iname", "petsc_options_value");
556 :
557 : // A reference to the PetscOptions object that contains the settings that will be used in the
558 : // solve
559 77633 : auto & po = fe_problem.getPetscOptions();
560 :
561 : // First process the single petsc options/flags
562 77633 : addPetscFlagsToPetscOptions(petsc_options, prefix, param_object, po);
563 :
564 : // Then process the option-value pairs
565 155266 : addPetscPairsToPetscOptions(
566 77633 : petsc_pair_options, fe_problem.mesh().dimension(), prefix, param_object, po);
567 77633 : }
568 :
569 : void
570 92266 : setSolveTypeFromParams(FEProblemBase & fe_problem, const InputParameters & params)
571 : {
572 : // Note: Options set in the Preconditioner block will override those set in the Executioner block
573 353366 : if (params.isParamValid("solve_type") && !params.isParamValid("_use_eigen_value"))
574 : {
575 : // Extract the solve type
576 37696 : const std::string & solve_type = params.get<MooseEnum>("solve_type");
577 75838 : for (const auto i : make_range(fe_problem.numNonlinearSystems()))
578 38142 : fe_problem.solverParams(i)._type = Moose::stringToEnum<Moose::SolveType>(solve_type);
579 37696 : }
580 92266 : }
581 :
582 : void
583 77633 : setLineSearchFromParams(FEProblemBase & fe_problem, const InputParameters & params)
584 : {
585 : // Note: Options set in the Preconditioner block will override those set in the Executioner block
586 155266 : if (params.isParamValid("line_search"))
587 : {
588 62870 : const auto & line_search = params.get<MooseEnum>("line_search");
589 124922 : for (const auto i : make_range(fe_problem.numNonlinearSystems()))
590 62052 : if (fe_problem.solverParams(i)._line_search == Moose::LS_INVALID || line_search != "default")
591 : {
592 : Moose::LineSearchType enum_line_search =
593 61666 : Moose::stringToEnum<Moose::LineSearchType>(line_search);
594 61666 : fe_problem.solverParams(i)._line_search = enum_line_search;
595 61666 : if (enum_line_search == LS_CONTACT || enum_line_search == LS_PROJECT)
596 : {
597 0 : NonlinearImplicitSystem * nl_system = dynamic_cast<NonlinearImplicitSystem *>(
598 0 : &fe_problem.getNonlinearSystemBase(i).system());
599 0 : if (!nl_system)
600 0 : mooseError("You've requested a line search but you must be solving an EigenProblem. "
601 : "These two things are not consistent.");
602 : PetscNonlinearSolver<Real> * petsc_nonlinear_solver =
603 0 : dynamic_cast<PetscNonlinearSolver<Real> *>(nl_system->nonlinear_solver.get());
604 0 : if (!petsc_nonlinear_solver)
605 0 : mooseError("Currently the MOOSE line searches all use Petsc, so you "
606 : "must use Petsc as your non-linear solver.");
607 : petsc_nonlinear_solver->linesearch_object =
608 0 : std::make_unique<ComputeLineSearchObjectWrapper>(fe_problem);
609 : }
610 : }
611 : }
612 77633 : }
613 :
614 : void
615 77633 : setMFFDTypeFromParams(FEProblemBase & fe_problem, const InputParameters & params)
616 : {
617 155266 : if (params.isParamValid("mffd_type"))
618 : {
619 77053 : const auto & mffd_type = params.get<MooseEnum>("mffd_type");
620 153556 : for (const auto i : make_range(fe_problem.numNonlinearSystems()))
621 76503 : fe_problem.solverParams(i)._mffd_type = Moose::stringToEnum<Moose::MffdType>(mffd_type);
622 : }
623 77633 : }
624 :
625 : template <typename T>
626 : void
627 59460 : checkUserProvidedPetscOption(const T & option, const ParallelParamObject & param_object)
628 : {
629 59460 : const auto & string_option = static_cast<const std::string &>(option);
630 59460 : if (string_option[0] != '-')
631 0 : param_object.mooseError("PETSc option '", string_option, "' does not begin with '-'");
632 59460 : }
633 :
634 : void
635 77633 : addPetscFlagsToPetscOptions(const MultiMooseEnum & petsc_flags,
636 : std::string prefix,
637 : const ParallelParamObject & param_object,
638 : PetscOptions & po)
639 : {
640 77633 : prefix.insert(prefix.begin(), '-');
641 : checkPrefix(prefix);
642 :
643 : // Update the PETSc single flags
644 80360 : for (const auto & option : petsc_flags)
645 : {
646 2727 : checkUserProvidedPetscOption(option, param_object);
647 :
648 2727 : const std::string & string_option = option.name();
649 :
650 : /**
651 : * "-log_summary" cannot be used in the input file. This option needs to be set when PETSc is
652 : * initialized
653 : * which happens before the parser is even created. We'll throw an error if somebody attempts
654 : * to add this option later.
655 : */
656 2727 : if (option == "-log_summary" || option == "-log_view")
657 0 : mooseError("The PETSc option \"-log_summary\" or \"-log_view\" can only be used on the "
658 : "command line. Please "
659 : "remove it from the input file");
660 :
661 : // Update the stored items, but do not create duplicates
662 2727 : const std::string prefixed_option = prefix + string_option.substr(1);
663 2727 : if (!po.flags.isValueSet(prefixed_option))
664 : {
665 2727 : po.flags.setAdditionalValue(prefixed_option);
666 2727 : po.user_set_options.setAdditionalValue(prefixed_option);
667 : }
668 2727 : }
669 77633 : }
670 :
671 : void
672 62870 : setConvergedReasonFlags(FEProblemBase & fe_problem, std::string prefix)
673 : {
674 62870 : prefix.insert(prefix.begin(), '-');
675 : checkPrefix(prefix);
676 62870 : libmesh_ignore(fe_problem); // avoid unused warnings for old PETSc
677 :
678 : #if !PETSC_VERSION_LESS_THAN(3, 14, 0)
679 : // the boolean in these pairs denote whether the user has specified any of the reason flags in the
680 : // input file
681 125740 : std::array<std::string, 2> reason_flags = {{"snes_converged_reason", "ksp_converged_reason"}};
682 :
683 62870 : auto & po = fe_problem.getPetscOptions();
684 :
685 188610 : for (const auto & reason_flag : reason_flags)
686 376363 : if (!po.flags.isValueSet(prefix + reason_flag) &&
687 124883 : (std::find_if(po.pairs.begin(),
688 : po.pairs.end(),
689 176405 : [&reason_flag, &prefix](auto & pair)
690 551911 : { return pair.first == (prefix + reason_flag); }) == po.pairs.end()))
691 124883 : po.pairs.emplace_back(prefix + reason_flag, "::failed");
692 : #endif
693 62870 : }
694 :
695 : void
696 77875 : addPetscPairsToPetscOptions(
697 : const std::vector<std::pair<MooseEnumItem, std::string>> & petsc_pair_options,
698 : const unsigned int mesh_dimension,
699 : std::string prefix,
700 : const ParallelParamObject & param_object,
701 : PetscOptions & po)
702 : {
703 77875 : prefix.insert(prefix.begin(), '-');
704 : checkPrefix(prefix);
705 :
706 : // Setup the name value pairs
707 77875 : bool boomeramg_found = false;
708 77875 : bool strong_threshold_found = false;
709 : #if !PETSC_VERSION_LESS_THAN(3, 7, 0)
710 77875 : bool superlu_dist_found = false;
711 77875 : bool fact_pattern_found = false;
712 77875 : bool tiny_pivot_found = false;
713 : #endif
714 77875 : std::string pc_description = "";
715 : #if !PETSC_VERSION_LESS_THAN(3, 12, 0)
716 : // If users use HMG, we would like to set
717 77875 : bool hmg_found = false;
718 77875 : bool matptap_found = false;
719 77875 : bool hmg_strong_threshold_found = false;
720 : #endif
721 77875 : std::vector<std::pair<std::string, std::string>> new_options;
722 :
723 134608 : for (const auto & [option_name, option_value] : petsc_pair_options)
724 : {
725 56733 : checkUserProvidedPetscOption(option_name, param_object);
726 :
727 56733 : new_options.clear();
728 : const std::string prefixed_option_name =
729 56733 : prefix + static_cast<const std::string &>(option_name).substr(1);
730 :
731 : // Do not add duplicate settings
732 56733 : if (auto it =
733 56733 : MooseUtils::findPair(po.pairs, po.pairs.begin(), prefixed_option_name, MooseUtils::Any);
734 56733 : it == po.pairs.end())
735 : {
736 : #if !PETSC_VERSION_LESS_THAN(3, 9, 0)
737 56613 : if (option_name == "-pc_factor_mat_solver_package")
738 156 : new_options.emplace_back(prefix + "pc_factor_mat_solver_type", option_value);
739 : #else
740 : if (option_name == "-pc_factor_mat_solver_type")
741 : new_options.push_back(prefix + "pc_factor_mat_solver_package", option_value);
742 : #endif
743 :
744 : // Look for a pc description
745 85193 : if (option_name == "-pc_type" || option_name == "-sub_pc_type" ||
746 28580 : option_name == "-pc_hypre_type")
747 49830 : pc_description += option_value + ' ';
748 :
749 : #if !PETSC_VERSION_LESS_THAN(3, 12, 0)
750 56613 : if (option_name == "-pc_type" && option_value == "hmg")
751 42 : hmg_found = true;
752 :
753 : // MPIAIJ for PETSc 3.12.0: -matptap_via
754 : // MAIJ for PETSc 3.12.0: -matmaijptap_via
755 : // MPIAIJ for PETSc 3.13 to 3.16: -matptap_via, -matproduct_ptap_via
756 : // MAIJ for PETSc 3.13 to 3.16: -matproduct_ptap_via
757 : // MPIAIJ for PETSc 3.17 and higher: -matptap_via, -mat_product_algorithm
758 : // MAIJ for PETSc 3.17 and higher: -mat_product_algorithm
759 : #if !PETSC_VERSION_LESS_THAN(3, 17, 0)
760 56739 : if (hmg_found && (option_name == "-matptap_via" || option_name == "-matmaijptap_via" ||
761 126 : option_name == "-matproduct_ptap_via"))
762 0 : new_options.emplace_back(prefix + "mat_product_algorithm", option_value);
763 : #elif !PETSC_VERSION_LESS_THAN(3, 13, 0)
764 : if (hmg_found && (option_name == "-matptap_via" || option_name == "-matmaijptap_via"))
765 : new_options.emplace_back(prefix + "matproduct_ptap_via", option_value);
766 : #else
767 : if (hmg_found && (option_name == "-matproduct_ptap_via"))
768 : {
769 : new_options.emplace_back(prefix + "matptap_via", option_value);
770 : new_options.emplace_back(prefix + "matmaijptap_via", option_value);
771 : }
772 : #endif
773 :
774 113226 : if (option_name == "-matptap_via" || option_name == "-matmaijptap_via" ||
775 113226 : option_name == "-matproduct_ptap_via" || option_name == "-mat_product_algorithm")
776 0 : matptap_found = true;
777 :
778 : // For 3D problems, we need to set this 0.7
779 56613 : if (option_name == "-hmg_inner_pc_hypre_boomeramg_strong_threshold")
780 0 : hmg_strong_threshold_found = true;
781 : #endif
782 : // This special case is common enough that we'd like to handle it for the user.
783 56613 : if (option_name == "-pc_hypre_type" && option_value == "boomeramg")
784 21797 : boomeramg_found = true;
785 56613 : if (option_name == "-pc_hypre_boomeramg_strong_threshold")
786 0 : strong_threshold_found = true;
787 : #if !PETSC_VERSION_LESS_THAN(3, 7, 0)
788 56613 : if ((option_name == "-pc_factor_mat_solver_package" ||
789 57575 : option_name == "-pc_factor_mat_solver_type") &&
790 962 : option_value == "superlu_dist")
791 311 : superlu_dist_found = true;
792 56613 : if (option_name == "-mat_superlu_dist_fact")
793 0 : fact_pattern_found = true;
794 56613 : if (option_name == "-mat_superlu_dist_replacetinypivot")
795 0 : tiny_pivot_found = true;
796 : #endif
797 :
798 56613 : if (!new_options.empty())
799 : {
800 156 : std::copy(new_options.begin(), new_options.end(), std::back_inserter(po.pairs));
801 312 : for (const auto & option : new_options)
802 156 : po.user_set_options.setAdditionalValue(option.first);
803 : }
804 : else
805 : {
806 56457 : po.pairs.push_back(std::make_pair(prefixed_option_name, option_value));
807 56457 : po.user_set_options.setAdditionalValue(prefixed_option_name);
808 : }
809 : }
810 : else
811 : {
812 : do
813 : {
814 120 : it->second = option_value;
815 120 : it = MooseUtils::findPair(po.pairs, std::next(it), prefixed_option_name, MooseUtils::Any);
816 120 : } while (it != po.pairs.end());
817 : }
818 56733 : }
819 :
820 : // When running a 3D mesh with boomeramg, it is almost always best to supply a strong threshold
821 : // value. We will provide that for the user here if they haven't supplied it themselves.
822 77875 : if (boomeramg_found && !strong_threshold_found && mesh_dimension == 3)
823 : {
824 984 : po.pairs.emplace_back(prefix + "pc_hypre_boomeramg_strong_threshold", "0.7");
825 984 : pc_description += "strong_threshold: 0.7 (auto)";
826 : }
827 :
828 : #if !PETSC_VERSION_LESS_THAN(3, 12, 0)
829 77875 : if (hmg_found && !hmg_strong_threshold_found && mesh_dimension == 3)
830 : {
831 28 : po.pairs.emplace_back(prefix + "hmg_inner_pc_hypre_boomeramg_strong_threshold", "0.7");
832 28 : pc_description += "strong_threshold: 0.7 (auto)";
833 : }
834 :
835 : // Default PETSc PtAP takes too much memory, and it is not quite useful
836 : // Let us switch to use new algorithm
837 77875 : if (hmg_found && !matptap_found)
838 : {
839 : #if !PETSC_VERSION_LESS_THAN(3, 17, 0)
840 42 : po.pairs.emplace_back(prefix + "mat_product_algorithm", "allatonce");
841 : #elif !PETSC_VERSION_LESS_THAN(3, 13, 0)
842 : po.pairs.emplace_back(prefix + "matproduct_ptap_via", "allatonce");
843 : #else
844 : po.pairs.emplace_back(prefix + "matptap_via", "allatonce");
845 : po.pairs.emplace_back(prefix + "matmaijptap_via", "allatonce");
846 : #endif
847 : }
848 : #endif
849 :
850 : #if !PETSC_VERSION_LESS_THAN(3, 7, 0)
851 : // In PETSc-3.7.{0--4}, there is a bug when using superlu_dist, and we have to use
852 : // SamePattern_SameRowPerm, otherwise we use whatever we have in PETSc
853 77875 : if (superlu_dist_found && !fact_pattern_found)
854 : {
855 311 : po.pairs.emplace_back(prefix + "mat_superlu_dist_fact",
856 : #if PETSC_VERSION_LESS_THAN(3, 7, 5)
857 : "SamePattern_SameRowPerm");
858 : pc_description += "mat_superlu_dist_fact: SamePattern_SameRowPerm ";
859 : #else
860 : "SamePattern");
861 311 : pc_description += "mat_superlu_dist_fact: SamePattern ";
862 : #endif
863 : }
864 :
865 : // restore this superlu option
866 77875 : if (superlu_dist_found && !tiny_pivot_found)
867 : {
868 311 : po.pairs.emplace_back(prefix + "mat_superlu_dist_replacetinypivot", "1");
869 311 : pc_description += " mat_superlu_dist_replacetinypivot: true ";
870 : }
871 : #endif
872 : // Set Preconditioner description
873 77875 : if (!pc_description.empty() && prefix.size() > 1)
874 810 : po.pc_description += "[" + prefix.substr(1, prefix.size() - 2) + "]: ";
875 77875 : po.pc_description += pc_description;
876 77875 : }
877 :
878 : std::set<std::string>
879 297171 : getPetscValidLineSearches()
880 : {
881 297171 : return {"default", "shell", "none", "basic", "l2", "bt", "cp"};
882 : }
883 :
884 : InputParameters
885 426320 : getPetscValidParams()
886 : {
887 426320 : InputParameters params = emptyInputParameters();
888 :
889 1705280 : MooseEnum solve_type("PJFNK JFNK NEWTON FD LINEAR");
890 1705280 : params.addParam<MooseEnum>("solve_type",
891 : solve_type,
892 : "PJFNK: Preconditioned Jacobian-Free Newton Krylov "
893 : "JFNK: Jacobian-Free Newton Krylov "
894 : "NEWTON: Full Newton Solve "
895 : "FD: Use finite differences to compute Jacobian "
896 : "LINEAR: Solving a linear problem");
897 :
898 1705280 : MooseEnum mffd_type("wp ds", "wp");
899 1705280 : params.addParam<MooseEnum>("mffd_type",
900 : mffd_type,
901 : "Specifies the finite differencing type for "
902 : "Jacobian-free solve types. Note that the "
903 : "default is wp (for Walker and Pernice).");
904 :
905 1278960 : params.addParam<MultiMooseEnum>(
906 852640 : "petsc_options", getCommonPetscFlags(), "Singleton PETSc options");
907 1278960 : params.addParam<MultiMooseEnum>(
908 852640 : "petsc_options_iname", getCommonPetscKeys(), "Names of PETSc name/value pairs");
909 1705280 : params.addParam<std::vector<std::string>>(
910 : "petsc_options_value",
911 : "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\"");
912 1278960 : params.addParamNamesToGroup("solve_type petsc_options petsc_options_iname petsc_options_value "
913 : "mffd_type",
914 : "PETSc");
915 :
916 852640 : return params;
917 426320 : }
918 :
919 : MultiMooseEnum
920 442907 : getCommonSNESFlags()
921 : {
922 : return MultiMooseEnum(
923 : "-ksp_monitor_snes_lg -snes_ksp_ew -ksp_snes_ew -snes_converged_reason "
924 : "-snes_ksp -snes_linesearch_monitor -snes_mf -snes_mf_operator -snes_monitor "
925 : "-snes_test_display -snes_view -snes_monitor_cancel",
926 : "",
927 885814 : true);
928 : }
929 :
930 : MultiMooseEnum
931 442243 : getCommonKSPFlags()
932 : {
933 : return MultiMooseEnum(
934 884486 : "-ksp_converged_reason -ksp_gmres_modifiedgramschmidt -ksp_monitor", "", true);
935 : }
936 :
937 : MultiMooseEnum
938 442243 : getCommonPetscFlags()
939 : {
940 442243 : auto options = MultiMooseEnum("-dm_moose_print_embedding -dm_view", "", true);
941 442243 : options.addValidName(getCommonKSPFlags());
942 442243 : options.addValidName(getCommonSNESFlags());
943 442243 : return options;
944 0 : }
945 :
946 : MultiMooseEnum
947 442907 : getCommonSNESKeys()
948 : {
949 : return MultiMooseEnum("-snes_atol -snes_linesearch_type -snes_ls -snes_max_it -snes_rtol "
950 : "-snes_divergence_tolerance -snes_type",
951 : "",
952 885814 : true);
953 : }
954 :
955 : MultiMooseEnum
956 442243 : getCommonKSPKeys()
957 : {
958 : return MultiMooseEnum("-ksp_atol -ksp_gmres_restart -ksp_max_it -ksp_pc_side -ksp_rtol "
959 : "-ksp_type -sub_ksp_type",
960 : "",
961 884486 : true);
962 : }
963 : MultiMooseEnum
964 442243 : getCommonPetscKeys()
965 : {
966 : auto options = MultiMooseEnum("-mat_fd_coloring_err -mat_fd_type -mat_mffd_type "
967 : "-pc_asm_overlap -pc_factor_levels "
968 : "-pc_factor_mat_ordering_type -pc_hypre_boomeramg_grid_sweeps_all "
969 : "-pc_hypre_boomeramg_max_iter "
970 : "-pc_hypre_boomeramg_strong_threshold -pc_hypre_type -pc_type "
971 : "-sub_pc_type",
972 : "",
973 442243 : true);
974 442243 : options.addValidName(getCommonKSPKeys());
975 442243 : options.addValidName(getCommonSNESKeys());
976 442243 : return options;
977 0 : }
978 :
979 : bool
980 604 : isSNESVI(FEProblemBase & fe_problem)
981 : {
982 604 : const PetscOptions & petsc = fe_problem.getPetscOptions();
983 :
984 : int argc;
985 : char ** args;
986 604 : LibmeshPetscCallA(fe_problem.comm().get(), PetscGetArgs(&argc, &args));
987 :
988 604 : std::vector<std::string> cml_arg;
989 7504 : for (int i = 0; i < argc; i++)
990 13800 : cml_arg.push_back(args[i]);
991 :
992 604 : if (MooseUtils::findPair(petsc.pairs, petsc.pairs.begin(), MooseUtils::Any, "vinewtonssls") ==
993 1208 : petsc.pairs.end() &&
994 604 : MooseUtils::findPair(petsc.pairs, petsc.pairs.begin(), MooseUtils::Any, "vinewtonrsls") ==
995 648 : petsc.pairs.end() &&
996 1252 : std::find(cml_arg.begin(), cml_arg.end(), "vinewtonssls") == cml_arg.end() &&
997 648 : std::find(cml_arg.begin(), cml_arg.end(), "vinewtonrsls") == cml_arg.end())
998 44 : return false;
999 :
1000 560 : return true;
1001 604 : }
1002 :
1003 : void
1004 235360 : setSinglePetscOption(const std::string & name,
1005 : const std::string & value /*=""*/,
1006 : FEProblemBase * const problem /*=nullptr*/)
1007 : {
1008 235360 : static const TIMPI::Communicator comm_world(PETSC_COMM_WORLD);
1009 235360 : const TIMPI::Communicator & comm = problem ? problem->comm() : comm_world;
1010 235360 : LibmeshPetscCallA(comm.get(),
1011 : PetscOptionsSetValue(LIBMESH_PETSC_NULLPTR,
1012 : name.c_str(),
1013 : value == "" ? LIBMESH_PETSC_NULLPTR : value.c_str()));
1014 235360 : const auto lower_case_name = MooseUtils::toLower(name);
1015 20 : auto check_problem = [problem, &lower_case_name]()
1016 : {
1017 20 : if (!problem)
1018 0 : mooseError(
1019 : "Setting the option '",
1020 : lower_case_name,
1021 : "' requires passing a 'problem' parameter. Contact a developer of your application "
1022 : "to have them update their code. If in doubt, reach out to the MOOSE team on Github "
1023 : "discussions");
1024 235380 : };
1025 :
1026 : // Select vector type from user-passed PETSc options
1027 235360 : if (lower_case_name.find("-vec_type") != std::string::npos)
1028 : {
1029 10 : check_problem();
1030 20 : for (auto & solver_system : problem->_solver_systems)
1031 : {
1032 10 : auto & lm_sys = solver_system->system();
1033 30 : for (auto & [vec_name, vec] : as_range(lm_sys.vectors_begin(), lm_sys.vectors_end()))
1034 : {
1035 20 : libmesh_ignore(vec_name);
1036 20 : auto * const petsc_vec = cast_ptr<PetscVector<Number> *>(vec.get());
1037 20 : LibmeshPetscCallA(comm.get(), VecSetFromOptions(petsc_vec->vec()));
1038 : }
1039 : // The solution vectors aren't included in the system vectors storage
1040 10 : auto * petsc_vec = cast_ptr<PetscVector<Number> *>(lm_sys.solution.get());
1041 10 : LibmeshPetscCallA(comm.get(), VecSetFromOptions(petsc_vec->vec()));
1042 10 : petsc_vec = cast_ptr<PetscVector<Number> *>(lm_sys.current_local_solution.get());
1043 10 : LibmeshPetscCallA(comm.get(), VecSetFromOptions(petsc_vec->vec()));
1044 : }
1045 : }
1046 : // Select matrix type from user-passed PETSc options
1047 235350 : else if (lower_case_name.find("mat_type") != std::string::npos)
1048 : {
1049 10 : check_problem();
1050 :
1051 10 : bool found_matching_prefix = false;
1052 10 : for (const auto i : index_range(problem->_solver_systems))
1053 : {
1054 10 : const auto & solver_sys_name = problem->_solver_sys_names[i];
1055 10 : if (lower_case_name.find("-" + MooseUtils::toLower(solver_sys_name) + "_mat_type") ==
1056 : std::string::npos)
1057 0 : continue;
1058 :
1059 10 : if (problem->solverParams(i)._type == Moose::ST_JFNK)
1060 0 : mooseError(
1061 : "Setting option '", lower_case_name, "' is incompatible with a JFNK 'solve_type'");
1062 :
1063 10 : auto & lm_sys = problem->_solver_systems[i]->system();
1064 20 : for (auto & [mat_name, mat] : as_range(lm_sys.matrices_begin(), lm_sys.matrices_end()))
1065 : {
1066 10 : libmesh_ignore(mat_name);
1067 10 : if (auto * const petsc_mat = dynamic_cast<PetscMatrixBase<Number> *>(mat.get()); petsc_mat)
1068 : {
1069 : #ifdef DEBUG
1070 : const char * mat_prefix;
1071 : LibmeshPetscCallA(comm.get(), MatGetOptionsPrefix(petsc_mat->mat(), &mat_prefix));
1072 : mooseAssert(strcmp(mat_prefix, (solver_sys_name + "_").c_str()) == 0,
1073 : "We should have prefixed these matrices previously");
1074 : #endif
1075 10 : LibmeshPetscCallA(comm.get(), MatSetFromOptions(petsc_mat->mat()));
1076 : }
1077 : }
1078 10 : found_matching_prefix = true;
1079 10 : break;
1080 : }
1081 :
1082 10 : if (!found_matching_prefix)
1083 0 : mooseError("We did not find a matching solver system name for the petsc option '",
1084 : lower_case_name,
1085 : "'");
1086 : }
1087 235360 : }
1088 :
1089 : void
1090 72306 : setSinglePetscOptionIfAppropriate(const MultiMooseEnum & dont_add_these_options,
1091 : const std::string & name,
1092 : const std::string & value /*=""*/,
1093 : FEProblemBase * const problem /*=nullptr*/)
1094 : {
1095 72306 : if (!dont_add_these_options.contains(name))
1096 70869 : setSinglePetscOption(name, value, problem);
1097 72306 : }
1098 :
1099 : void
1100 0 : colorAdjacencyMatrix(PetscScalar * adjacency_matrix,
1101 : unsigned int size,
1102 : unsigned int colors,
1103 : std::vector<unsigned int> & vertex_colors,
1104 : const char * coloring_algorithm)
1105 : {
1106 : // Mat A will be a dense matrix from the incoming data structure
1107 : Mat A;
1108 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatCreate(PETSC_COMM_SELF, &A));
1109 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatSetSizes(A, size, size, size, size));
1110 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatSetType(A, MATSEQDENSE));
1111 : // PETSc requires a non-const data array to populate the matrix
1112 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatSeqDenseSetPreallocation(A, adjacency_matrix));
1113 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
1114 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
1115 :
1116 : // Convert A to a sparse matrix
1117 : #if PETSC_VERSION_LESS_THAN(3, 7, 0)
1118 : LibmeshPetscCallA(PETSC_COMM_SELF, MatConvert(A, MATAIJ, MAT_REUSE_MATRIX, &A));
1119 : #else
1120 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatConvert(A, MATAIJ, MAT_INPLACE_MATRIX, &A));
1121 : #endif
1122 :
1123 : ISColoring iscoloring;
1124 : MatColoring mc;
1125 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringCreate(A, &mc));
1126 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringSetType(mc, coloring_algorithm));
1127 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringSetMaxColors(mc, static_cast<PetscInt>(colors)));
1128 :
1129 : // Petsc normally colors by distance two (neighbors of neighbors), we just want one
1130 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringSetDistance(mc, 1));
1131 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringSetFromOptions(mc));
1132 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringApply(mc, &iscoloring));
1133 :
1134 : PetscInt nn;
1135 : IS * is;
1136 : #if PETSC_RELEASE_LESS_THAN(3, 12, 0)
1137 : LibmeshPetscCallA(PETSC_COMM_SELF, ISColoringGetIS(iscoloring, &nn, &is));
1138 : #else
1139 0 : LibmeshPetscCallA(PETSC_COMM_SELF, ISColoringGetIS(iscoloring, PETSC_USE_POINTER, &nn, &is));
1140 : #endif
1141 :
1142 0 : if (nn > static_cast<PetscInt>(colors))
1143 0 : throw std::runtime_error("Not able to color with designated number of colors");
1144 :
1145 0 : for (int i = 0; i < nn; i++)
1146 : {
1147 : PetscInt isize;
1148 : const PetscInt * indices;
1149 0 : LibmeshPetscCallA(PETSC_COMM_SELF, ISGetLocalSize(is[i], &isize));
1150 0 : LibmeshPetscCallA(PETSC_COMM_SELF, ISGetIndices(is[i], &indices));
1151 0 : for (int j = 0; j < isize; j++)
1152 : {
1153 : mooseAssert(indices[j] < static_cast<PetscInt>(vertex_colors.size()), "Index out of bounds");
1154 0 : vertex_colors[indices[j]] = i;
1155 : }
1156 0 : LibmeshPetscCallA(PETSC_COMM_SELF, ISRestoreIndices(is[i], &indices));
1157 : }
1158 :
1159 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatDestroy(&A));
1160 0 : LibmeshPetscCallA(PETSC_COMM_SELF, MatColoringDestroy(&mc));
1161 0 : LibmeshPetscCallA(PETSC_COMM_SELF, ISColoringDestroy(&iscoloring));
1162 0 : }
1163 :
1164 : void
1165 3748 : dontAddPetscFlag(const std::string & flag, PetscOptions & petsc_options)
1166 : {
1167 3748 : if (!petsc_options.dont_add_these_options.contains(flag))
1168 3748 : petsc_options.dont_add_these_options.setAdditionalValue(flag);
1169 3748 : }
1170 :
1171 : void
1172 656 : dontAddNonlinearConvergedReason(FEProblemBase & fe_problem)
1173 : {
1174 1312 : dontAddPetscFlag("-snes_converged_reason", fe_problem.getPetscOptions());
1175 656 : }
1176 :
1177 : void
1178 666 : dontAddLinearConvergedReason(FEProblemBase & fe_problem)
1179 : {
1180 1332 : dontAddPetscFlag("-ksp_converged_reason", fe_problem.getPetscOptions());
1181 666 : }
1182 :
1183 : void
1184 0 : dontAddCommonKSPOptions(FEProblemBase & fe_problem)
1185 : {
1186 0 : auto & petsc_options = fe_problem.getPetscOptions();
1187 0 : for (const auto & flag : getCommonKSPFlags().getNames())
1188 0 : dontAddPetscFlag(flag, petsc_options);
1189 0 : for (const auto & key : getCommonKSPKeys().getNames())
1190 0 : dontAddPetscFlag(key, petsc_options);
1191 0 : }
1192 :
1193 : void
1194 664 : dontAddCommonSNESOptions(FEProblemBase & fe_problem)
1195 : {
1196 664 : auto & petsc_options = fe_problem.getPetscOptions();
1197 8632 : for (const auto & flag : getCommonSNESFlags().getNames())
1198 7968 : if (!petsc_options.dont_add_these_options.contains(flag))
1199 4648 : petsc_options.dont_add_these_options.setAdditionalValue(flag);
1200 5312 : for (const auto & key : getCommonSNESKeys().getNames())
1201 4648 : if (!petsc_options.dont_add_these_options.contains(key))
1202 2988 : petsc_options.dont_add_these_options.setAdditionalValue(key);
1203 664 : }
1204 :
1205 : std::unique_ptr<PetscMatrix<Number>>
1206 48 : createMatrixFromFile(const libMesh::Parallel::Communicator & comm,
1207 : Mat & mat,
1208 : const std::string & binary_mat_file,
1209 : const unsigned int mat_number_to_load)
1210 : {
1211 48 : LibmeshPetscCallA(comm.get(), MatCreate(comm.get(), &mat));
1212 : PetscViewer matviewer;
1213 48 : LibmeshPetscCallA(
1214 : comm.get(),
1215 : PetscViewerBinaryOpen(comm.get(), binary_mat_file.c_str(), FILE_MODE_READ, &matviewer));
1216 96 : for (unsigned int i = 0; i < mat_number_to_load; ++i)
1217 48 : LibmeshPetscCallA(comm.get(), MatLoad(mat, matviewer));
1218 48 : LibmeshPetscCallA(comm.get(), PetscViewerDestroy(&matviewer));
1219 :
1220 96 : return std::make_unique<PetscMatrix<Number>>(mat, comm);
1221 : }
1222 :
1223 : } // Namespace PetscSupport
1224 : } // Namespace MOOSE
|
