www.mooseframework.org
Public Member Functions | Protected Attributes | List of all members
ReferenceResidualProblem Class Reference

FEProblemBase derived class to enable convergence checking relative to a user-specified postprocessor. More...

#include <ReferenceResidualProblem.h>

Inheritance diagram for ReferenceResidualProblem:
[legend]

Public Member Functions

 ReferenceResidualProblem (const InputParameters &params)
 
virtual ~ReferenceResidualProblem ()
 
virtual void initialSetup () override
 
virtual void timestepSetup () override
 
void updateReferenceResidual ()
 
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 PetscBool force_iteration, const Real ref_resid, const Real div_threshold) override
 
bool checkConvergenceIndividVars (const Real fnorm, const Real abstol, const Real rtol, const Real ref_resid)
 
void addSolutionVariables (std::set< std::string > &sol_vars)
 Add solution variables to ReferenceResidualProblem. More...
 
void addReferenceResidualVariables (std::set< std::string > &ref_vars)
 Add reference residual variables to ReferenceResidualProblem. More...
 
void addGroupVariables (std::set< std::string > &group_vars)
 Add a set of variables that need to be grouped together. More...
 

Protected Attributes

std::vector< unsigned int > _variable_group_num_index
 Group number index for each variable. More...
 
std::vector< Real > _scaling_factors
 Local storage for the scaling factors applied to each of the variables to apply to _ref_resid_vars. More...
 
std::vector< std::vector< std::string > > _group_variables
 Name of variables that are grouped together to check convergence. More...
 
bool _use_group_variables
 True if any variables are grouped. More...
 
std::vector< std::string > _soln_var_names
 
std::vector< std::string > _ref_resid_var_names
 
std::vector< std::string > _group_soln_var_names
 
std::vector< std::string > _group_ref_resid_var_names
 
std::vector< unsigned int > _soln_vars
 
std::vector< unsigned int > _ref_resid_vars
 
Real _accept_mult
 
int _accept_iters
 
std::vector< Real > _ref_resid
 
std::vector< Real > _resid
 
std::vector< Real > _group_ref_resid
 
std::vector< Real > _group_resid
 

Detailed Description

FEProblemBase derived class to enable convergence checking relative to a user-specified postprocessor.

Definition at line 24 of file ReferenceResidualProblem.h.

Constructor & Destructor Documentation

◆ ReferenceResidualProblem()

ReferenceResidualProblem::ReferenceResidualProblem ( const InputParameters &  params)

Definition at line 52 of file ReferenceResidualProblem.C.

53  : FEProblem(params), _use_group_variables(false)
54 {
55  if (params.isParamValid("solution_variables"))
56  _soln_var_names = params.get<std::vector<std::string>>("solution_variables");
57  if (params.isParamValid("reference_residual_variables"))
58  _ref_resid_var_names = params.get<std::vector<std::string>>("reference_residual_variables");
59  if (_soln_var_names.size() != _ref_resid_var_names.size())
60  mooseError("In ReferenceResidualProblem, size of solution_variables (",
61  _soln_var_names.size(),
62  ") != size of reference_residual_variables (",
63  _ref_resid_var_names.size(),
64  ")");
65 
66  if (params.isParamValid("group_variables"))
67  {
68  _group_variables = params.get<std::vector<std::vector<std::string>>>("group_variables");
69  _use_group_variables = true;
70  }
71 
72  _accept_mult = params.get<Real>("acceptable_multiplier");
73  _accept_iters = params.get<int>("acceptable_iterations");
74 }
bool _use_group_variables
True if any variables are grouped.
std::vector< std::vector< std::string > > _group_variables
Name of variables that are grouped together to check convergence.
std::vector< std::string > _soln_var_names
std::vector< std::string > _ref_resid_var_names

◆ ~ReferenceResidualProblem()

ReferenceResidualProblem::~ReferenceResidualProblem ( )
virtual

Definition at line 76 of file ReferenceResidualProblem.C.

76 {}

Member Function Documentation

◆ addGroupVariables()

