https://mooseframework.inl.gov
FEProblemSolve.C
Go to the documentation of this file.
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 "FEProblemSolve.h"
11 
12 #include "FEProblem.h"
13 #include "NonlinearSystemBase.h"
14 #include "LinearSystem.h"
15 #include "Convergence.h"
16 #include "Executioner.h"
17 
18 std::set<std::string> const FEProblemSolve::_moose_line_searches = {"contact", "project"};
19 
20 const std::set<std::string> &
22 {
23  return _moose_line_searches;
24 }
25 
28 {
30 
31  params.addParam<unsigned int>("nl_max_its", 50, "Max Nonlinear Iterations");
32  params.addParam<unsigned int>("nl_forced_its", 0, "The Number of Forced Nonlinear Iterations");
33  params.addParam<unsigned int>("nl_max_funcs", 10000, "Max Nonlinear solver function evaluations");
34  params.addParam<Real>("nl_abs_tol", 1.0e-50, "Nonlinear Absolute Tolerance");
35  params.addParam<Real>("nl_rel_tol", 1.0e-8, "Nonlinear Relative Tolerance");
36  params.addParam<Real>(
37  "nl_div_tol",
38  1.0e10,
39  "Nonlinear Relative Divergence Tolerance. A negative value disables this check.");
40  params.addParam<Real>(
41  "nl_abs_div_tol",
42  1.0e50,
43  "Nonlinear Absolute Divergence Tolerance. A negative value disables this check.");
44  params.addParam<Real>("nl_abs_step_tol", 0., "Nonlinear Absolute step Tolerance");
45  params.addParam<Real>("nl_rel_step_tol", 0., "Nonlinear Relative step Tolerance");
46  params.addParam<unsigned int>("n_max_nonlinear_pingpong",
47  100,
48  "The maximum number of times the nonlinear residual can ping pong "
49  "before requesting halting the current evaluation and requesting "
50  "timestep cut for transient simulations");
51 
52  params.addParamNamesToGroup(
53  "nl_max_its nl_forced_its nl_max_funcs nl_abs_tol nl_rel_tol "
54  "nl_rel_step_tol nl_abs_step_tol nl_div_tol nl_abs_div_tol n_max_nonlinear_pingpong",
55  "Nonlinear Solver");
56 
57  return params;
58 }
59 
62 {
65 
66  params.addParam<std::vector<std::vector<std::string>>>(
67  "splitting",
68  {},
69  "Top-level splitting defining a hierarchical decomposition into "
70  "subsystems to help the solver. Outer-vector of this vector-of-vector parameter correspond "
71  "to each nonlinear system.");
72 
73  std::set<std::string> line_searches = mooseLineSearches();
74 
75  std::set<std::string> alias_line_searches = {"default", "none", "basic"};
76  line_searches.insert(alias_line_searches.begin(), alias_line_searches.end());
77  std::set<std::string> petsc_line_searches = Moose::PetscSupport::getPetscValidLineSearches();
78  line_searches.insert(petsc_line_searches.begin(), petsc_line_searches.end());
79  std::string line_search_string = Moose::stringify(line_searches, " ");
80  MooseEnum line_search(line_search_string, "default");
81  std::string addtl_doc_str(" (Note: none = basic)");
82  params.addParam<MooseEnum>(
83  "line_search", line_search, "Specifies the line search type" + addtl_doc_str);
84  MooseEnum line_search_package("petsc moose", "petsc");
85  params.addParam<MooseEnum>("line_search_package",
86  line_search_package,
87  "The solver package to use to conduct the line-search");
88 
89  params.addParam<unsigned>("contact_line_search_allowed_lambda_cuts",
90  2,
91  "The number of times lambda is allowed to be cut in half in the "
92  "contact line search. We recommend this number be roughly bounded by 0 "
93  "<= allowed_lambda_cuts <= 3");
94  params.addParam<Real>("contact_line_search_ltol",
95  "The linear relative tolerance to be used while the contact state is "
96  "changing between non-linear iterations. We recommend that this tolerance "
97  "be looser than the standard linear tolerance");
98 
100  params.addParam<Real>("l_tol", 1.0e-5, "Linear Relative Tolerance");
101  params.addParam<Real>("l_abs_tol", 1.0e-50, "Linear Absolute Tolerance");
102  params.addParam<unsigned int>("l_max_its", 10000, "Max Linear Iterations");
103  params.addParam<std::vector<ConvergenceName>>(
104  "nonlinear_convergence",
105  "Name of the Convergence object(s) to use to assess convergence of the "
106  "nonlinear system(s) solve. If not provided, the default Convergence "
107  "associated with the Problem will be constructed internally.");
108  params.addParam<std::vector<ConvergenceName>>(
109  "linear_convergence",
110  "Name of the Convergence object(s) to use to assess convergence of the "
111  "linear system(s) solve. If not provided, the linear solver tolerance parameters are used");
112  params.addParam<bool>(
113  "snesmf_reuse_base",
114  true,
115  "Specifies whether or not to reuse the base vector for matrix-free calculation");
116  params.addParam<bool>(
117  "skip_exception_check", false, "Specifies whether or not to skip exception check");
118  params.addParam<bool>("compute_initial_residual_before_preset_bcs",
119  false,
120  "Use the residual norm computed *before* solution modifying objects like "
121  "preset BCs are imposed in relative convergence check.");
122  params.deprecateParam(
123  "compute_initial_residual_before_preset_bcs", "use_pre_SMO_residual", "12/31/2024");
124  params.addParam<bool>(
125  "use_pre_SMO_residual",
126  false,
127  "Compute the pre-SMO residual norm and use it in the relative convergence check. The "
128  "pre-SMO residual is computed at the begining of the time step before solution-modifying "
129  "objects are executed. Solution-modifying objects include preset BCs, constraints, "
130  "predictors, etc.");
131  params.addParam<bool>("automatic_scaling", "Whether to use automatic scaling for the variables.");
132  params.addParam<std::vector<bool>>(
133  "compute_scaling_once",
134  {true},
135  "Whether the scaling factors should only be computed once at the beginning of the simulation "
136  "through an extra Jacobian evaluation. If this is set to false, then the scaling factors "
137  "will be computed during an extra Jacobian evaluation at the beginning of every time step. "
138  "Vector entries correspond to each nonlinear system.");
139  params.addParam<std::vector<bool>>(
140  "off_diagonals_in_auto_scaling",
141  {false},
142  "Whether to consider off-diagonals when determining automatic scaling factors. Vector "
143  "entries correspond to each nonlinear system.");
144  params.addRangeCheckedParam<std::vector<Real>>(
145  "resid_vs_jac_scaling_param",
146  {0},
147  "0<=resid_vs_jac_scaling_param<=1",
148  "A parameter that indicates the weighting of the residual vs the Jacobian in determining "
149  "variable scaling parameters. A value of 1 indicates pure residual-based scaling. A value of "
150  "0 indicates pure Jacobian-based scaling. Vector entries correspond to each nonlinear "
151  "system.");
152  params.addParam<std::vector<std::vector<std::vector<std::string>>>>(
153  "scaling_group_variables",
154  "Name of variables that are grouped together for determining scale factors. (Multiple "
155  "groups can be provided, separated by semicolon). Vector entries correspond to each "
156  "nonlinear system.");
157  params.addParam<std::vector<std::vector<std::string>>>(
158  "ignore_variables_for_autoscaling",
159  "List of variables that do not participate in autoscaling. Vector entries correspond to each "
160  "nonlinear system.");
161  params.addRangeCheckedParam<unsigned int>(
162  "num_grids",
163  1,
164  "num_grids>0",
165  "The number of grids to use for a grid sequencing algorithm. This includes the final grid, "
166  "so num_grids = 1 indicates just one solve in a time-step");
167  params.addParam<std::vector<bool>>("residual_and_jacobian_together",
168  {false},
169  "Whether to compute the residual and Jacobian together. "
170  "Vector entries correspond to each nonlinear system.");
171 
172  params.addParam<bool>("reuse_preconditioner",
173  false,
174  "If true reuse the previously calculated "
175  "preconditioner for the linearized "
176  "system across multiple solves "
177  "spanning nonlinear iterations and time steps. "
178  "The preconditioner resets as controlled by "
179  "reuse_preconditioner_max_linear_its");
180  params.addParam<unsigned int>("reuse_preconditioner_max_linear_its",
181  25,
182  "Reuse the previously calculated "
183  "preconditioner for the linear system "
184  "until the number of linear iterations "
185  "exceeds this number");
186 
187  // Multi-system fixed point
188  // Defaults to false because of the difficulty of defining a good multi-system convergence
189  // criterion, unless we add a default one to the simulation?
190  params.addParam<bool>(
191  "multi_system_fixed_point",
192  false,
193  "Whether to perform fixed point (Picard) iterations between the nonlinear systems.");
194  params.addParam<ConvergenceName>(
195  "multi_system_fixed_point_convergence",
196  "Convergence object to determine the convergence of the multi-system fixed point iteration. "
197  "If unspecified, defaults to checking that every system is converged (based on their own "
198  "convergence criterion)");
199 
200  params.addParamNamesToGroup("l_tol l_abs_tol l_max_its reuse_preconditioner "
201  "reuse_preconditioner_max_linear_its",
202  "Linear Solver");
203  params.addParamNamesToGroup(
204  "solve_type snesmf_reuse_base use_pre_SMO_residual "
205  "num_grids residual_and_jacobian_together splitting nonlinear_convergence linear_convergence",
206  "Nonlinear Solver");
207  params.addParamNamesToGroup(
208  "automatic_scaling compute_scaling_once off_diagonals_in_auto_scaling "
209  "scaling_group_variables resid_vs_jac_scaling_param ignore_variables_for_autoscaling",
210  "Solver variable scaling");
211  params.addParamNamesToGroup("line_search line_search_package contact_line_search_ltol "
212  "contact_line_search_allowed_lambda_cuts",
213  "Solver line search");
214  params.addParamNamesToGroup("multi_system_fixed_point multi_system_fixed_point_convergence",
215  "Multiple solver system");
216  params.addParamNamesToGroup("skip_exception_check", "Advanced");
217 
218  return params;
219 }
220 
223  _num_grid_steps(cast_int<unsigned int>(getParam<unsigned int>("num_grids") - 1)),
224  _using_multi_sys_fp_iterations(getParam<bool>("multi_system_fixed_point")),
225  _multi_sys_fp_convergence(nullptr) // has not been created yet
226 {
227  if (_moose_line_searches.find(getParam<MooseEnum>("line_search").operator std::string()) !=
228  _moose_line_searches.end())
230 
231  auto set_solver_params = [this, &ex](const SolverSystem & sys)
232  {
233  const auto prefix = sys.prefix();
236 
237  // Set solver parameter prefix and system number
238  auto & solver_params = _problem.solverParams(sys.number());
239  solver_params._prefix = prefix;
240  solver_params._solver_sys_num = sys.number();
241  };
242 
243  // Extract and store PETSc related settings on FEProblemBase
244  for (const auto * const sys : _systems)
245  set_solver_params(*sys);
246 
247  // Set linear solve parameters in the equation system
248  // Nonlinear solve parameters are added in the DefaultNonlinearConvergence
249  EquationSystems & es = _problem.es();
250  es.parameters.set<Real>("linear solver tolerance") = getParam<Real>("l_tol");
251  es.parameters.set<Real>("linear solver absolute tolerance") = getParam<Real>("l_abs_tol");
252  es.parameters.set<unsigned int>("linear solver maximum iterations") =
253  getParam<unsigned int>("l_max_its");
254  es.parameters.set<bool>("reuse preconditioner") = getParam<bool>("reuse_preconditioner");
255  es.parameters.set<unsigned int>("reuse preconditioner maximum linear iterations") =
256  getParam<unsigned int>("reuse_preconditioner_max_linear_its");
257 
258  // Transfer to the Problem misc nonlinear solve optimization parameters
259  _problem.setSNESMFReuseBase(getParam<bool>("snesmf_reuse_base"),
260  _pars.isParamSetByUser("snesmf_reuse_base"));
261  _problem.skipExceptionCheck(getParam<bool>("skip_exception_check"));
262 
263  if (isParamValid("nonlinear_convergence"))
264  {
266  mooseError("The selected problem does not allow 'nonlinear_convergence' to be set.");
268  getParam<std::vector<ConvergenceName>>("nonlinear_convergence"));
269  }
270  else
272  if (isParamValid("linear_convergence"))
273  {
274  if (_problem.numLinearSystems() == 0)
275  paramError(
276  "linear_convergence",
277  "Setting 'linear_convergence' is currently only possible for solving linear systems");
279  getParam<std::vector<ConvergenceName>>("linear_convergence"));
280  }
281 
282  // Check whether the user has explicitly requested automatic scaling and is using a solve type
283  // without a matrix. If so, then we warn them
284  if ((_pars.isParamSetByUser("automatic_scaling") && getParam<bool>("automatic_scaling")) &&
285  std::all_of(_systems.begin(),
286  _systems.end(),
287  [this](const auto & solver_sys)
288  { return _problem.solverParams(solver_sys->number())._type == Moose::ST_JFNK; }))
289  {
290  paramWarning("automatic_scaling",
291  "Automatic scaling isn't implemented for the case where you do not have a "
292  "preconditioning matrix. No scaling will be applied");
293  _problem.automaticScaling(false);
294  }
295  else
296  // Check to see whether automatic_scaling has been specified anywhere, including at the
297  // application level. No matter what: if we don't have a matrix, we don't do scaling
299  isParamValid("automatic_scaling")
300  ? getParam<bool>("automatic_scaling")
301  : (getMooseApp().defaultAutomaticScaling() &&
302  std::any_of(_systems.begin(),
303  _systems.end(),
304  [this](const auto & solver_sys) {
305  return _problem.solverParams(solver_sys->number())._type !=
307  })));
308 
309  if (!_using_multi_sys_fp_iterations && isParamValid("multi_system_fixed_point_convergence"))
310  paramError("multi_system_fixed_point_convergence",
311  "Cannot set a convergence object for multi-system fixed point iterations if "
312  "'multi_system_fixed_point' is set to false");
313  if (_using_multi_sys_fp_iterations && !isParamValid("multi_system_fixed_point_convergence"))
314  paramError("multi_system_fixed_point_convergence",
315  "Must set a convergence object for multi-system fixed point iterations if using "
316  "multi-system fixed point iterations");
317 
318  // Set the same parameters to every nonlinear system by default
319  int i_nl_sys = -1;
320  for (const auto i_sys : index_range(_systems))
321  {
322  auto nl_ptr = dynamic_cast<NonlinearSystemBase *>(_systems[i_sys]);
323  // Linear systems have very different parameters at the moment
324  if (!nl_ptr)
325  continue;
326  auto & nl = *nl_ptr;
327  i_nl_sys++;
328 
329  nl.setPreSMOResidual(getParam<bool>("use_pre_SMO_residual"));
330 
331  const auto & all_splittings = getParam<std::vector<std::vector<std::string>>>("splitting");
332  if (all_splittings.size())
333  nl.setDecomposition(
334  getParamFromNonlinearSystemVectorParam<std::vector<std::string>>("splitting", i_nl_sys));
335  else
336  nl.setDecomposition({});
337 
338  const auto res_and_jac =
339  getParamFromNonlinearSystemVectorParam<bool>("residual_and_jacobian_together", i_nl_sys);
340  if (res_and_jac)
341  nl.residualAndJacobianTogether();
342 
343  // Automatic scaling parameters
344  nl.computeScalingOnce(
345  getParamFromNonlinearSystemVectorParam<bool>("compute_scaling_once", i_nl_sys));
346  nl.autoScalingParam(
347  getParamFromNonlinearSystemVectorParam<Real>("resid_vs_jac_scaling_param", i_nl_sys));
348  nl.offDiagonalsInAutoScaling(
349  getParamFromNonlinearSystemVectorParam<bool>("off_diagonals_in_auto_scaling", i_nl_sys));
350  if (isParamValid("scaling_group_variables"))
351  nl.scalingGroupVariables(
352  getParamFromNonlinearSystemVectorParam<std::vector<std::vector<std::string>>>(
353  "scaling_group_variables", i_nl_sys));
354  if (isParamValid("ignore_variables_for_autoscaling"))
355  {
356  // Before setting ignore_variables_for_autoscaling, check that they are not present in
357  // scaling_group_variables
358  if (isParamValid("scaling_group_variables"))
359  {
360  const auto & ignore_variables_for_autoscaling =
361  getParamFromNonlinearSystemVectorParam<std::vector<std::string>>(
362  "ignore_variables_for_autoscaling", i_nl_sys);
363  const auto & scaling_group_variables =
364  getParamFromNonlinearSystemVectorParam<std::vector<std::vector<std::string>>>(
365  "scaling_group_variables", i_nl_sys);
366  for (const auto & group : scaling_group_variables)
367  for (const auto & var_name : group)
368  if (std::find(ignore_variables_for_autoscaling.begin(),
369  ignore_variables_for_autoscaling.end(),
370  var_name) != ignore_variables_for_autoscaling.end())
371  paramError("ignore_variables_for_autoscaling",
372  "Variables cannot be in a scaling grouping and also be ignored");
373  }
374  nl.ignoreVariablesForAutoscaling(
375  getParamFromNonlinearSystemVectorParam<std::vector<std::string>>(
376  "ignore_variables_for_autoscaling", i_nl_sys));
377  }
378  }
379 
380  // Multi-grid options
382 }
383 
384 template <typename T>
385 T
387  unsigned int index) const
388 {
389  const auto & param_vec = getParam<std::vector<T>>(param_name);
390  if (index > _num_nl_systems)
391  paramError(param_name,
392  "Vector parameter is requested at index (" + std::to_string(index) +
393  ") which is larger than number of nonlinear systems (" +
394  std::to_string(_num_nl_systems) + ").");
395  if (param_vec.size() == 0)
396  paramError(
397  param_name,
398  "This parameter was passed to a routine which cannot handle empty vector parameters");
399  if (param_vec.size() != 1 && param_vec.size() != _num_nl_systems)
400  paramError(param_name,
401  "Vector parameter size (" + std::to_string(param_vec.size()) +
402  ") is different than the number of nonlinear systems (" +
403  std::to_string(_num_nl_systems) + ").");
404 
405  // User passed only one parameter, assume it applies to all nonlinear systems
406  if (param_vec.size() == 1)
407  return param_vec[0];
408  else
409  return param_vec[index];
410 }
411 
412 bool
414 {
415  // This should be late enough to retrieve the convergence object.
416  // TODO: Move this to a setup phase, which does not exist for SolveObjects
417  if (isParamValid("multi_system_fixed_point_convergence"))
419  &_problem.getConvergence(getParam<ConvergenceName>("multi_system_fixed_point_convergence"));
420 
421  // Outer loop for multi-grid convergence
422  bool converged = false;
423  unsigned int num_fp_multisys_iters = 0;
424  for (MooseIndex(_num_grid_steps) grid_step = 0; grid_step <= _num_grid_steps; ++grid_step)
425  {
426  // Multi-system fixed point loop
427  // Use a convergence object if provided, if not, use a reasonable default of every nested system
428  // being converged
429  num_fp_multisys_iters = 0;
430  converged = false;
431  while (!converged)
432  {
433  // Loop over each system
434  for (const auto sys : _systems)
435  {
436  const bool is_nonlinear = (dynamic_cast<NonlinearSystemBase *>(sys) != nullptr);
437 
438  // Call solve on the problem for that system
439  if (is_nonlinear)
440  _problem.solve(sys->number());
441  else
442  {
443  const auto linear_sys_number =
444  cast_int<unsigned int>(sys->number() - _problem.numNonlinearSystems());
445  _problem.solveLinearSystem(linear_sys_number, &_problem.getPetscOptions());
446 
447  // This is for postprocessing purposes in case none of the objects
448  // request the gradients.
449  // TODO: Somehow collect information if the postprocessors
450  // need gradients and if nothing needs this, just skip it
451  _problem.getLinearSystem(linear_sys_number).computeGradients();
452  }
453 
454  // Check convergence
455  const auto solve_name =
456  _systems.size() == 1 ? " Solve" : "System " + sys->name() + ": Solve";
457  if (_problem.shouldSolve())
458  {
459  if (_problem.converged(sys->number()))
460  _console << COLOR_GREEN << solve_name << " Converged!" << COLOR_DEFAULT << std::endl;
461  else
462  {
463  _console << COLOR_RED << solve_name << " Did NOT Converge!" << COLOR_DEFAULT
464  << std::endl;
465  return false;
466  }
467  }
468  else
469  _console << COLOR_GREEN << solve_name << " Skipped!" << COLOR_DEFAULT << std::endl;
470  }
471 
472  // Assess convergence of the multi-system fixed point iteration
474  converged = true;
475  else
476  {
477  converged = _multi_sys_fp_convergence->checkConvergence(num_fp_multisys_iters) ==
479  if (_multi_sys_fp_convergence->checkConvergence(num_fp_multisys_iters) ==
481  break;
482  }
483  num_fp_multisys_iters++;
484  }
485 
486  if (grid_step != _num_grid_steps)
488  }
489 
491  return (_multi_sys_fp_convergence->checkConvergence(num_fp_multisys_iters) ==
493  else
494  return converged;
495 }
static InputParameters validParams()
bool shouldSolve() const
FEProblemBase & _problem
Reference to FEProblem.
Definition: SolveObject.h:47
Moose::PetscSupport::PetscOptions & getPetscOptions()
Retrieve a writable reference the PETSc options (used by PetscSupport)
void computeGradients()
Compute the Green-Gauss gradients.
Definition: LinearSystem.C:142
std::string _prefix
Definition: SolverParams.h:35
virtual std::size_t numNonlinearSystems() const override
std::set< std::string > getPetscValidLineSearches()
Returns the valid petsc line search options as a set of strings.
Definition: PetscSupport.C:919
virtual MooseConvergenceStatus checkConvergence(unsigned int iter)=0
Returns convergence status.
static std::set< std::string > const _moose_line_searches
Moose provided line searches.
virtual bool onlyAllowDefaultNonlinearConvergence() const
Returns true if an error will result if the user supplies &#39;nonlinear_convergence&#39;.
virtual bool solve() override
Picard solve the FEProblem.
FEProblemSolve(Executioner &ex)
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
const unsigned int _num_grid_steps
The number of steps to perform in a grid sequencing algorithm.
virtual void solve(const unsigned int nl_sys_num)
MooseApp & getMooseApp() const
Get the MooseApp this class is associated with.
Definition: MooseBase.h:45
InputParameters emptyInputParameters()
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Nonlinear system to be solved.
void skipExceptionCheck(bool skip_exception_check)
Set a flag that indicates if we want to skip exception and stop solve.
virtual void addLineSearch(const InputParameters &)
add a MOOSE line search
void uniformRefine()
uniformly refine the problem mesh(es).
const std::string _type
The type of this class.
Definition: MooseBase.h:87
void numGridSteps(unsigned int num_grid_steps)
Set the number of steps in a grid sequences.
Convergence * _multi_sys_fp_convergence
Convergence object to assess the convergence of the multi-system fixed point iteration.
virtual Convergence & getConvergence(const std::string &name, const THREAD_ID tid=0) const
Gets a Convergence object.
Jacobian-Free Newton Krylov.
Definition: MooseTypes.h:845
virtual libMesh::EquationSystems & es() override
virtual bool converged(const unsigned int sys_num)
Eventually we want to convert this virtual over to taking a solver system number argument.
Definition: SubProblem.h:113
static InputParameters feProblemDefaultConvergenceParams()
const T & getParam(const std::string &name) const
Retrieve a parameter for the object.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:33
void setLinearConvergenceNames(const std::vector< ConvergenceName > &convergence_names)
Sets the linear convergence object name(s) if there is one.
Executioners are objects that do the actual work of solving your problem.
Definition: Executioner.h:30
void paramError(const std::string &param, Args... args) const
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
void setNonlinearConvergenceNames(const std::vector< ConvergenceName > &convergence_names)
Sets the nonlinear convergence object name(s) if there is one.
static InputParameters validParams()
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:64
void setNeedToAddDefaultNonlinearConvergence()
Sets _need_to_add_default_nonlinear_convergence to true.
LinearSystem & getLinearSystem(unsigned int sys_num)
Get non-constant reference to a linear system.
void setConvergedReasonFlags(FEProblemBase &fe_problem, const std::string &prefix)
Set flags that will instruct the user on the reason their simulation diverged from PETSc&#39;s perspectiv...
Definition: PetscSupport.C:714
bool isParamSetByUser(const std::string &name) const
Method returns true if the parameter was set by the user.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
InputParameters getPetscValidParams()
Returns the PETSc options that are common between Executioners and Preconditioners.
Definition: PetscSupport.C:925
T getParamFromNonlinearSystemVectorParam(const std::string &param_name, unsigned int index) const
Helper routine to get the nonlinear system parameter at the right index.
T & set(const std::string &)
unsigned int _num_nl_systems
Number of nonlinear systems.
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
const InputParameters & _pars
Parameters of this object, references the InputParameters stored in the InputParametersWarehouse.
SolverParams & solverParams(unsigned int solver_sys_num=0)
Get the solver parameters.
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an optional parameter and a documentation string to the InputParameters object...
const bool _using_multi_sys_fp_iterations
Whether we are using fixed point iterations for multi-system.
virtual std::size_t numLinearSystems() const override
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
void paramWarning(const std::string &param, Args... args) const
Emits a warning prefixed with the file and line number of the given param (from the input file) along...
void automaticScaling(bool automatic_scaling) override
Automatic scaling setter.
void storePetscOptions(FEProblemBase &fe_problem, const std::string &prefix, const ParallelParamObject &param_object)
Stores the PETSc options supplied from the parameter object on the problem.
Definition: PetscSupport.C:588
void setSNESMFReuseBase(bool reuse, bool set_by_user)
If or not to reuse the base vector for matrix-free calculation.
virtual void solveLinearSystem(const unsigned int linear_sys_num, const Moose::PetscSupport::PetscOptions *po=nullptr)
Build and solve a linear system.
void ErrorVector unsigned int
auto index_range(const T &sizable)
std::vector< SolverSystem * > _systems
Vector of pointers to the systems.
void setPreSMOResidual(bool use)
Set whether to evaluate the pre-SMO residual and use it in the subsequent relative convergence checks...
Tnew cast_int(Told oldvar)
static const std::set< std::string > & mooseLineSearches()
A solve object for use when wanting to solve multiple systems.
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)
This method takes a space delimited list of parameter names and adds them to the specified group name...