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

Class to manage nested solution for augmented Lagrange contact. More...

#include <AugmentedLagrangianContactProblem.h>

Inheritance diagram for AugmentedLagrangianContactProblem:
[legend]

Public Member Functions

 AugmentedLagrangianContactProblem (const InputParameters &params)
 
virtual ~AugmentedLagrangianContactProblem ()
 
virtual void timestepSetup () override
 
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
 
virtual void initialSetup () override
 
void updateReferenceResidual ()
 
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
 

Private Attributes

int _num_lagmul_iterations
 
int _max_lagmul_iters
 

Detailed Description

Class to manage nested solution for augmented Lagrange contact.

The AugmentedLagrangianContactProblem manages the nested solution procedure, repeating the solution until convergence has been achieved, checking for convergence, and updating the Lagrangian multipliers.

Definition at line 29 of file AugmentedLagrangianContactProblem.h.

Constructor & Destructor Documentation

◆ AugmentedLagrangianContactProblem()

AugmentedLagrangianContactProblem::AugmentedLagrangianContactProblem ( const InputParameters &  params)

Definition at line 41 of file AugmentedLagrangianContactProblem.C.

42  : ReferenceResidualProblem(params),
44  _max_lagmul_iters(getParam<int>("maximum_lagrangian_update_iterations"))
45 {
46 }
ReferenceResidualProblem(const InputParameters &params)

◆ ~AugmentedLagrangianContactProblem()

virtual AugmentedLagrangianContactProblem::~AugmentedLagrangianContactProblem ( )
inlinevirtual

Definition at line 33 of file AugmentedLagrangianContactProblem.h.

33 {}

Member Function Documentation

◆ addGroupVariables()

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

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)
inherited

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)
inherited

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 
)
inherited

Definition at line 437 of file ReferenceResidualProblem.C.

Referenced by ReferenceResidualProblem::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 AugmentedLagrangianContactProblem::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 from ReferenceResidualProblem.

Definition at line 56 of file AugmentedLagrangianContactProblem.C.

69 {
70 
71  Real my_max_funcs = std::numeric_limits<int>::max();
72  Real my_div_threshold = std::numeric_limits<Real>::max();
73 
74  MooseNonlinearConvergenceReason reason =
76  it,
77  xnorm,
78  snorm,
79  fnorm,
80  rtol,
81  stol,
82  abstol,
83  nfuncs,
84  my_max_funcs,
85  force_iteration,
86  ref_resid,
87  my_div_threshold);
88 
89  _console << "Augmented Lagrangian contact iteration " << _num_lagmul_iterations << "\n";
90 
91  bool _augLM_repeat_step;
92 
93  if (reason > 0)
94  {
96  {
97  NonlinearSystemBase & nonlinear_sys = getNonlinearSystemBase();
98  nonlinear_sys.update();
99 
100  const ConstraintWarehouse & constraints = nonlinear_sys.getConstraintWarehouse();
101 
102  std::map<std::pair<unsigned int, unsigned int>, PenetrationLocator *> * penetration_locators =
103  NULL;
104 
105  bool displaced = false;
106  _augLM_repeat_step = false;
107  if (getDisplacedProblem() == NULL)
108  {
109  GeometricSearchData & geom_search_data = geomSearchData();
110  penetration_locators = &geom_search_data._penetration_locators;
111  }
112  else
113  {
114  GeometricSearchData & displaced_geom_search_data = getDisplacedProblem()->geomSearchData();
115  penetration_locators = &displaced_geom_search_data._penetration_locators;
116  displaced = true;
117  }
118 
119  for (const auto & it : *penetration_locators)
120  {
121  PenetrationLocator & pen_loc = *(it.second);
122 
123  BoundaryID slave_boundary = pen_loc._slave_boundary;
124 
125  if (constraints.hasActiveNodeFaceConstraints(slave_boundary, displaced))
126  {
127  const auto & ncs = constraints.getActiveNodeFaceConstraints(slave_boundary, displaced);
128 
129  for (const auto & nc : ncs)
130  {
131  if (std::dynamic_pointer_cast<MechanicalContactConstraint>(nc) == NULL)
132  mooseError("AugmentedLagrangianContactProblem: dynamic cast of "
133  "MechanicalContactConstraint object failed.");
134 
135  if (!(std::dynamic_pointer_cast<MechanicalContactConstraint>(nc))
136  ->AugmentedLagrangianContactConverged())
137  {
138  (std::dynamic_pointer_cast<MechanicalContactConstraint>(nc))
139  ->updateAugmentedLagrangianMultiplier(false);
140  _augLM_repeat_step = true;
141  break;
142  }
143  }
144  }
145  }
146 
147  if (_augLM_repeat_step)
148  {
149  // force it to keep iterating
150  reason = MOOSE_NONLINEAR_ITERATING;
151  _console << "Augmented Lagrangian Multiplier needs updating." << std::endl;
153  }
154  else
155  _console << "Augmented Lagrangian contact constraint enforcement is satisfied."
156  << std::endl;
157  }
158  else
159  {
160  // maxed out
161  _console << "Maximum Augmented Lagrangian contact iterations have been reached." << std::endl;
162  reason = MOOSE_DIVERGED_FUNCTION_COUNT;
163  }
164  }
165 
166  return reason;
167 }
A MechanicalContactConstraint forces the value of a variable to be the same on both sides of an inter...
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