void ReferenceResidualProblem::addGroupVariables ( std::set< std::string > &  group_vars)

Add a set of variables that need to be grouped together.

Parameters
group_varsA set of solution variables that need to be grouped.

Definition at line 93 of file ReferenceResidualProblem.C.

94 {
95  std::vector<std::string> group_vars_vector;
96 
97  for (const auto & var : group_vars)
98  group_vars_vector.push_back(var);
99 
100  _group_variables.push_back(group_vars_vector);
101 
102  _use_group_variables = true;
103 }
bool _use_group_variables
True if any variables are grouped.
std::vector< std::vector< std::string > > _group_variables
Name of variables that are grouped together to check convergence.

◆ addReferenceResidualVariables()

void ReferenceResidualProblem::addReferenceResidualVariables ( std::set< std::string > &  ref_vars)

Add reference residual variables to ReferenceResidualProblem.

Parameters
ref_varsA set of reference residual variables that need to be added to ReferenceResidualProblem.

Definition at line 86 of file ReferenceResidualProblem.C.

87 {
88  for (const auto & var : ref_vars)
89  _ref_resid_var_names.push_back(var);
90 }
std::vector< std::string > _ref_resid_var_names

◆ addSolutionVariables()

void ReferenceResidualProblem::addSolutionVariables ( std::set< std::string > &  sol_vars)

Add solution variables to ReferenceResidualProblem.

Parameters
sol_varsA set of solution variables that need to be added to ReferenceResidualProblem.

Definition at line 79 of file ReferenceResidualProblem.C.

80 {
81  for (const auto & var : sol_vars)
82  _soln_var_names.push_back(var);
83 }
std::vector< std::string > _soln_var_names

◆ checkConvergenceIndividVars()

bool ReferenceResidualProblem::checkConvergenceIndividVars ( const Real  fnorm,
const Real  abstol,
const Real  rtol,
const Real  ref_resid 
)

Definition at line 437 of file ReferenceResidualProblem.C.

Referenced by checkNonlinearConvergence().

441 {
442  bool convergedRelative = true;
443  if (_group_resid.size() > 0)
444  {
445  for (unsigned int i = 0; i < _group_resid.size(); ++i)
446  convergedRelative &=
447  ((_group_resid[i] < _group_ref_resid[i] * rtol) || (_group_resid[i] < abstol));
448  }
449 
450  else if (fnorm > ref_resid * rtol)
451  convergedRelative = false;
452 
453  return convergedRelative;
454 }
std::vector< Real > _group_ref_resid

◆ checkNonlinearConvergence()

MooseNonlinearConvergenceReason ReferenceResidualProblem::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 PetscBool  force_iteration,
const Real  ref_resid,
const Real  div_threshold 
)
overridevirtual

Reimplemented in AugmentedLagrangianContactProblem.

Definition at line 334 of file ReferenceResidualProblem.C.

Referenced by AugmentedLagrangianContactProblem::checkNonlinearConvergence().

