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 template <>
18 {
20 
21  params.addParam<std::vector<std::string>>("splitting",
22  "Top-level splitting defining a "
23  "hierarchical decomposition into "
24  "subsystems to help the solver.");
25 
26  std::set<std::string> line_searches = {"contact", "default", "none", "basic"};
27 #ifdef LIBMESH_HAVE_PETSC
28  std::set<std::string> petsc_line_searches = Moose::PetscSupport::getPetscValidLineSearches();
29  line_searches.insert(petsc_line_searches.begin(), petsc_line_searches.end());
30 #endif // LIBMESH_HAVE_PETSC
31  std::string line_search_string = Moose::stringify(line_searches, " ");
32  MooseEnum line_search(line_search_string, "default");
33  std::string addtl_doc_str(" (Note: none = basic)");
34  params.addParam<MooseEnum>(
35  "line_search", line_search, "Specifies the line search type" + addtl_doc_str);
36  MooseEnum line_search_package("petsc moose", "petsc");
37  params.addParam<MooseEnum>("line_search_package",
38  line_search_package,
39  "The solver package to use to conduct the line-search");
40 
41  params.addParam<unsigned>("contact_line_search_allowed_lambda_cuts",
42  2,
43  "The number of times lambda is allowed to be cut in half in the "
44  "contact line search. We recommend this number be roughly bounded by 0 "
45  "<= allowed_lambda_cuts <= 3");
46  params.addParam<Real>("contact_line_search_ltol",
47  "The linear relative tolerance to be used while the contact state is "
48  "changing between non-linear iterations. We recommend that this tolerance "
49  "be looser than the standard linear tolerance");
50 
51  // Default Solver Behavior
52 #ifdef LIBMESH_HAVE_PETSC
54 #endif // LIBMESH_HAVE_PETSC
55  params.addParam<Real>("l_tol", 1.0e-5, "Linear Tolerance");
56  params.addParam<Real>("l_abs_step_tol", -1, "Linear Absolute Step Tolerance");
57  params.addParam<unsigned int>("l_max_its", 10000, "Max Linear Iterations");
58  params.addParam<unsigned int>("nl_max_its", 50, "Max Nonlinear Iterations");
59  params.addParam<unsigned int>("nl_max_funcs", 10000, "Max Nonlinear solver function evaluations");
60  params.addParam<Real>("nl_abs_tol", 1.0e-50, "Nonlinear Absolute Tolerance");
61  params.addParam<Real>("nl_rel_tol", 1.0e-8, "Nonlinear Relative Tolerance");
62  params.addParam<Real>("nl_abs_step_tol", 1.0e-50, "Nonlinear Absolute step Tolerance");
63  params.addParam<Real>("nl_rel_step_tol", 1.0e-50, "Nonlinear Relative step Tolerance");
64  params.addParam<bool>(
65  "snesmf_reuse_base",
66  true,
67  "Specifies whether or not to reuse the base vector for matrix-free calculation");
68  params.addParam<bool>("compute_initial_residual_before_preset_bcs",
69  false,
70  "Use the residual norm computed *before* PresetBCs are imposed in relative "
71  "convergence check");
72 
73  params.addParamNamesToGroup("l_tol l_abs_step_tol l_max_its nl_max_its nl_max_funcs "
74  "nl_abs_tol nl_rel_tol nl_abs_step_tol nl_rel_step_tol "
75  "snesmf_reuse_base compute_initial_residual_before_preset_bcs",
76  "Solver");
77  return params;
78 }
79 
81  : SolveObject(ex), _splitting(getParam<std::vector<std::string>>("splitting"))
82 {
83  if (_pars.isParamSetByUser("line_search"))
85 
86 // Extract and store PETSc related settings on FEProblemBase
87 #ifdef LIBMESH_HAVE_PETSC
89 #endif // LIBMESH_HAVE_PETSC
90 
91  EquationSystems & es = _problem.es();
92  es.parameters.set<Real>("linear solver tolerance") = getParam<Real>("l_tol");
93 
94  es.parameters.set<Real>("linear solver absolute step tolerance") =
95  getParam<Real>("l_abs_step_tol");
96 
97  es.parameters.set<unsigned int>("linear solver maximum iterations") =
98  getParam<unsigned int>("l_max_its");
99 
100  es.parameters.set<unsigned int>("nonlinear solver maximum iterations") =
101  getParam<unsigned int>("nl_max_its");
102 
103  es.parameters.set<unsigned int>("nonlinear solver maximum function evaluations") =
104  getParam<unsigned int>("nl_max_funcs");
105 
106  es.parameters.set<Real>("nonlinear solver absolute residual tolerance") =
107  getParam<Real>("nl_abs_tol");
108 
109  es.parameters.set<Real>("nonlinear solver relative residual tolerance") =
110  getParam<Real>("nl_rel_tol");
111 
112  es.parameters.set<Real>("nonlinear solver absolute step tolerance") =
113  getParam<Real>("nl_abs_step_tol");
114 
115  es.parameters.set<Real>("nonlinear solver relative step tolerance") =
116  getParam<Real>("nl_rel_step_tol");
117 
119  getParam<bool>("compute_initial_residual_before_preset_bcs");
120 
121  _nl._l_abs_step_tol = getParam<Real>("l_abs_step_tol");
122 
123  _problem.setSNESMFReuseBase(getParam<bool>("snesmf_reuse_base"),
124  _pars.isParamSetByUser("snesmf_reuse_base"));
125 
127 }
128 
129 bool
131 {
132  _problem.solve();
133  return _problem.converged();
134 }
FEProblemBase & _problem
Reference to FEProblem.
Definition: SolveObject.h:41
void storePetscOptions(FEProblemBase &fe_problem, const InputParameters &params)
Stores the PETSc options supplied from the InputParameters with MOOSE.
Definition: PetscSupport.C:634
std::set< std::string > getPetscValidLineSearches()
Returns the valid petsc line search options as a set of strings.
Definition: PetscSupport.C:819
FEProblemSolve(Executioner *ex)
virtual bool solve() override
Picard solve the FEProblem.
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.
InputParameters emptyInputParameters()
virtual EquationSystems & es() override
virtual void addLineSearch(const InputParameters &)
add a MOOSE line search
virtual bool converged() override
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
Executioners are objects that do the actual work of solving your problem.
Definition: Executioner.h:32
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:60
const InputParameters & _pars
Parameters of this object, references the InputParameters stored in the InputParametersWarehouse.
Definition: MooseObject.h:174
bool _compute_initial_residual_before_preset_bcs
bool isParamSetByUser(const std::string &name) const
Method returns true if the parameter was by the user.
InputParameters getPetscValidParams()
Returns the PETSc options that are common between Executioners and Preconditioners.
Definition: PetscSupport.C:829
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...
void setSNESMFReuseBase(bool reuse, bool set_by_user)
If or not to reuse the base vector for matrix-free calculation.
virtual void solve() override
NonlinearSystemBase & _nl
Reference to nonlinear system base for faster access.
Definition: SolveObject.h:49
InputParameters validParams< FEProblemSolve >()