www.mooseframework.org
FEProblemSolve.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 "FEProblemSolve.h"
11 
12 #include "FEProblem.h"
13 #include "NonlinearSystemBase.h"
14 
15 std::set<std::string> const FEProblemSolve::_moose_line_searches = {"contact", "project"};
16 
17 const std::set<std::string> &
19 {
20  return _moose_line_searches;
21 }
22 
25 {
27 
28  params.addParam<std::vector<std::string>>("splitting",
29  {},
30  "Top-level splitting defining a "
31  "hierarchical decomposition into "
32  "subsystems to help the solver.");
33 
34  std::set<std::string> line_searches = mooseLineSearches();
35 
36  std::set<std::string> alias_line_searches = {"default", "none", "basic"};
37  line_searches.insert(alias_line_searches.begin(), alias_line_searches.end());
38  std::set<std::string> petsc_line_searches = Moose::PetscSupport::getPetscValidLineSearches();
39  line_searches.insert(petsc_line_searches.begin(), petsc_line_searches.end());
40  std::string line_search_string = Moose::stringify(line_searches, " ");
41  MooseEnum line_search(line_search_string, "default");
42  std::string addtl_doc_str(" (Note: none = basic)");
43  params.addParam<MooseEnum>(
44  "line_search", line_search, "Specifies the line search type" + addtl_doc_str);
45  MooseEnum line_search_package("petsc moose", "petsc");
46  params.addParam<MooseEnum>("line_search_package",
47  line_search_package,
48  "The solver package to use to conduct the line-search");
49 
50  params.addParam<unsigned>("contact_line_search_allowed_lambda_cuts",
51  2,
52  "The number of times lambda is allowed to be cut in half in the "
53  "contact line search. We recommend this number be roughly bounded by 0 "
54  "<= allowed_lambda_cuts <= 3");
55  params.addParam<Real>("contact_line_search_ltol",
56  "The linear relative tolerance to be used while the contact state is "
57  "changing between non-linear iterations. We recommend that this tolerance "
58  "be looser than the standard linear tolerance");
59 
61  params.addParam<Real>("l_tol", 1.0e-5, "Linear Relative Tolerance");
62  params.addParam<Real>("l_abs_tol", 1.0e-50, "Linear Absolute Tolerance");
63  params.addParam<unsigned int>("l_max_its", 10000, "Max Linear Iterations");
64  params.addParam<unsigned int>("nl_max_its", 50, "Max Nonlinear Iterations");
65  params.addParam<unsigned int>("nl_forced_its", 0, "The Number of Forced Nonlinear Iterations");
66  params.addParam<unsigned int>("nl_max_funcs", 10000, "Max Nonlinear solver function evaluations");
67  params.addParam<Real>("nl_abs_tol", 1.0e-50, "Nonlinear Absolute Tolerance");
68  params.addParam<Real>("nl_rel_tol", 1.0e-8, "Nonlinear Relative Tolerance");
69  params.addParam<Real>(
70  "nl_div_tol",
71  1.0e10,
72  "Nonlinear Relative Divergence Tolerance. A negative value disables this check.");
73  params.addParam<Real>(
74  "nl_abs_div_tol",
75  1.0e50,
76  "Nonlinear Absolute Divergence Tolerance. A negative value disables this check.");
77  params.addParam<Real>("nl_abs_step_tol", 0., "Nonlinear Absolute step Tolerance");
78  params.addParam<Real>("nl_rel_step_tol", 0., "Nonlinear Relative step Tolerance");
79  params.addParam<unsigned int>(
80  "n_max_nonlinear_pingpong",
81  100,
82  "The maximum number of times the nonlinear residual can ping pong "
83  "before requesting halting the current evaluation and requesting timestep cut");
84  params.addParam<bool>(
85  "snesmf_reuse_base",
86  true,
87  "Specifies whether or not to reuse the base vector for matrix-free calculation");
88  params.addParam<bool>(
89  "skip_exception_check", false, "Specifies whether or not to skip exception check");
90  params.addParam<bool>(
91  "compute_initial_residual_before_preset_bcs",
92  false,
93  "Use the residual norm computed *before* preset BCs are imposed in relative "
94  "convergence check");
95  params.addParam<bool>("automatic_scaling", "Whether to use automatic scaling for the variables.");
96  params.addParam<bool>(
97  "compute_scaling_once",
98  true,
99  "Whether the scaling factors should only be computed once at the beginning of the simulation "
100  "through an extra Jacobian evaluation. If this is set to false, then the scaling factors "
101  "will be computed during an extra Jacobian evaluation at the beginning of every time step.");
102  params.addParam<bool>(
103  "off_diagonals_in_auto_scaling",
104  false,
105  "Whether to consider off-diagonals when determining automatic scaling factors.");
106  params.addRangeCheckedParam<Real>(
107  "resid_vs_jac_scaling_param",
108  0,
109  "0<=resid_vs_jac_scaling_param<=1",
110  "A parameter that indicates the weighting of the residual vs the Jacobian in determining "
111  "variable scaling parameters. A value of 1 indicates pure residual-based scaling. A value of "
112  "0 indicates pure Jacobian-based scaling");
113  params.addParam<std::vector<std::vector<std::string>>>(
114  "scaling_group_variables",
115  "Name of variables that are grouped together for determining scale factors. (Multiple "
116  "groups can be provided, separated by semicolon)");
117  params.addParam<std::vector<std::string>>(
118  "ignore_variables_for_autoscaling",
119  "List of variables that do not participate in autoscaling.");
120  params.addRangeCheckedParam<unsigned int>(
121  "num_grids",
122  1,
123  "num_grids>0",
124  "The number of grids to use for a grid sequencing algorithm. This includes the final grid, "
125  "so num_grids = 1 indicates just one solve in a time-step");
126  params.addParam<bool>("residual_and_jacobian_together",
127  false,
128  "Whether to compute the residual and Jacobian together.");
129 
130  params.addParam<bool>("reuse_preconditioner",
131  false,
132  "If true reuse the previously calculated "
133  "preconditioner for the linearized "
134  "system across multiple solves "
135  "spanning nonlinear iterations and time steps. "
136  "The preconditioner resets as controlled by "
137  "reuse_preconditioner_max_linear_its");
138  params.addParam<unsigned int>("reuse_preconditioner_max_linear_its",
139  25,
140  "Reuse the previously calculated "
141  "preconditioner for the linear system "
142  "until the number of linear iterations "
143  "exceeds this number");
144 
145  params.addParamNamesToGroup("l_tol l_abs_tol l_max_its reuse_preconditioner "
146  "reuse_preconditioner_max_linear_its",
147  "Linear Solver");
148  params.addParamNamesToGroup("solve_type nl_max_its nl_forced_its nl_max_funcs "
149  "nl_abs_tol nl_rel_tol nl_abs_step_tol nl_rel_step_tol "
150  "snesmf_reuse_base compute_initial_residual_before_preset_bcs "
151  "num_grids nl_div_tol nl_abs_div_tol residual_and_jacobian_together "
152  "n_max_nonlinear_pingpong",
153  "Nonlinear Solver");
154  params.addParamNamesToGroup(
155  "automatic_scaling compute_scaling_once off_diagonals_in_auto_scaling "
156  "scaling_group_variables resid_vs_jac_scaling_param ignore_variables_for_autoscaling",
157  "Solver variable scaling");
158  params.addParamNamesToGroup("line_search line_search_package contact_line_search_ltol "
159  "contact_line_search_allowed_lambda_cuts",
160  "Solver line search");
161  params.addParamNamesToGroup("skip_exception_check", "Advanced");
162 
163  return params;
164 }
165 
167  : NonlinearSolveObject(ex),
168  _splitting(getParam<std::vector<std::string>>("splitting")),
169  _num_grid_steps(getParam<unsigned int>("num_grids") - 1)
170 {
171  if (_moose_line_searches.find(getParam<MooseEnum>("line_search").operator std::string()) !=
172  _moose_line_searches.end())
174 
175  // Extract and store PETSc related settings on FEProblemBase
177 
178  EquationSystems & es = _problem.es();
179  es.parameters.set<Real>("linear solver tolerance") = getParam<Real>("l_tol");
180 
181  es.parameters.set<Real>("linear solver absolute tolerance") = getParam<Real>("l_abs_tol");
182 
183  es.parameters.set<unsigned int>("linear solver maximum iterations") =
184  getParam<unsigned int>("l_max_its");
185 
186  es.parameters.set<unsigned int>("nonlinear solver maximum iterations") =
187  getParam<unsigned int>("nl_max_its");
188 
189  es.parameters.set<unsigned int>("nonlinear solver maximum function evaluations") =
190  getParam<unsigned int>("nl_max_funcs");
191 
192  es.parameters.set<Real>("nonlinear solver absolute residual tolerance") =
193  getParam<Real>("nl_abs_tol");
194 
195  es.parameters.set<Real>("nonlinear solver relative residual tolerance") =
196  getParam<Real>("nl_rel_tol");
197 
198  es.parameters.set<Real>("nonlinear solver divergence tolerance") = getParam<Real>("nl_div_tol");
199 
200  es.parameters.set<Real>("nonlinear solver absolute step tolerance") =
201  getParam<Real>("nl_abs_step_tol");
202 
203  es.parameters.set<Real>("nonlinear solver relative step tolerance") =
204  getParam<Real>("nl_rel_step_tol");
205 
206  es.parameters.set<bool>("reuse preconditioner") = getParam<bool>("reuse_preconditioner");
207 
208  es.parameters.set<unsigned int>("reuse preconditioner maximum linear iterations") =
209  getParam<unsigned int>("reuse_preconditioner_max_linear_its");
210 
212  getParam<bool>("compute_initial_residual_before_preset_bcs");
213 
214  _problem.setSNESMFReuseBase(getParam<bool>("snesmf_reuse_base"),
215  _pars.isParamSetByUser("snesmf_reuse_base"));
216 
217  _problem.skipExceptionCheck(getParam<bool>("skip_exception_check"));
218 
219  _problem.setMaxNLPingPong(getParam<unsigned int>("n_max_nonlinear_pingpong"));
220 
221  _problem.setNonlinearForcedIterations(getParam<unsigned int>("nl_forced_its"));
222 
223  _problem.setNonlinearAbsoluteDivergenceTolerance(getParam<Real>("nl_abs_div_tol"));
224 
226 
227  if (getParam<bool>("residual_and_jacobian_together"))
229 
230  // Check whether the user has explicitly requested automatic scaling and is using a solve type
231  // without a matrix. If so, then we warn them
232  if ((_pars.isParamSetByUser("automatic_scaling") && getParam<bool>("automatic_scaling")) &&
234  {
235  paramWarning("automatic_scaling",
236  "Automatic scaling isn't implemented for the case where you do not have a "
237  "preconditioning matrix. No scaling will be applied");
238  _problem.automaticScaling(false);
239  }
240  else
241  // Check to see whether automatic_scaling has been specified anywhere, including at the
242  // application level. No matter what: if we don't have a matrix, we don't do scaling
243  _problem.automaticScaling((isParamValid("automatic_scaling")
244  ? getParam<bool>("automatic_scaling")
245  : getMooseApp().defaultAutomaticScaling()) &&
247 
248  _nl.computeScalingOnce(getParam<bool>("compute_scaling_once"));
249  _nl.autoScalingParam(getParam<Real>("resid_vs_jac_scaling_param"));
250  _nl.offDiagonalsInAutoScaling(getParam<bool>("off_diagonals_in_auto_scaling"));
251  if (isParamValid("scaling_group_variables"))
253  getParam<std::vector<std::vector<std::string>>>("scaling_group_variables"));
254  if (isParamValid("ignore_variables_for_autoscaling"))
255  {
256  // Before setting ignore_variables_for_autoscaling, check that they are not present in
257  // scaling_group_variables
258  if (isParamValid("scaling_group_variables"))
259  {
260  const auto & ignore_variables_for_autoscaling =
261  getParam<std::vector<std::string>>("ignore_variables_for_autoscaling");
262  const auto & scaling_group_variables =
263  getParam<std::vector<std::vector<std::string>>>("scaling_group_variables");
264  for (const auto & group : scaling_group_variables)
265  for (const auto & var_name : group)
266  if (std::find(ignore_variables_for_autoscaling.begin(),
267  ignore_variables_for_autoscaling.end(),
268  var_name) != ignore_variables_for_autoscaling.end())
269  paramError("ignore_variables_for_autoscaling",
270  "Variables cannot be in a scaling grouping and also be ignored");
271  }
273  getParam<std::vector<std::string>>("ignore_variables_for_autoscaling"));
274  }
275 
277 }
278 
279 bool
281 {
282  // This loop is for nonlinear multigrids (developed by Alex)
283  for (MooseIndex(_num_grid_steps) grid_step = 0; grid_step <= _num_grid_steps; ++grid_step)
284  {
286 
287  if (_problem.shouldSolve())
288  {
289  if (_problem.converged(_nl.number()))
290  _console << COLOR_GREEN << " Solve Converged!" << COLOR_DEFAULT << std::endl;
291  else
292  {
293  _console << COLOR_RED << " Solve Did NOT Converge!" << COLOR_DEFAULT << std::endl;
294  return false;
295  }
296  }
297  else
298  _console << COLOR_GREEN << " Solve Skipped!" << COLOR_DEFAULT << std::endl;
299 
300  if (grid_step != _num_grid_steps)
302  }
303  return _problem.converged(_nl.number());
304 }
bool shouldSolve() const
FEProblemBase & _problem
Reference to FEProblem.
Definition: SolveObject.h:42
void setMaxNLPingPong(const unsigned int n_max_nl_pingpong)
method setting the maximum number of allowable non linear residual pingpong
void storePetscOptions(FEProblemBase &fe_problem, const InputParameters &params)
Stores the PETSc options supplied from the InputParameters with MOOSE.
Definition: PetscSupport.C:534
SolverParams & solverParams()
Get the solver parameters.
virtual bool converged(const unsigned int nl_sys_num)
Eventually we want to convert this virtual over to taking a nonlinear system number argument...
Definition: SubProblem.h:101
bool computeScalingOnce() const
std::set< std::string > getPetscValidLineSearches()
Returns the valid petsc line search options as a set of strings.
Definition: PetscSupport.C:833
static std::set< std::string > const _moose_line_searches
Moose provided line searches.
virtual bool solve() override
Picard solve the FEProblem.
FEProblemSolve(Executioner &ex)
void setDecomposition(const std::vector< std::string > &decomposition)
If called with a single string, it is used as the name of a the top-level decomposition split...
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
std::vector< std::string > _splitting
Splitting.
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:44
InputParameters emptyInputParameters()
virtual EquationSystems & es() override
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
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).
void scalingGroupVariables(const std::vector< std::vector< std::string >> &scaling_group_variables)
void numGridSteps(unsigned int num_grid_steps)
Set the number of steps in a grid sequences.
Jacobian-Free Newton Krylov.
Definition: MooseTypes.h:758
void ignoreVariablesForAutoscaling(const std::vector< std::string > &ignore_variables_for_autoscaling)
A solve object for use with a nonlinear system solver.
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:31
Moose::SolveType _type
Definition: SolverParams.h:19
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 ...
unsigned int number() const
Gets the number of this system.
Definition: SystemBase.C:1125
static InputParameters validParams()
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:62
void setNonlinearForcedIterations(const unsigned int nl_forced_its)
method setting the minimum number of nonlinear iterations before performing divergence checks ...
bool offDiagonalsInAutoScaling() const
bool _compute_initial_residual_before_preset_bcs
bool isParamSetByUser(const std::string &name) const
Method returns true if the parameter was 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:839
void setNonlinearAbsoluteDivergenceTolerance(const Real nl_abs_div_tol)
method setting the absolute divergence tolerance
const InputParameters & _pars
Parameters of this object, references the InputParameters stored in the InputParametersWarehouse.
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...
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.
virtual void residualAndJacobianTogether()=0
Call this method if you want the residual and Jacobian to be computed simultaneously.
void setSNESMFReuseBase(bool reuse, bool set_by_user)
If or not to reuse the base vector for matrix-free calculation.
NonlinearSystemBase & _nl
Reference to nonlinear system base for faster access.
Definition: SolveObject.h:50
void autoScalingParam(Real resid_vs_jac_scaling_param)
Sets the param that indicates the weighting of the residual vs the Jacobian in determining variable s...
void ErrorVector unsigned int
static const std::set< std::string > & mooseLineSearches()