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