347 {
349 
350  if (_group_soln_var_names.size() > 0)
351  {
352  _console << "Solution, reference convergence variable norms:" << std::endl;
353  unsigned int maxwsv = 0;
354  unsigned int maxwrv = 0;
355  for (unsigned int i = 0; i < _group_soln_var_names.size(); ++i)
356  {
357  if (_group_soln_var_names[i].size() > maxwsv)
358  maxwsv = _group_soln_var_names[i].size();
359  if (_group_ref_resid_var_names[i].size() > maxwrv)
360  maxwrv = _group_ref_resid_var_names[i].size();
361  }
362 
363  for (unsigned int i = 0; i < _group_soln_var_names.size(); ++i)
364  _console << std::setw(maxwsv + 2) << std::left << _group_soln_var_names[i] + ":"
365  << _group_resid[i] << " " << std::setw(maxwrv + 2)
366  << _group_ref_resid_var_names[i] + ":" << _group_ref_resid[i] << std::endl;
367  }
368 
369  NonlinearSystemBase & system = getNonlinearSystemBase();
370  MooseNonlinearConvergenceReason reason = MOOSE_NONLINEAR_ITERATING;
371  std::stringstream oss;
372 
373  if (fnorm != fnorm)
374  {
375  oss << "Failed to converge, function norm is NaN\n";
376  reason = MOOSE_DIVERGED_FNORM_NAN;
377  }
378  else if (fnorm < abstol && (it || !force_iteration))
379  {
380  oss << "Converged due to function norm " << fnorm << " < " << abstol << std::endl;
381  reason = MOOSE_CONVERGED_FNORM_ABS;
382  }
383  else if (nfuncs >= max_funcs)
384  {
385  oss << "Exceeded maximum number of function evaluations: " << nfuncs << " > " << max_funcs
386  << std::endl;
387  reason = MOOSE_DIVERGED_FUNCTION_COUNT;
388  }
389 
390  if (it && !reason)
391  {
392  if (checkConvergenceIndividVars(fnorm, abstol, rtol, ref_resid))
393  {
394  if (_group_resid.size() > 0)
395  oss << "Converged due to function norm "
396  << " < "
397  << " (relative tolerance) or (absolute tolerance) for all solution variables"
398  << std::endl;
399  else
400  oss << "Converged due to function norm " << fnorm << " < "
401  << " (relative tolerance)" << std::endl;
402  reason = MOOSE_CONVERGED_FNORM_RELATIVE;
403  }
404  else if (it >= _accept_iters &&
406  fnorm, abstol * _accept_mult, rtol * _accept_mult, ref_resid))
407  {
408  if (_group_resid.size() > 0)
409  oss << "Converged due to function norm "
410  << " < "
411  << " (acceptable relative tolerance) or (acceptable absolute tolerance) for all "
412  "solution variables"
413  << std::endl;
414  else
415  oss << "Converged due to function norm " << fnorm << " < "
416  << " (acceptable relative tolerance)" << std::endl;
417  _console << "ACCEPTABLE" << std::endl;
418  reason = MOOSE_CONVERGED_FNORM_RELATIVE;
419  }
420 
421  else if (snorm < stol * xnorm)
422  {
423  oss << "Converged due to small update length: " << snorm << " < " << stol << " * " << xnorm
424  << std::endl;
425  reason = MOOSE_CONVERGED_SNORM_RELATIVE;
426  }
427  }
428 
429  system._last_nl_rnorm = fnorm;
430  system._current_nl_its = static_cast<unsigned int>(it);
431 
432  msg = oss.str();
433  return reason;
434 }
std::vector< std::string > _group_soln_var_names
bool checkConvergenceIndividVars(const Real fnorm, const Real abstol, const Real rtol, const Real ref_resid)
std::vector< Real > _group_ref_resid
std::vector< std::string > _group_ref_resid_var_names

◆ initialSetup()

void ReferenceResidualProblem::initialSetup ( )
overridevirtual

Definition at line 106 of file ReferenceResidualProblem.C.

