www.mooseframework.org
ReferenceResidualProblem.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
8 // MOOSE includes
10 
11 #include "AuxiliarySystem.h"
12 #include "MooseApp.h"
13 #include "MooseMesh.h"
14 #include "NonlinearSystem.h"
15 
16 template <>
17 InputParameters
19 {
20  InputParameters params = validParams<FEProblem>();
21  params.addClassDescription("Problem that checks for convergence relative to "
22  "a user-supplied reference quantity rather than "
23  "the initial residual");
24  params.addParam<std::vector<std::string>>(
25  "solution_variables", "Set of solution variables to be checked for relative convergence");
26  params.addParam<std::vector<std::string>>(
27  "reference_residual_variables",
28  "Set of variables that provide reference residuals for relative convergence check");
29  params.addParam<Real>("acceptable_multiplier",
30  1.0,
31  "Multiplier applied to relative tolerance for acceptable limit");
32  params.addParam<int>("acceptable_iterations",
33  0,
34  "Iterations after which convergence to acceptable limits is accepted");
35  return params;
36 }
37 
39  : FEProblem(params)
40 {
41  if (params.isParamValid("solution_variables"))
42  _solnVarNames = params.get<std::vector<std::string>>("solution_variables");
43  if (params.isParamValid("reference_residual_variables"))
44  _refResidVarNames = params.get<std::vector<std::string>>("reference_residual_variables");
45  if (_solnVarNames.size() != _refResidVarNames.size())
46  mooseError("In ReferenceResidualProblem, size of solution_variables (",
47  _solnVarNames.size(),
48  ") != size of reference_residual_variables (",
49  _refResidVarNames.size(),
50  ")");
51 
52  _accept_mult = params.get<Real>("acceptable_multiplier");
53  _accept_iters = params.get<int>("acceptable_iterations");
54 }
55 
57 
58 void
60 {
61  NonlinearSystemBase & nonlinear_sys = getNonlinearSystemBase();
62  AuxiliarySystem & aux_sys = getAuxiliarySystem();
63  System & s = nonlinear_sys.system();
64  TransientExplicitSystem & as = aux_sys.sys();
65 
66  if (_solnVarNames.size() > 0 && _solnVarNames.size() != s.n_vars())
67  mooseError("In ReferenceResidualProblem, size of solution_variables (",
68  _solnVarNames.size(),
69  ") != number of variables in system (",
70  s.n_vars(),
71  ")");
72 
73  _solnVars.clear();
74  for (unsigned int i = 0; i < _solnVarNames.size(); ++i)
75  {
76  bool foundMatch = false;
77  for (unsigned int var_num = 0; var_num < s.n_vars(); var_num++)
78  {
79  if (_solnVarNames[i] == s.variable_name(var_num))
80  {
81  _solnVars.push_back(var_num);
82  _resid.push_back(0.0);
83  foundMatch = true;
84  break;
85  }
86  }
87  if (!foundMatch)
88  mooseError("Could not find solution variable '", _solnVarNames[i], "' in system");
89  }
90 
91  _refResidVars.clear();
92  for (unsigned int i = 0; i < _refResidVarNames.size(); ++i)
93  {
94  bool foundMatch = false;
95  for (unsigned int var_num = 0; var_num < as.n_vars(); var_num++)
96  {
97  if (_refResidVarNames[i] == as.variable_name(var_num))
98  {
99  _refResidVars.push_back(var_num);
100  _refResid.push_back(0.0);
101  foundMatch = true;
102  break;
103  }
104  }
105  if (!foundMatch)
106  mooseError("Could not find variable '", _refResidVarNames[i], "' in auxiliary system");
107  }
108 
109  FEProblemBase::initialSetup();
110 }
111 
112 void
114 {
115  for (unsigned int i = 0; i < _refResid.size(); ++i)
116  {
117  _refResid[i] = 0.0;
118  _resid[i] = 0.0;
119  }
120  FEProblemBase::timestepSetup();
121 }
122 
123 void
125 {
126  NonlinearSystemBase & nonlinear_sys = getNonlinearSystemBase();
127  AuxiliarySystem & aux_sys = getAuxiliarySystem();
128  System & s = nonlinear_sys.system();
129  TransientExplicitSystem & as = aux_sys.sys();
130 
131  for (unsigned int i = 0; i < _solnVars.size(); ++i)
132  _resid[i] = s.calculate_norm(nonlinear_sys.RHS(), _solnVars[i], DISCRETE_L2);
133 
134  for (unsigned int i = 0; i < _refResidVars.size(); ++i)
135  _refResid[i] = as.calculate_norm(*as.current_local_solution, _refResidVars[i], DISCRETE_L2);
136 }
137 
138 MooseNonlinearConvergenceReason
140  const PetscInt it,
141  const Real xnorm,
142  const Real snorm,
143  const Real fnorm,
144  const Real rtol,
145  const Real stol,
146  const Real abstol,
147  const PetscInt nfuncs,
148  const PetscInt max_funcs,
149  const Real ref_resid,
150  const Real /*div_threshold*/)
151 {
153 
154  if (_solnVars.size() > 0)
155  {
156  _console << "Solution, reference convergence variable norms:" << std::endl;
157  unsigned int maxwsv = 0;
158  unsigned int maxwrv = 0;
159  for (unsigned int i = 0; i < _solnVars.size(); ++i)
160  {
161  if (_solnVarNames[i].size() > maxwsv)
162  maxwsv = _solnVarNames[i].size();
163  if (_refResidVarNames[i].size() > maxwrv)
164  maxwrv = _refResidVarNames[i].size();
165  }
166 
167  for (unsigned int i = 0; i < _solnVars.size(); ++i)
168  _console << std::setw(maxwsv + 2) << std::left << _solnVarNames[i] + ":" << _resid[i] << " "
169  << std::setw(maxwrv + 2) << _refResidVarNames[i] + ":" << _refResid[i] << std::endl;
170  }
171 
172  NonlinearSystemBase & system = getNonlinearSystemBase();
173  MooseNonlinearConvergenceReason reason = MOOSE_NONLINEAR_ITERATING;
174  std::stringstream oss;
175 
176  if (fnorm != fnorm)
177  {
178  oss << "Failed to converge, function norm is NaN\n";
179  reason = MOOSE_DIVERGED_FNORM_NAN;
180  }
181  else if (fnorm < abstol)
182  {
183  oss << "Converged due to function norm " << fnorm << " < " << abstol << std::endl;
184  reason = MOOSE_CONVERGED_FNORM_ABS;
185  }
186  else if (nfuncs >= max_funcs)
187  {
188  oss << "Exceeded maximum number of function evaluations: " << nfuncs << " > " << max_funcs
189  << std::endl;
190  reason = MOOSE_DIVERGED_FUNCTION_COUNT;
191  }
192 
193  if (it && !reason)
194  {
195  if (checkConvergenceIndividVars(fnorm, abstol, rtol, ref_resid))
196  {
197  if (_resid.size() > 0)
198  oss << "Converged due to function norm "
199  << " < "
200  << " (relative tolerance) or (absolute tolerance) for all solution variables"
201  << std::endl;
202  else
203  oss << "Converged due to function norm " << fnorm << " < "
204  << " (relative tolerance)" << std::endl;
205  reason = MOOSE_CONVERGED_FNORM_RELATIVE;
206  }
207  else if (it >= _accept_iters &&
209  fnorm, abstol * _accept_mult, rtol * _accept_mult, ref_resid))
210  {
211  if (_resid.size() > 0)
212  oss << "Converged due to function norm "
213  << " < "
214  << " (acceptable relative tolerance) or (acceptable absolute tolerance) for all "
215  "solution variables"
216  << std::endl;
217  else
218  oss << "Converged due to function norm " << fnorm << " < "
219  << " (acceptable relative tolerance)" << std::endl;
220  _console << "ACCEPTABLE" << std::endl;
221  reason = MOOSE_CONVERGED_FNORM_RELATIVE;
222  }
223 
224  else if (snorm < stol * xnorm)
225  {
226  oss << "Converged due to small update length: " << snorm << " < " << stol << " * " << xnorm
227  << std::endl;
228  reason = MOOSE_CONVERGED_SNORM_RELATIVE;
229  }
230  }
231 
232  system._last_nl_rnorm = fnorm;
233  system._current_nl_its = static_cast<unsigned int>(it);
234 
235  msg = oss.str();
236  return reason;
237 }
238 
239 bool
241  const Real abstol,
242  const Real rtol,
243  const Real ref_resid)
244 {
245  bool convergedRelative = true;
246  if (_resid.size() > 0)
247  {
248  for (unsigned int i = 0; i < _resid.size(); ++i)
249  convergedRelative &= ((_resid[i] < _refResid[i] * rtol) || (_resid[i] < abstol));
250  }
251 
252  else if (fnorm > ref_resid * rtol)
253  convergedRelative = false;
254 
255  return convergedRelative;
256 }
std::vector< unsigned int > _solnVars
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 ref_resid, const Real div_threshold)
std::vector< std::string > _solnVarNames
bool checkConvergenceIndividVars(const Real fnorm, const Real abstol, const Real rtol, const Real ref_resid)
std::vector< unsigned int > _refResidVars
InputParameters validParams< ReferenceResidualProblem >()
std::vector< std::string > _refResidVarNames
ReferenceResidualProblem(const InputParameters &params)