www.mooseframework.org
PetscSupport.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
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 #ifdef LIBMESH_HAVE_PETSC
13 
14 // MOOSE includes
15 #include "MooseApp.h"
16 #include "FEProblem.h"
17 #include "DisplacedProblem.h"
18 #include "NonlinearSystem.h"
19 #include "DisplacedProblem.h"
20 #include "PenetrationLocator.h"
21 #include "NearestNodeLocator.h"
22 #include "MooseTypes.h"
23 #include "MooseUtils.h"
24 #include "CommandLine.h"
25 #include "Console.h"
26 #include "MultiMooseEnum.h"
27 #include "Conversion.h"
28 #include "Executioner.h"
29 #include "MooseMesh.h"
31 
32 #include "libmesh/equation_systems.h"
33 #include "libmesh/linear_implicit_system.h"
34 #include "libmesh/nonlinear_implicit_system.h"
35 #include "libmesh/petsc_linear_solver.h"
36 #include "libmesh/petsc_matrix.h"
37 #include "libmesh/petsc_nonlinear_solver.h"
38 #include "libmesh/petsc_preconditioner.h"
39 #include "libmesh/petsc_vector.h"
40 #include "libmesh/sparse_matrix.h"
41 
42 // PETSc includes
43 #include <petsc.h>
44 #include <petscsnes.h>
45 #include <petscksp.h>
46 
47 // For graph coloring
48 #include <petscmat.h>
49 #include <petscis.h>
50 
51 #if PETSC_VERSION_LESS_THAN(3, 3, 0)
52 // PETSc 3.2.x and lower
53 #include <private/kspimpl.h>
54 #include <private/snesimpl.h>
55 #else
56 // PETSc 3.3.0+
57 #include <petscdm.h>
58 #endif
59 
60 // PetscDMMoose include
61 #include "PetscDMMoose.h"
62 
63 // Standard includes
64 #include <ostream>
65 #include <fstream>
66 #include <string>
67 
68 namespace Moose
69 {
70 namespace PetscSupport
71 {
72 
73 std::string
75 {
76  switch (t)
77  {
78  case LS_BASIC:
79  return "basic";
80  case LS_DEFAULT:
81  return "default";
82  case LS_NONE:
83  return "none";
84 #if PETSC_VERSION_LESS_THAN(3, 3, 0)
85  case LS_CUBIC:
86  return "cubic";
87  case LS_QUADRATIC:
88  return "quadratic";
89  case LS_BASICNONORMS:
90  return "basicnonorms";
91 #else
92  case LS_SHELL:
93  return "shell";
94  case LS_L2:
95  return "l2";
96  case LS_BT:
97  return "bt";
98  case LS_CP:
99  return "cp";
100  case LS_CONTACT:
101  return "contact";
102 #endif
103  case LS_INVALID:
104  mooseError("Invalid LineSearchType");
105  }
106  return "";
107 }
108 
109 std::string
110 stringify(const MffdType & t)
111 {
112  switch (t)
113  {
114  case MFFD_WP:
115  return "wp";
116  case MFFD_DS:
117  return "ds";
118  case MFFD_INVALID:
119  mooseError("Invalid MffdType");
120  }
121  return "";
122 }
123 
124 void
126 {
127  // set PETSc options implied by a solve type
128  switch (solver_params._type)
129  {
130  case Moose::ST_PJFNK:
131  setSinglePetscOption("-snes_mf_operator");
132  setSinglePetscOption("-mat_mffd_type", stringify(solver_params._mffd_type));
133  break;
134 
135  case Moose::ST_JFNK:
136  setSinglePetscOption("-snes_mf");
137  setSinglePetscOption("-mat_mffd_type", stringify(solver_params._mffd_type));
138  break;
139 
140  case Moose::ST_NEWTON:
141  break;
142 
143  case Moose::ST_FD:
144  setSinglePetscOption("-snes_fd");
145  break;
146 
147  case Moose::ST_LINEAR:
148  setSinglePetscOption("-snes_type", "ksponly");
149  setSinglePetscOption("-snes_monitor_cancel");
150  break;
151  }
152 
153  Moose::LineSearchType ls_type = solver_params._line_search;
154  if (ls_type == Moose::LS_NONE)
155  ls_type = Moose::LS_BASIC;
156 
157  if (ls_type != Moose::LS_DEFAULT && ls_type != Moose::LS_CONTACT)
158  {
159 #if PETSC_VERSION_LESS_THAN(3, 3, 0)
160  setSinglePetscOption("-snes_type", "ls");
161  setSinglePetscOption("-snes_ls", stringify(ls_type));
162 #else
163  setSinglePetscOption("-snes_linesearch_type", stringify(ls_type));
164 #endif
165  }
166 }
167 
168 void
170 {
171 #if !PETSC_VERSION_LESS_THAN(3, 3, 0)
172  PetscErrorCode ierr;
173  PetscBool ismoose;
174  DM dm = PETSC_NULL;
175 
176  // Initialize the part of the DM package that's packaged with Moose; in the PETSc source tree this
177  // call would be in DMInitializePackage()
179  CHKERRABORT(nl.comm().get(), ierr);
180  // Create and set up the DM that will consume the split options and deal with block matrices.
181  PetscNonlinearSolver<Number> * petsc_solver =
182  dynamic_cast<PetscNonlinearSolver<Number> *>(nl.nonlinearSolver());
183  SNES snes = petsc_solver->snes();
184  // if there exists a DMMoose object, not to recreate a new one
185  ierr = SNESGetDM(snes, &dm);
186  CHKERRABORT(nl.comm().get(), ierr);
187  if (dm)
188  {
189  ierr = PetscObjectTypeCompare((PetscObject)dm, DMMOOSE, &ismoose);
190  CHKERRABORT(nl.comm().get(), ierr);
191  if (ismoose)
192  return;
193  }
194  ierr = DMCreateMoose(nl.comm().get(), nl, &dm);
195  CHKERRABORT(nl.comm().get(), ierr);
196  ierr = DMSetFromOptions(dm);
197  CHKERRABORT(nl.comm().get(), ierr);
198  ierr = DMSetUp(dm);
199  CHKERRABORT(nl.comm().get(), ierr);
200  ierr = SNESSetDM(snes, dm);
201  CHKERRABORT(nl.comm().get(), ierr);
202  ierr = DMDestroy(&dm);
203  CHKERRABORT(nl.comm().get(), ierr);
204 // We temporarily comment out this updating function because
205 // we lack an approach to check if the problem
206 // structure has been changed from the last iteration.
207 // The indices will be rebuilt for every timestep.
208 // TODO: figure out a way to check the structure changes of the
209 // matrix
210 // ierr = SNESSetUpdate(snes,SNESUpdateDMMoose);
211 // CHKERRABORT(nl.comm().get(),ierr);
212 #endif
213 }
214 
215 void
217 {
218  // commandline options always win
219  // the options from a user commandline will overwrite the existing ones if any conflicts
220  { // Get any options specified on the command-line
221  int argc;
222  char ** args;
223 
224  PetscGetArgs(&argc, &args);
225 #if PETSC_VERSION_LESS_THAN(3, 7, 0)
226  PetscOptionsInsert(&argc, &args, NULL);
227 #else
228  PetscOptionsInsert(PETSC_NULL, &argc, &args, NULL);
229 #endif
230  }
231 }
232 
233 void
235 {
236  // Reference to the options stored in FEPRoblem
237  PetscOptions & petsc = problem.getPetscOptions();
238 
239  if (petsc.inames.size() != petsc.values.size())
240  mooseError("PETSc names and options are not the same length");
241 
242 #if PETSC_VERSION_LESS_THAN(3, 7, 0)
243  PetscOptionsClear();
244 #else
245  PetscOptionsClear(PETSC_NULL);
246 #endif
247 
248  setSolverOptions(problem.solverParams());
249 
250  // Add any additional options specified in the input file
251  for (const auto & flag : petsc.flags)
252  setSinglePetscOption(flag.rawName().c_str());
253  for (unsigned int i = 0; i < petsc.inames.size(); ++i)
254  setSinglePetscOption(petsc.inames[i], petsc.values[i]);
255 
256  // set up DM which is required if use a field split preconditioner
259 
261 }
262 
263 PetscErrorCode
265 {
266  char code[10] = {45, 45, 109, 111, 111, 115, 101};
267  const std::vector<std::string> argv = cmd_line->getArguments();
268  for (const auto & arg : argv)
269  {
270  if (arg == std::string(code, 10))
271  {
273  break;
274  }
275  }
276  return 0;
277 }
278 
279 PetscErrorCode
280 petscConverged(KSP ksp, PetscInt n, PetscReal rnorm, KSPConvergedReason * reason, void * ctx)
281 {
282  // Cast the context pointer coming from PETSc to an FEProblemBase& and
283  // get a reference to the System from it.
284  FEProblemBase & problem = *static_cast<FEProblemBase *>(ctx);
285 
286  // Let's be nice and always check PETSc error codes.
287  PetscErrorCode ierr = 0;
288 
289  // We want the default behavior of the KSPDefaultConverged test, but
290  // we don't want PETSc to die in that function with a CHKERRQ
291  // call... that is probably extremely unlikely/impossible, but just
292  // to be on the safe side, we push a different error handler before
293  // calling KSPDefaultConverged().
294  ierr = PetscPushErrorHandler(PetscReturnErrorHandler, /*void* ctx=*/PETSC_NULL);
295  CHKERRABORT(problem.comm().get(), ierr);
296 
297 #if PETSC_VERSION_LESS_THAN(3, 0, 0)
298  // Prior to PETSc 3.0.0, you could call KSPDefaultConverged with a NULL context
299  // pointer, as it was unused.
300  KSPDefaultConverged(ksp, n, rnorm, reason, PETSC_NULL);
301 #elif PETSC_VERSION_LESS_THAN(3, 5, 0)
302  // As of PETSc 3.0.0, you must call KSPDefaultConverged with a
303  // non-NULL context pointer which must be created with
304  // KSPDefaultConvergedCreate(), and destroyed with
305  // KSPDefaultConvergedDestroy().
306  void * default_ctx = NULL;
307  KSPDefaultConvergedCreate(&default_ctx);
308  KSPDefaultConverged(ksp, n, rnorm, reason, default_ctx);
309  KSPDefaultConvergedDestroy(default_ctx);
310 #else
311  // As of PETSc 3.5.0, use KSPConvergedDefaultXXX
312  void * default_ctx = NULL;
313  KSPConvergedDefaultCreate(&default_ctx);
314  KSPConvergedDefault(ksp, n, rnorm, reason, default_ctx);
315  KSPConvergedDefaultDestroy(default_ctx);
316 #endif
317 
318  // Pop the Error handler we pushed on the stack to go back
319  // to default PETSc error handling behavior.
320  ierr = PetscPopErrorHandler();
321  CHKERRABORT(problem.comm().get(), ierr);
322 
323  // Get tolerances from the KSP object
324  PetscReal rtol = 0.;
325  PetscReal atol = 0.;
326  PetscReal dtol = 0.;
327  PetscInt maxits = 0;
328  ierr = KSPGetTolerances(ksp, &rtol, &atol, &dtol, &maxits);
329  CHKERRABORT(problem.comm().get(), ierr);
330 
331  // Now do some additional MOOSE-specific tests...
332  std::string msg;
333  MooseLinearConvergenceReason moose_reason =
334  problem.checkLinearConvergence(msg, n, rnorm, rtol, atol, dtol, maxits);
335 
336  switch (moose_reason)
337  {
339  *reason = KSP_CONVERGED_RTOL;
340  break;
341 
343  *reason = KSP_CONVERGED_ITS;
344  break;
345 
347 #if PETSC_VERSION_LESS_THAN(3, 4, 0)
348  // Report divergence due to exceeding the divergence tolerance.
349  *reason = KSP_DIVERGED_DTOL;
350 #else
351  // KSP_DIVERGED_NANORINF was added in PETSc 3.4.0.
352  *reason = KSP_DIVERGED_NANORINF;
353 #endif
354  break;
355 #if !PETSC_VERSION_LESS_THAN(3, 6, 0) // A new convergence enum in PETSc 3.6
357 #if PETSC_VERSION_LESS_THAN(3, 11, 0) && PETSC_VERSION_RELEASE
358  *reason = KSP_DIVERGED_PCSETUP_FAILED;
359 #else
360  *reason = KSP_DIVERGED_PC_FAILED;
361 #endif
362  break;
363 #endif
364  default:
365  {
366  // If it's not either of the two specific cases we handle, just go
367  // with whatever PETSc decided in KSPDefaultConverged.
368  break;
369  }
370  }
371 
372  return 0;
373 }
374 
375 PetscErrorCode
377  PetscInt it,
378  PetscReal xnorm,
379  PetscReal snorm,
380  PetscReal fnorm,
381  SNESConvergedReason * reason,
382  void * ctx)
383 {
384  FEProblemBase & problem = *static_cast<FEProblemBase *>(ctx);
386 
387  // Let's be nice and always check PETSc error codes.
388  PetscErrorCode ierr = 0;
389 
390  // Temporary variables to store SNES tolerances. Usual C-style would be to declare
391  // but not initialize these... but it bothers me to leave anything uninitialized.
392  PetscReal atol = 0.; // absolute convergence tolerance
393  PetscReal rtol = 0.; // relative convergence tolerance
394  PetscReal stol = 0.; // convergence (step) tolerance in terms of the norm of the change in the
395  // solution between steps
396  PetscInt maxit = 0; // maximum number of iterations
397  PetscInt maxf = 0; // maximum number of function evaluations
398 
399  // Ask the SNES object about its tolerances.
400  ierr = SNESGetTolerances(snes, &atol, &rtol, &stol, &maxit, &maxf);
401  CHKERRABORT(problem.comm().get(), ierr);
402 
403  // Get current number of function evaluations done by SNES.
404  PetscInt nfuncs = 0;
405  ierr = SNESGetNumberFunctionEvals(snes, &nfuncs);
406  CHKERRABORT(problem.comm().get(), ierr);
407 
408  // Whether or not to force SNESSolve() take at least one iteration regardless of the initial
409  // residual norm
410  PetscBool force_iteration = PETSC_FALSE;
411 #if !PETSC_VERSION_LESS_THAN(3, 8, 4)
412  ierr = SNESGetForceIteration(snes, &force_iteration);
413  CHKERRABORT(problem.comm().get(), ierr);
414 #endif
415 
416 // See if SNESSetFunctionDomainError() has been called. Note:
417 // SNESSetFunctionDomainError() and SNESGetFunctionDomainError()
418 // were added in different releases of PETSc.
419 #if !PETSC_VERSION_LESS_THAN(3, 3, 0)
420  PetscBool domainerror;
421  ierr = SNESGetFunctionDomainError(snes, &domainerror);
422  CHKERRABORT(problem.comm().get(), ierr);
423  if (domainerror)
424  {
425  *reason = SNES_DIVERGED_FUNCTION_DOMAIN;
426  return 0;
427  }
428 #endif
429 
430  // Error message that will be set by the FEProblemBase.
431  std::string msg;
432 
433  // xnorm: 2-norm of current iterate
434  // snorm: 2-norm of current step
435  // fnorm: 2-norm of function at current iterate
436  MooseNonlinearConvergenceReason moose_reason =
437  problem.checkNonlinearConvergence(msg,
438  it,
439  xnorm,
440  snorm,
441  fnorm,
442  rtol,
443  stol,
444  atol,
445  nfuncs,
446  maxf,
447  force_iteration,
448  system._initial_residual_before_preset_bcs,
449  std::numeric_limits<Real>::max());
450 
451  if (msg.length() > 0)
452  PetscInfo(snes, msg.c_str());
453 
454  switch (moose_reason)
455  {
457  *reason = SNES_CONVERGED_ITERATING;
458  break;
459 
461  *reason = SNES_CONVERGED_FNORM_ABS;
462  break;
463 
465  *reason = SNES_CONVERGED_FNORM_RELATIVE;
466  break;
467 
469 #if PETSC_VERSION_LESS_THAN(3, 3, 0)
470  *reason = SNES_CONVERGED_PNORM_RELATIVE;
471 #else
472  *reason = SNES_CONVERGED_SNORM_RELATIVE;
473 #endif
474  break;
475 
477  *reason = SNES_DIVERGED_FUNCTION_COUNT;
478  break;
479 
481  *reason = SNES_DIVERGED_FNORM_NAN;
482  break;
483 
485 #if PETSC_VERSION_LESS_THAN(3, 2, 0)
486  *reason = SNES_DIVERGED_LS_FAILURE;
487 #else
488  *reason = SNES_DIVERGED_LINE_SEARCH;
489 #endif
490  break;
491  }
492 
493  return 0;
494 }
495 
496 PCSide
498 {
499  switch (pcs)
500  {
501  case Moose::PCS_LEFT:
502  return PC_LEFT;
503  case Moose::PCS_RIGHT:
504  return PC_RIGHT;
506  return PC_SYMMETRIC;
507  default:
508  mooseError("Unknown PC side requested.");
509  break;
510  }
511 }
512 
513 KSPNormType
515 {
516  switch (kspnorm)
517  {
518  case Moose::KSPN_NONE:
519  return KSP_NORM_NONE;
521  return KSP_NORM_PRECONDITIONED;
523  return KSP_NORM_UNPRECONDITIONED;
524  case Moose::KSPN_NATURAL:
525  return KSP_NORM_NATURAL;
526  case Moose::KSPN_DEFAULT:
527  return KSP_NORM_DEFAULT;
528  default:
529  mooseError("Unknown KSP norm type requested.");
530  break;
531  }
532 }
533 
534 void
536 {
538 
539  KSPSetNormType(ksp, getPetscKSPNormType(nl.getMooseKSPNormType()));
540 }
541 
542 void
544 {
546 
547 #if PETSC_VERSION_LESS_THAN(3, 2, 0)
548  // pc_side is NOT set, PETSc will make the decision
549  // PETSc 3.1.x-
551  KSPSetPreconditionerSide(ksp, getPetscPCSide(nl.getPCSide()));
552 #else
553  // PETSc 3.2.x+
555  KSPSetPCSide(ksp, getPetscPCSide(nl.getPCSide()));
556 #endif
557 }
558 
559 void
561 {
563 
564 #if PETSC_VERSION_LESS_THAN(3, 0, 0)
565  // PETSc 2.3.3-
566  KSPSetConvergenceTest(ksp, petscConverged, &problem);
567 #else
568  // PETSc 3.0.0+
569 
570  // In 3.0.0, the context pointer must actually be used, and the
571  // final argument to KSPSetConvergenceTest() is a pointer to a
572  // routine for destroying said private data context. In this case,
573  // we use the default context provided by PETSc in addition to
574  // a few other tests.
575  {
576  PetscErrorCode ierr = KSPSetConvergenceTest(ksp, petscConverged, &problem, PETSC_NULL);
577  CHKERRABORT(nl.comm().get(), ierr);
578  }
579 #endif
580 
581  auto & es = problem.es();
582 
583  PetscReal rtol = es.parameters.get<Real>("linear solver tolerance");
584  PetscReal atol = es.parameters.get<Real>("linear solver absolute step tolerance");
585 
586  // MOOSE defaults this to -1 for some dumb reason
587  if (atol < 0)
588  atol = 1e-50;
589 
590  PetscReal maxits = es.parameters.get<unsigned int>("linear solver maximum iterations");
591 
592  // 1e100 is because we don't use divtol currently
593  KSPSetTolerances(ksp, rtol, atol, 1e100, maxits);
594 
595  petscSetDefaultPCSide(problem, ksp);
596 
597  petscSetDefaultKSPNormType(problem, ksp);
598 }
599 
600 void
602 {
603  // dig out Petsc solver
605  PetscNonlinearSolver<Number> * petsc_solver =
606  dynamic_cast<PetscNonlinearSolver<Number> *>(nl.nonlinearSolver());
607  SNES snes = petsc_solver->snes();
608  KSP ksp;
609  SNESGetKSP(snes, &ksp);
610 
611  SNESSetMaxLinearSolveFailures(snes, 1000000);
612 
613 #if PETSC_VERSION_LESS_THAN(3, 0, 0)
614  // PETSc 2.3.3-
615  SNESSetConvergenceTest(snes, petscNonlinearConverged, &problem);
616 #else
617  // PETSc 3.0.0+
618 
619  // In 3.0.0, the context pointer must actually be used, and the
620  // final argument to KSPSetConvergenceTest() is a pointer to a
621  // routine for destroying said private data context. In this case,
622  // we use the default context provided by PETSc in addition to
623  // a few other tests.
624  {
625  auto ierr = SNESSetConvergenceTest(snes, petscNonlinearConverged, &problem, PETSC_NULL);
626  CHKERRABORT(nl.comm().get(), ierr);
627  }
628 #endif
629 
630  petscSetKSPDefaults(problem, ksp);
631 }
632 
633 void
634 storePetscOptions(FEProblemBase & fe_problem, const InputParameters & params)
635 {
636  // Note: Options set in the Preconditioner block will override those set in the Executioner block
637  if (params.isParamValid("solve_type") && !params.isParamValid("_use_eigen_value"))
638  {
639  // Extract the solve type
640  const std::string & solve_type = params.get<MooseEnum>("solve_type");
641  fe_problem.solverParams()._type = Moose::stringToEnum<Moose::SolveType>(solve_type);
642  }
643 
644  if (params.isParamValid("line_search"))
645  {
646  MooseEnum line_search = params.get<MooseEnum>("line_search");
647  if (fe_problem.solverParams()._line_search == Moose::LS_INVALID || line_search != "default")
648  {
649  Moose::LineSearchType enum_line_search =
650  Moose::stringToEnum<Moose::LineSearchType>(line_search);
651  fe_problem.solverParams()._line_search = enum_line_search;
652  if (enum_line_search == LS_CONTACT)
653  {
654  NonlinearImplicitSystem * nl_system =
655  dynamic_cast<NonlinearImplicitSystem *>(&fe_problem.getNonlinearSystemBase().system());
656  if (!nl_system)
657  mooseError("You've requested a line search but you must be solving an EigenProblem. "
658  "These two things are not consistent.");
659  PetscNonlinearSolver<Real> * petsc_nonlinear_solver =
660  dynamic_cast<PetscNonlinearSolver<Real> *>(nl_system->nonlinear_solver.get());
661  if (!petsc_nonlinear_solver)
662  mooseError("Currently the contact line search is only implemented through Petsc, so you "
663  "must use Petsc as your non-linear solver.");
664  petsc_nonlinear_solver->linesearch_object =
665  libmesh_make_unique<ComputeLineSearchObjectWrapper>(fe_problem);
666  }
667  }
668  }
669 
670  if (params.isParamValid("mffd_type"))
671  {
672  MooseEnum mffd_type = params.get<MooseEnum>("mffd_type");
673  fe_problem.solverParams()._mffd_type = Moose::stringToEnum<Moose::MffdType>(mffd_type);
674  }
675 
676  // The parameters contained in the Action
677  const MultiMooseEnum & petsc_options = params.get<MultiMooseEnum>("petsc_options");
678  const MultiMooseEnum & petsc_options_inames = params.get<MultiMooseEnum>("petsc_options_iname");
679  const std::vector<std::string> & petsc_options_values =
680  params.get<std::vector<std::string>>("petsc_options_value");
681 
682  // A reference to the PetscOptions object that contains the settings that will be used in the
683  // solve
685 
686  // Update the PETSc single flags
687  for (const auto & option : petsc_options)
688  {
695  if (option == "-log_summary")
696  mooseError("The PETSc option \"-log_summary\" can only be used on the command line. Please "
697  "remove it from the input file");
698 
699  // Warn about superseded PETSc options (Note: -snes is not a REAL option, but people used it in
700  // their input files)
701  else
702  {
703  std::string help_string;
704  if (option == "-snes" || option == "-snes_mf" || option == "-snes_mf_operator")
705  help_string = "Please set the solver type through \"solve_type\".";
706  else if (option == "-ksp_monitor")
707  help_string = "Please use \"Outputs/print_linear_residuals=true\"";
708 
709  if (help_string != "")
710  mooseWarning("The PETSc option ",
711  std::string(option),
712  " should not be used directly in a MOOSE input file. ",
713  help_string);
714  }
715 
716  // Update the stored items, but do not create duplicates
717  if (!po.flags.contains(option))
718  po.flags.push_back(option);
719  }
720 
721  // Check that the name value pairs are sized correctly
722  if (petsc_options_inames.size() != petsc_options_values.size())
723  mooseError("PETSc names and options are not the same length");
724 
725  // Setup the name value pairs
726  bool boomeramg_found = false;
727  bool strong_threshold_found = false;
728 #if !PETSC_VERSION_LESS_THAN(3, 7, 0)
729  bool superlu_dist_found = false;
730  bool fact_pattern_found = false;
731  bool tiny_pivot_found = false;
732 #endif
733  std::string pc_description = "";
734  for (unsigned int i = 0; i < petsc_options_inames.size(); i++)
735  {
736  // Do not add duplicate settings
737  if (find(po.inames.begin(), po.inames.end(), petsc_options_inames[i]) == po.inames.end())
738  {
739 #if !PETSC_VERSION_LESS_THAN(3, 9, 0)
740  if (petsc_options_inames[i] == "-pc_factor_mat_solver_package")
741  po.inames.push_back("-pc_factor_mat_solver_type");
742  else
743  po.inames.push_back(petsc_options_inames[i]);
744 #else
745  if (petsc_options_inames[i] == "-pc_factor_mat_solver_type")
746  po.inames.push_back("-pc_factor_mat_solver_package");
747  else
748  po.inames.push_back(petsc_options_inames[i]);
749 #endif
750  po.values.push_back(petsc_options_values[i]);
751 
752  // Look for a pc description
753  if (petsc_options_inames[i] == "-pc_type" || petsc_options_inames[i] == "-pc_sub_type" ||
754  petsc_options_inames[i] == "-pc_hypre_type")
755  pc_description += petsc_options_values[i] + ' ';
756 
757  // This special case is common enough that we'd like to handle it for the user.
758  if (petsc_options_inames[i] == "-pc_hypre_type" && petsc_options_values[i] == "boomeramg")
759  boomeramg_found = true;
760  if (petsc_options_inames[i] == "-pc_hypre_boomeramg_strong_threshold")
761  strong_threshold_found = true;
762 #if !PETSC_VERSION_LESS_THAN(3, 7, 0)
763  if ((petsc_options_inames[i] == "-pc_factor_mat_solver_package" ||
764  petsc_options_inames[i] == "-pc_factor_mat_solver_type") &&
765  petsc_options_values[i] == "superlu_dist")
766  superlu_dist_found = true;
767  if (petsc_options_inames[i] == "-mat_superlu_dist_fact")
768  fact_pattern_found = true;
769  if (petsc_options_inames[i] == "-mat_superlu_dist_replacetinypivot")
770  tiny_pivot_found = true;
771 #endif
772  }
773  else
774  {
775  for (unsigned int j = 0; j < po.inames.size(); j++)
776  if (po.inames[j] == petsc_options_inames[i])
777  po.values[j] = petsc_options_values[i];
778  }
779  }
780 
781  // When running a 3D mesh with boomeramg, it is almost always best to supply a strong threshold
782  // value
783  // We will provide that for the user here if they haven't supplied it themselves.
784  if (boomeramg_found && !strong_threshold_found && fe_problem.mesh().dimension() == 3)
785  {
786  po.inames.push_back("-pc_hypre_boomeramg_strong_threshold");
787  po.values.push_back("0.7");
788  pc_description += "strong_threshold: 0.7 (auto)";
789  }
790 
791 #if !PETSC_VERSION_LESS_THAN(3, 7, 0)
792  // In PETSc-3.7.{0--4}, there is a bug when using superlu_dist, and we have to use
793  // SamePattern_SameRowPerm, otherwise we use whatever we have in PETSc
794  if (superlu_dist_found && !fact_pattern_found)
795  {
796  po.inames.push_back("-mat_superlu_dist_fact");
797 #if PETSC_VERSION_LESS_THAN(3, 7, 5)
798  po.values.push_back("SamePattern_SameRowPerm");
799  pc_description += "mat_superlu_dist_fact: SamePattern_SameRowPerm ";
800 #else
801  po.values.push_back("SamePattern");
802  pc_description += "mat_superlu_dist_fact: SamePattern ";
803 #endif
804  }
805 
806  // restore this superlu option
807  if (superlu_dist_found && !tiny_pivot_found)
808  {
809  po.inames.push_back("-mat_superlu_dist_replacetinypivot");
810  po.values.push_back("1");
811  pc_description += " mat_superlu_dist_replacetinypivot: true ";
812  }
813 #endif
814  // Set Preconditioner description
815  po.pc_description = pc_description;
816 }
817 
818 std::set<std::string>
820 {
821 #if PETSC_VERSION_LESS_THAN(3, 3, 0)
822  return {"default", "cubic", "quadratic", "none", "basic", "basicnonorms"};
823 #else
824  return {"default", "shell", "none", "basic", "l2", "bt", "cp"};
825 #endif
826 }
827 
830 {
832 
833  MooseEnum solve_type("PJFNK JFNK NEWTON FD LINEAR");
834  params.addParam<MooseEnum>("solve_type",
835  solve_type,
836  "PJFNK: Preconditioned Jacobian-Free Newton Krylov "
837  "JFNK: Jacobian-Free Newton Krylov "
838  "NEWTON: Full Newton Solve "
839  "FD: Use finite differences to compute Jacobian "
840  "LINEAR: Solving a linear problem");
841 
842  MooseEnum mffd_type("wp ds", "wp");
843  params.addParam<MooseEnum>("mffd_type",
844  mffd_type,
845  "Specifies the finite differencing type for "
846  "Jacobian-free solve types. Note that the "
847  "default is wp (for Walker and Pernice).");
848 
849  params.addParam<MultiMooseEnum>(
850  "petsc_options", getCommonPetscFlags(), "Singleton PETSc options");
851  params.addParam<MultiMooseEnum>(
852  "petsc_options_iname", getCommonPetscKeys(), "Names of PETSc name/value pairs");
853  params.addParam<std::vector<std::string>>(
854  "petsc_options_value",
855  "Values of PETSc name/value pairs (must correspond with \"petsc_options_iname\"");
856  return params;
857 }
858 
861 {
862  return MultiMooseEnum(
863  "-dm_moose_print_embedding -dm_view -ksp_converged_reason -ksp_gmres_modifiedgramschmidt "
864  "-ksp_monitor -ksp_monitor_snes_lg-snes_ksp_ew -ksp_snes_ew -snes_converged_reason "
865  "-snes_ksp -snes_ksp_ew -snes_linesearch_monitor -snes_mf -snes_mf_operator -snes_monitor "
866  "-snes_test_display -snes_view",
867  "",
868  true);
869 }
870 
873 {
874  return MultiMooseEnum("-ksp_atol -ksp_gmres_restart -ksp_max_it -ksp_pc_side -ksp_rtol "
875  "-ksp_type -mat_fd_coloring_err -mat_fd_type -mat_mffd_type "
876  "-pc_asm_overlap -pc_factor_levels "
877  "-pc_factor_mat_ordering_type -pc_hypre_boomeramg_grid_sweeps_all "
878  "-pc_hypre_boomeramg_max_iter "
879  "-pc_hypre_boomeramg_strong_threshold -pc_hypre_type -pc_type -snes_atol "
880  "-snes_linesearch_type "
881  "-snes_ls -snes_max_it -snes_rtol -snes_type -sub_ksp_type -sub_pc_type",
882  "",
883  true);
884 }
885 
886 void
887 setSinglePetscOption(const std::string & name, const std::string & value)
888 {
889  PetscErrorCode ierr;
890 
891 #if PETSC_VERSION_LESS_THAN(3, 7, 0)
892  ierr = PetscOptionsSetValue(name.c_str(), value == "" ? PETSC_NULL : value.c_str());
893 #else
894  // PETSc 3.7.0 and later version. First argument is the options
895  // database to use, NULL indicates the default global database.
896  ierr = PetscOptionsSetValue(PETSC_NULL, name.c_str(), value == "" ? PETSC_NULL : value.c_str());
897 #endif
898 
899  // Not convenient to use the usual error checking macro, because we
900  // don't have a specific communicator in this helper function.
901  if (ierr)
902  mooseError("Error setting PETSc option: ", name);
903 }
904 
905 void
906 colorAdjacencyMatrix(PetscScalar * adjacency_matrix,
907  unsigned int size,
908  unsigned int colors,
909  std::vector<unsigned int> & vertex_colors,
910  const char * coloring_algorithm)
911 {
912  // Mat A will be a dense matrix from the incoming data structure
913  Mat A;
914  MatCreate(MPI_COMM_SELF, &A);
915  MatSetSizes(A, size, size, size, size);
916  MatSetType(A, MATSEQDENSE);
917  // PETSc requires a non-const data array to populate the matrix
918  MatSeqDenseSetPreallocation(A, adjacency_matrix);
919  MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
920  MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
921 
922  // Convert A to a sparse matrix
923  MatConvert(A,
924  MATAIJ,
925 #if PETSC_VERSION_LESS_THAN(3, 7, 0)
926  MAT_REUSE_MATRIX,
927 #else
928  MAT_INPLACE_MATRIX,
929 #endif
930  &A);
931 
932  ISColoring iscoloring;
933 #if PETSC_VERSION_LESS_THAN(3, 5, 0)
934  MatGetColoring(A, coloring_algorithm, &iscoloring);
935 #else
936  MatColoring mc;
937  MatColoringCreate(A, &mc);
938  MatColoringSetType(mc, coloring_algorithm);
939  MatColoringSetMaxColors(mc, static_cast<PetscInt>(colors));
940 
941  // Petsc normally colors by distance two (neighbors of neighbors), we just want one
942  MatColoringSetDistance(mc, 1);
943  MatColoringSetFromOptions(mc);
944  MatColoringApply(mc, &iscoloring);
945 #endif
946 
947  PetscInt nn;
948  IS * is;
949  ISColoringGetIS(iscoloring, &nn, &is);
950 
951  if (nn > static_cast<PetscInt>(colors))
952  throw std::runtime_error("Not able to color with designated number of colors");
953 
954  for (int i = 0; i < nn; i++)
955  {
956  PetscInt isize;
957  const PetscInt * indices;
958  ISGetLocalSize(is[i], &isize);
959  ISGetIndices(is[i], &indices);
960  for (int j = 0; j < isize; j++)
961  {
962  mooseAssert(indices[j] < static_cast<PetscInt>(vertex_colors.size()), "Index out of bounds");
963  vertex_colors[indices[j]] = i;
964  }
965  ISRestoreIndices(is[i], &indices);
966  }
967 
968  MatDestroy(&A);
969 #if !PETSC_VERSION_LESS_THAN(3, 5, 0)
970  MatColoringDestroy(&mc);
971 #endif
972  ISColoringDestroy(&iscoloring);
973 }
974 
975 } // Namespace PetscSupport
976 } // Namespace MOOSE
977 
978 #endif // LIBMESH_HAVE_PETSC
MultiMooseEnum getCommonPetscKeys()
A helper function to produce a MultiMooseEnum with commonly used PETSc iname options (keys in key-val...
Definition: PetscSupport.C:872
std::vector< std::string > values
Values for PETSc key-value pairs.
Definition: PetscSupport.h:48
Moose::PetscSupport::PetscOptions & getPetscOptions()
Retrieve a writable reference the PETSc options (used by PetscSupport)
KSPNormType getPetscKSPNormType(Moose::MooseKSPNormType kspnorm)
Definition: PetscSupport.C:514
void addPetscOptionsFromCommandline()
Definition: PetscSupport.C:216
void storePetscOptions(FEProblemBase &fe_problem, const InputParameters &params)
Stores the PETSc options supplied from the InputParameters with MOOSE.
Definition: PetscSupport.C:634
Full Newton Solve.
Definition: MooseTypes.h:593
SolverParams & solverParams()
Get the solver parameters.
static void petscSetupOutput()
Output string for setting up PETSC output.
Definition: Console.C:674
NonlinearSystemBase & getNonlinearSystemBase()
virtual NonlinearSolver< Number > * nonlinearSolver()=0
std::set< std::string > getPetscValidLineSearches()
Returns the valid petsc line search options as a set of strings.
Definition: PetscSupport.C:819
Moose::PCSideType getPCSide()
bool haveFieldSplitPreconditioner() const
Moose::LineSearchType _line_search
Definition: SolverParams.h:20
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:207
void petscSetDefaults(FEProblemBase &problem)
Sets the default options for PETSc.
Definition: PetscSupport.C:601
void mooseWarning(Args &&... args)
Emit a warning message with the given stringified, concatenated args.
Definition: MooseError.h:219
A struct for storing the various types of petsc options and values.
Definition: PetscSupport.h:39
MultiMooseEnum flags
Single value PETSc options (flags)
Definition: PetscSupport.h:51
PetscBool ismoose
PetscErrorCode DMCreateMoose(MPI_Comm, NonlinearSystemBase &, DM *)
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
void petscSetDefaultKSPNormType(FEProblemBase &problem, KSP ksp)
Definition: PetscSupport.C:535
Solving a linear problem.
Definition: MooseTypes.h:595
void petscSetKSPDefaults(FEProblemBase &problem, KSP ksp)
Set the default options for a KSP.
Definition: PetscSupport.C:560
PetscErrorCode petscConverged(KSP ksp, PetscInt n, PetscReal rnorm, KSPConvergedReason *reason, void *ctx)
Definition: PetscSupport.C:280
MffdType
Type of the matrix-free finite-differencing parameter.
Definition: MooseTypes.h:695
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
InputParameters emptyInputParameters()
This class wraps provides and tracks access to command line parameters.
Definition: CommandLine.h:30
virtual EquationSystems & es() override
MooseKSPNormType
Norm type for converge test.
Definition: MooseTypes.h:577
NonlinearSystemBase * nl
Nonlinear system to be solved.
PCSide getPetscPCSide(Moose::PCSideType pcs)
Definition: PetscSupport.C:497
MooseLinearConvergenceReason
virtual MooseLinearConvergenceReason checkLinearConvergence(std::string &msg, const PetscInt n, const Real rnorm, const Real rtol, const Real atol, const Real dtol, const PetscInt maxits)
Check for convergence of the linear solution.
Use whatever we have in PETSc.
Definition: MooseTypes.h:571
nl system()
MultiMooseEnum getCommonPetscFlags()
A helper function to produce a MultiMooseEnum with commonly used PETSc single options (flags) ...
Definition: PetscSupport.C:860
Moose::MffdType _mffd_type
Definition: SolverParams.h:21
LineSearchType
Type of the line search.
Definition: MooseTypes.h:671
static PetscErrorCode Mat * A
virtual unsigned int dimension() const
Returns MeshBase::mesh_dimsension(), (not MeshBase::spatial_dimension()!) of the underlying libMesh m...
Definition: MooseMesh.C:2133
Use whatever we have in PETSc.
Definition: MooseTypes.h:583
Jacobian-Free Newton Krylov.
Definition: MooseTypes.h:592
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
Moose::SolveType _type
Definition: SolverParams.h:19
Use finite differences to compute Jacobian.
Definition: MooseTypes.h:594
void setSolverOptions(SolverParams &solver_params)
Definition: PetscSupport.C:125
MooseNonlinearConvergenceReason
Enumeration for nonlinear convergence reasons.
Definition: FEProblemBase.h:93
InputParameters getPetscValidParams()
Returns the PETSc options that are common between Executioners and Preconditioners.
Definition: PetscSupport.C:829
PetscErrorCode petscSetupOutput(CommandLine *cmd_line)
Definition: PetscSupport.C:264
virtual System & system() override
Get the reference to the libMesh system.
PetscInt n
std::string stringify(const LineSearchType &t)
Definition: PetscSupport.C:74
means not set
Definition: MooseTypes.h:697
void colorAdjacencyMatrix(PetscScalar *adjacency_matrix, unsigned int size, unsigned int colors, std::vector< unsigned int > &vertex_colors, const char *coloring_algorithm)
This method takes an adjacency matrix, and a desired number of colors and applies a graph coloring al...
Definition: PetscSupport.C:906
virtual MooseMesh & mesh() override
PCSideType
Preconditioning side.
Definition: MooseTypes.h:566
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
Moose::MooseKSPNormType getMooseKSPNormType()
ierr
Preconditioned Jacobian-Free Newton Krylov.
Definition: MooseTypes.h:591
Definition: Moose.h:112
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
means not set
Definition: MooseTypes.h:673
PetscErrorCode petscNonlinearConverged(SNES snes, PetscInt it, PetscReal xnorm, PetscReal snorm, PetscReal fnorm, SNESConvergedReason *reason, void *ctx)
Definition: PetscSupport.C:376
void petscSetOptions(FEProblemBase &problem)
A function for setting the PETSc options in PETSc from the options supplied to MOOSE.
Definition: PetscSupport.C:234
PetscErrorCode DMMooseRegisterAll()
void petscSetupDM(NonlinearSystemBase &nl)
Definition: PetscSupport.C:169
static PetscErrorCode Vec Mat Mat void * ctx
void setSinglePetscOption(const std::string &name, const std::string &value="")
A wrapper function for dealing with different versions of PetscOptionsSetValue.
Definition: PetscSupport.C:887
const std::vector< std::string > & getArguments()
Return the raw argv arguments as a vector.
Definition: CommandLine.h:60
std::vector< std::string > inames
Keys for PETSc key-value pairs.
Definition: PetscSupport.h:45
void petscSetDefaultPCSide(FEProblemBase &problem, KSP ksp)
Definition: PetscSupport.C:543
virtual MooseNonlinearConvergenceReason checkNonlinearConvergence(std::string &msg, const PetscInt it, const Real xnorm, const Real snorm, const Real fnorm, const Real rtol, const Real stol, const Real abstol, const PetscInt nfuncs, const PetscInt max_funcs, const PetscBool force_iteration, const Real initial_residual_before_preset_bcs, const Real div_threshold)
Check for converence of the nonlinear solution.
bool isParamValid(const std::string &name) const
This method returns parameters that have been initialized in one fashion or another, i.e.