107 {
109 
110  unsigned int group_variable_num = 0;
112  {
113  for (unsigned int i = 0; i < _group_variables.size(); ++i)
114  {
115  group_variable_num += _group_variables[i].size();
116  if (_group_variables[i].size() == 1)
117  mooseError(" In the 'group_variables' parameter, variable ",
118  _group_variables[i][0],
119  " is not grouped with other variables.");
120  }
121 
122  unsigned int size = _soln_var_names.size() - group_variable_num + _group_variables.size();
123  _group_ref_resid.resize(size);
124  _group_resid.resize(size);
125  _group_soln_var_names.resize(size);
126  _group_ref_resid_var_names.resize(size);
127  }
128  else
129  {
130  _group_ref_resid.resize(_soln_var_names.size());
131  _group_resid.resize(_soln_var_names.size());
132  _group_soln_var_names.resize(_soln_var_names.size());
134  }
135 
136  std::set<std::string> check_duplicate;
138  {
139  for (unsigned int i = 0; i < _group_variables.size(); ++i)
140  for (unsigned int j = 0; j < _group_variables[i].size(); ++j)
141  check_duplicate.insert(_group_variables[i][j]);
142 
143  if (check_duplicate.size() != group_variable_num)
144  mooseError(
145  "A variable cannot be included in multiple groups in the 'group_variables' parameter.");
146  }
147 
148  NonlinearSystemBase & nonlinear_sys = getNonlinearSystemBase();
149  AuxiliarySystem & aux_sys = getAuxiliarySystem();
150  System & s = nonlinear_sys.system();
151  TransientExplicitSystem & as = aux_sys.sys();
152 
153  if (_soln_var_names.size() > 0 && _soln_var_names.size() != s.n_vars())
154  mooseError("In ReferenceResidualProblem, size of solution_variables (",
155  _soln_var_names.size(),
156  ") != number of variables in system (",
157  s.n_vars(),
158  ")");
159 
160  _soln_vars.clear();
161  for (unsigned int i = 0; i < _soln_var_names.size(); ++i)
162  {
163  bool foundMatch = false;
164  for (unsigned int var_num = 0; var_num < s.n_vars(); var_num++)
165  {
166  if (_soln_var_names[i] == s.variable_name(var_num))
167  {
168  _soln_vars.push_back(var_num);
169  _resid.push_back(0.0);
170  foundMatch = true;
171  break;
172  }
173  }
174  if (!foundMatch)
175  mooseError("Could not find solution variable '", _soln_var_names[i], "' in system");
176  }
177 
178  _ref_resid_vars.clear();
179  for (unsigned int i = 0; i < _ref_resid_var_names.size(); ++i)
180  {
181  bool foundMatch = false;
182  for (unsigned int var_num = 0; var_num < as.n_vars(); var_num++)
183  {
184  if (_ref_resid_var_names[i] == as.variable_name(var_num))
185  {
186  _ref_resid_vars.push_back(var_num);
187  _ref_resid.push_back(0.0);
188  foundMatch = true;
189  break;
190  }
191  }
192  if (!foundMatch)
193  mooseError("Could not find variable '", _ref_resid_var_names[i], "' in auxiliary system");
194  }
195 
196  unsigned int ungroup_index = 0;
198  ungroup_index = _group_variables.size();
199 
200  for (unsigned int i = 0; i < _soln_vars.size(); ++i)
201  {
202  bool find_group = false;
204  {
205  for (unsigned int j = 0; j < _group_variables.size(); ++j)
206  {
207  if (std::find(_group_variables[j].begin(),
208  _group_variables[j].end(),
209  s.variable_name(_soln_vars[i])) != _group_variables[j].end())
210  {
212  find_group = true;
213  break;
214  }
215  }
216  if (!find_group)
217  {
218  _variable_group_num_index[i] = ungroup_index;
219  ungroup_index++;
220  }
221  }
222  else
223  {
225  }
226  }
227 
229  {
230  for (unsigned int i = 0; i < _group_variables.size(); ++i)
231  {
232  bool same_variable = true;
233  if (_group_variables[i].size() > 1)
234  {
235  for (unsigned int j = 0; j < _group_variables[i].size(); ++j)
236  {
237  for (unsigned int var_num = 0; var_num < s.n_vars(); var_num++)
238  {
239  if (_group_variables[i][j] == s.variable_name(var_num))
240  {
241  if (j == 0)
242  same_variable = nonlinear_sys.isScalarVariable(_soln_vars[var_num]);
243  else
244  same_variable =
245  (same_variable == nonlinear_sys.isScalarVariable(_soln_vars[var_num]) ? true
246  : false);
247  break;
248  }
249  }
250  }
251  }
252  if (!same_variable)
253  mooseWarning("In the 'group_variables' parameter, standard variables and scalar variables "
254  "are grouped together in group ",
255  i);
256  }
257  }
258 
259  for (unsigned int i = 0; i < _soln_var_names.size(); ++i)
260  {
262  {
266  }
267 
269  {
273  }
274  }
275 
276  const unsigned int size_solnVars = _soln_vars.size();
277  _scaling_factors.resize(size_solnVars);
278  for (unsigned int i = 0; i < size_solnVars; ++i)
279  if (nonlinear_sys.isScalarVariable(_soln_vars[i]))
280  _scaling_factors[i] = nonlinear_sys.getScalarVariable(0, _soln_vars[i]).scalingFactor();
281  else
282  _scaling_factors[i] = nonlinear_sys.getVariable(/*tid*/ 0, _soln_vars[i]).scalingFactor();
283 
284  FEProblemBase::initialSetup();
285 }
std::vector< Real > _scaling_factors
Local storage for the scaling factors applied to each of the variables to apply to _ref_resid_vars...
bool _use_group_variables
True if any variables are grouped.
std::vector< std::string > _group_soln_var_names
std::vector< unsigned int > _soln_vars
std::vector< std::vector< std::string > > _group_variables
Name of variables that are grouped together to check convergence.
std::vector< std::string > _soln_var_names
std::vector< Real > _group_ref_resid
std::vector< unsigned int > _variable_group_num_index
Group number index for each variable.
std::vector< std::string > _ref_resid_var_names
std::vector< std::string > _group_ref_resid_var_names
std::vector< unsigned int > _ref_resid_vars