◆ initialSetup()

void ReferenceResidualProblem::initialSetup ( )
overridevirtualinherited

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 AugmentedLagrangianContactProblem::timestepSetup ( )
overridevirtual

◆ updateReferenceResidual()

void ReferenceResidualProblem::updateReferenceResidual ( )
inherited

Definition at line 299 of file ReferenceResidualProblem.C.

Referenced by ReferenceResidualProblem::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
protectedinherited

◆ _accept_mult

Real ReferenceResidualProblem::_accept_mult
protectedinherited

"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 ReferenceResidualProblem::checkNonlinearConvergence(), and ReferenceResidualProblem::ReferenceResidualProblem().

◆ _group_ref_resid

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

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 ReferenceResidualProblem::checkConvergenceIndividVars(), ReferenceResidualProblem::checkNonlinearConvergence(), ReferenceResidualProblem::initialSetup(), and ReferenceResidualProblem::updateReferenceResidual().

◆ _group_ref_resid_var_names

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

◆ _group_resid

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

◆ _group_soln_var_names

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

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 ReferenceResidualProblem::checkNonlinearConvergence(), and ReferenceResidualProblem::initialSetup().

◆ _group_variables

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

Name of variables that are grouped together to check convergence.

Definition at line 120 of file ReferenceResidualProblem.h.

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

◆ _max_lagmul_iters

int AugmentedLagrangianContactProblem::_max_lagmul_iters
private

Definition at line 54 of file AugmentedLagrangianContactProblem.h.

Referenced by checkNonlinearConvergence().

◆ _num_lagmul_iterations

int AugmentedLagrangianContactProblem::_num_lagmul_iterations
private

Definition at line 53 of file AugmentedLagrangianContactProblem.h.

Referenced by checkNonlinearConvergence(), and timestepSetup().

◆ _ref_resid

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

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 ReferenceResidualProblem::initialSetup(), ReferenceResidualProblem::timestepSetup(), and ReferenceResidualProblem::updateReferenceResidual().

◆ _ref_resid_var_names

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

◆ _ref_resid_vars

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

◆ _resid

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

◆ _scaling_factors

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

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 ReferenceResidualProblem::initialSetup(), and ReferenceResidualProblem::updateReferenceResidual().

◆ _soln_var_names

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

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 ReferenceResidualProblem::addSolutionVariables(), ReferenceResidualProblem::initialSetup(), and ReferenceResidualProblem::ReferenceResidualProblem().

◆ _soln_vars

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

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 ReferenceResidualProblem::initialSetup(), and ReferenceResidualProblem::updateReferenceResidual().

◆ _use_group_variables

bool ReferenceResidualProblem::_use_group_variables
protectedinherited

◆ _variable_group_num_index

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

Group number index for each variable.

Definition at line 114 of file ReferenceResidualProblem.h.

Referenced by ReferenceResidualProblem::initialSetup(), and ReferenceResidualProblem::updateReferenceResidual().


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