◆ timestepSetup()

void ReferenceResidualProblem::timestepSetup ( )
overridevirtual

Reimplemented in AugmentedLagrangianContactProblem.

Definition at line 288 of file ReferenceResidualProblem.C.

Referenced by AugmentedLagrangianContactProblem::timestepSetup().

289 {
290  for (unsigned int i = 0; i < _ref_resid.size(); ++i)
291  {
292  _ref_resid[i] = 0.0;
293  _resid[i] = 0.0;
294  }
295  FEProblemBase::timestepSetup();
296 }

◆ updateReferenceResidual()

void ReferenceResidualProblem::updateReferenceResidual ( )

Definition at line 299 of file ReferenceResidualProblem.C.

Referenced by checkNonlinearConvergence().

300 {
301  NonlinearSystemBase & nonlinear_sys = getNonlinearSystemBase();
302  AuxiliarySystem & aux_sys = getAuxiliarySystem();
303  System & s = nonlinear_sys.system();
304  TransientExplicitSystem & as = aux_sys.sys();
305 
306  for (unsigned int i = 0; i < _group_resid.size(); ++i)
307  {
308  _group_resid[i] = 0.0;
309  _group_ref_resid[i] = 0.0;
310  }
311 
312  for (unsigned int i = 0; i < _soln_vars.size(); ++i)
313  {
314  _resid[i] = s.calculate_norm(nonlinear_sys.RHS(), _soln_vars[i], DISCRETE_L2);
315  _group_resid[_variable_group_num_index[i]] += Utility::pow<2>(_resid[i]);
316  }
317 
318  for (unsigned int i = 0; i < _ref_resid_vars.size(); ++i)
319  {
320  const Real refResidual =
321  as.calculate_norm(*as.current_local_solution, _ref_resid_vars[i], DISCRETE_L2);
322  _ref_resid[i] = refResidual * _scaling_factors[i];
323  _group_ref_resid[_variable_group_num_index[i]] += Utility::pow<2>(_ref_resid[i]);
324  }
325 
326  for (unsigned int i = 0; i < _group_resid.size(); ++i)
327  {
328  _group_resid[i] = sqrt(_group_resid[i]);
329  _group_ref_resid[i] = sqrt(_group_ref_resid[i]);
330  }
331 }
std::vector< Real > _scaling_factors
Local storage for the scaling factors applied to each of the variables to apply to _ref_resid_vars...
std::vector< unsigned int > _soln_vars
std::vector< Real > _group_ref_resid
std::vector< unsigned int > _variable_group_num_index
Group number index for each variable.
std::vector< unsigned int > _ref_resid_vars

Member Data Documentation

◆ _accept_iters

int ReferenceResidualProblem::_accept_iters
protected

◆ _accept_mult

Real ReferenceResidualProblem::_accept_mult
protected

"Acceptable" absolute and relative tolerance multiplier and acceptable number of iterations. Used when checking the convergence of individual variables.

Definition at line 97 of file ReferenceResidualProblem.h.

Referenced by checkNonlinearConvergence(), and ReferenceResidualProblem().

◆ _group_ref_resid

std::vector<Real> ReferenceResidualProblem::_group_ref_resid
protected

Local storage for discrete L2 residual norms of the grouped variables listed in _group_ref_resid_var_names.

Definition at line 109 of file ReferenceResidualProblem.h.

Referenced by checkConvergenceIndividVars(), checkNonlinearConvergence(), initialSetup(), and updateReferenceResidual().

◆ _group_ref_resid_var_names

std::vector<std::string> ReferenceResidualProblem::_group_ref_resid_var_names
protected

Definition at line 84 of file ReferenceResidualProblem.h.

Referenced by checkNonlinearConvergence(), and initialSetup().

◆ _group_resid

std::vector<Real> ReferenceResidualProblem::_group_resid
protected

◆ _group_soln_var_names

std::vector<std::string> ReferenceResidualProblem::_group_soln_var_names
protected

List of grouped solution variable names whose reference residuals will be stored, and the residual variable names that will store them.

Definition at line 83 of file ReferenceResidualProblem.h.

Referenced by checkNonlinearConvergence(), and initialSetup().

◆ _group_variables

std::vector<std::vector<std::string> > ReferenceResidualProblem::_group_variables
protected

Name of variables that are grouped together to check convergence.

Definition at line 120 of file ReferenceResidualProblem.h.

Referenced by addGroupVariables(), initialSetup(), and ReferenceResidualProblem().

◆ _ref_resid

std::vector<Real> ReferenceResidualProblem::_ref_resid
protected

Local storage for discrete L2 residual norms of the variables listed in _ref_resid_var_names.

Definition at line 103 of file ReferenceResidualProblem.h.

Referenced by initialSetup(), timestepSetup(), and updateReferenceResidual().

◆ _ref_resid_var_names

std::vector<std::string> ReferenceResidualProblem::_ref_resid_var_names
protected

◆ _ref_resid_vars

std::vector<unsigned int> ReferenceResidualProblem::_ref_resid_vars
protected

Definition at line 90 of file ReferenceResidualProblem.h.

Referenced by initialSetup(), and updateReferenceResidual().

◆ _resid

std::vector<Real> ReferenceResidualProblem::_resid
protected

◆ _scaling_factors

std::vector<Real> ReferenceResidualProblem::_scaling_factors
protected

Local storage for the scaling factors applied to each of the variables to apply to _ref_resid_vars.

Definition at line 117 of file ReferenceResidualProblem.h.

Referenced by initialSetup(), and updateReferenceResidual().

◆ _soln_var_names

std::vector<std::string> ReferenceResidualProblem::_soln_var_names
protected

List of solution variable names whose reference residuals will be stored, and the residual variable names that will store them.

Definition at line 76 of file ReferenceResidualProblem.h.

Referenced by addSolutionVariables(), initialSetup(), and ReferenceResidualProblem().

◆ _soln_vars

std::vector<unsigned int> ReferenceResidualProblem::_soln_vars
protected

Variable numbers assoicated with the names in _soln_var_names and _ref_resid_var_names.

Definition at line 89 of file ReferenceResidualProblem.h.

Referenced by initialSetup(), and updateReferenceResidual().

◆ _use_group_variables

bool ReferenceResidualProblem::_use_group_variables
protected

True if any variables are grouped.

Definition at line 123 of file ReferenceResidualProblem.h.

Referenced by addGroupVariables(), initialSetup(), and ReferenceResidualProblem().

◆ _variable_group_num_index

std::vector<unsigned int> ReferenceResidualProblem::_variable_group_num_index
protected

Group number index for each variable.

Definition at line 114 of file ReferenceResidualProblem.h.

Referenced by initialSetup(), and updateReferenceResidual().


The documentation for this class was generated from the following files: