www.mooseframework.org
GluedContactConstraint.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 
9 
10 #include "SystemBase.h"
11 #include "PenetrationLocator.h"
12 #include "AddVariableAction.h"
13 
14 #include "libmesh/string_to_enum.h"
15 #include "libmesh/sparse_matrix.h"
16 
17 template <>
18 InputParameters
20 {
21  MooseEnum orders(AddVariableAction::getNonlinearVariableOrders());
22 
23  InputParameters params = validParams<SparsityBasedContactConstraint>();
24  params.addRequiredParam<BoundaryName>("boundary", "The master boundary");
25  params.addRequiredParam<BoundaryName>("slave", "The slave boundary");
26  params.addRequiredParam<unsigned int>("component",
27  "An integer corresponding to the direction "
28  "the variable this kernel acts in. (0 for x, "
29  "1 for y, 2 for z)");
30 
31  params.addCoupledVar("disp_x", "The x displacement");
32  params.addCoupledVar("disp_y", "The y displacement");
33  params.addCoupledVar("disp_z", "The z displacement");
34 
35  params.addCoupledVar(
36  "displacements",
37  "The displacements appropriate for the simulation geometry and coordinate system");
38 
39  params.addRequiredCoupledVar("nodal_area", "The nodal area");
40  params.addParam<std::string>("model", "frictionless", "The contact model to use");
41 
42  params.set<bool>("use_displaced_mesh") = true;
43  params.addParam<Real>(
44  "penalty",
45  1e8,
46  "The penalty to apply. This can vary depending on the stiffness of your materials");
47  params.addParam<Real>("friction_coefficient", 0.0, "The friction coefficient");
48  params.addParam<Real>("tangential_tolerance",
49  "Tangential distance to extend edges of contact surfaces");
50  params.addParam<Real>(
51  "normal_smoothing_distance",
52  "Distance from edge in parametric coordinates over which to smooth contact normal");
53  params.addParam<std::string>("normal_smoothing_method",
54  "Method to use to smooth normals (edge_based|nodal_normal_based)");
55  params.addParam<MooseEnum>("order", orders, "The finite element order");
56 
57  params.addParam<Real>("tension_release",
58  0.0,
59  "Tension release threshold. A node in contact "
60  "will not be released if its tensile load is below "
61  "this value. Must be positive.");
62 
63  params.addParam<std::string>("formulation", "default", "The contact formulation");
64  return params;
65 }
66 
67 GluedContactConstraint::GluedContactConstraint(const InputParameters & parameters)
68  : SparsityBasedContactConstraint(parameters),
69  _component(getParam<unsigned int>("component")),
70  _model(ContactMaster::contactModel(getParam<std::string>("model"))),
71  _formulation(ContactMaster::contactFormulation(getParam<std::string>("formulation"))),
72  _penalty(getParam<Real>("penalty")),
73  _friction_coefficient(getParam<Real>("friction_coefficient")),
74  _tension_release(getParam<Real>("tension_release")),
75  _updateContactSet(true),
76  _residual_copy(_sys.residualGhosted()),
77  _vars(3, libMesh::invalid_uint),
78  _nodal_area_var(getVar("nodal_area", 0)),
79  _aux_system(_nodal_area_var->sys()),
80  _aux_solution(_aux_system.currentSolution())
81 {
82  if (parameters.isParamValid("tangential_tolerance"))
83  _penetration_locator.setTangentialTolerance(getParam<Real>("tangential_tolerance"));
84 
85  if (parameters.isParamValid("normal_smoothing_distance"))
86  _penetration_locator.setNormalSmoothingDistance(getParam<Real>("normal_smoothing_distance"));
87 
88  if (parameters.isParamValid("normal_smoothing_method"))
89  _penetration_locator.setNormalSmoothingMethod(
90  parameters.get<std::string>("normal_smoothing_method"));
91 
92  if (isParamValid("displacements"))
93  {
94  // modern parameter scheme for displacements
95  for (unsigned int i = 0; i < coupledComponents("displacements"); ++i)
96  _vars[i] = coupled("displacements", i);
97  }
98  else
99  {
100  // Legacy parameter scheme for displacements
101  if (isParamValid("disp_x"))
102  _vars[0] = coupled("disp_x");
103  if (isParamValid("disp_y"))
104  _vars[1] = coupled("disp_y");
105  if (isParamValid("disp_z"))
106  _vars[2] = coupled("disp_z");
107 
108  mooseDeprecated("use the `displacements` parameter rather than the `disp_*` parameters (those "
109  "will go away with the deprecation of the Solid Mechanics module).");
110  }
111 
112  _penetration_locator.setUpdate(false);
113 }
114 
115 void
117 {
118  if (_component == 0)
119  {
120  updateContactSet(true);
121  _updateContactSet = false;
122  }
123 }
124 
125 void
127 {
128  if (_component == 0)
129  {
130  if (_updateContactSet)
132 
133  _updateContactSet = true;
134  }
135 }
136 
137 void
139 {
140  auto & has_penetrated = _penetration_locator._has_penetrated;
141 
142  for (auto & pinfo_pair : _penetration_locator._penetration_info)
143  {
144  PenetrationInfo * pinfo = pinfo_pair.second;
145 
146  // Skip this pinfo if there are no DOFs on this node.
147  if (!pinfo || pinfo->_node->n_comp(_sys.number(), _vars[_component]) < 1)
148  continue;
149 
150  const dof_id_type slave_node_num = pinfo_pair.first;
151  auto hpit = has_penetrated.find(slave_node_num);
152 
153  if (beginning_of_step)
154  {
155  pinfo->_starting_elem = pinfo->_elem;
156  pinfo->_starting_side_num = pinfo->_side_num;
157  pinfo->_starting_closest_point_ref = pinfo->_closest_point_ref;
158  }
159 
160  if (pinfo->_distance >= 0)
161  if (hpit == has_penetrated.end())
162  has_penetrated.insert(slave_node_num);
163  }
164 }
165 
166 bool
168 {
169  auto hpit = _penetration_locator._has_penetrated.find(_current_node->id());
170  return hpit != _penetration_locator._has_penetrated.end();
171 }
172 
173 Real
175 {
176  return _u_slave[_qp];
177 }
178 
179 Real
180 GluedContactConstraint::computeQpResidual(Moose::ConstraintType type)
181 {
182  switch (type)
183  {
184  case Moose::Slave:
185  {
186  PenetrationInfo * pinfo = _penetration_locator._penetration_info[_current_node->id()];
187  Real distance = (*_current_node)(_component)-pinfo->_closest_point(_component);
188  Real pen_force = _penalty * distance;
189 
190  Real resid = pen_force;
191  pinfo->_contact_force(_component) = resid;
192  pinfo->_mech_status = PenetrationInfo::MS_STICKING;
193 
194  return _test_slave[_i][_qp] * resid;
195  }
196 
197  case Moose::Master:
198  {
199  PenetrationInfo * pinfo = _penetration_locator._penetration_info[_current_node->id()];
200 
201  dof_id_type dof_number = _current_node->dof_number(0, _vars[_component], 0);
202  Real resid = _residual_copy(dof_number);
203 
204  pinfo->_contact_force(_component) = -resid;
205  pinfo->_mech_status = PenetrationInfo::MS_STICKING;
206 
207  return _test_master[_i][_qp] * resid;
208  }
209  }
210 
211  return 0.0;
212 }
213 
214 Real
215 GluedContactConstraint::computeQpJacobian(Moose::ConstraintJacobianType type)
216 {
217  switch (type)
218  {
219  case Moose::SlaveSlave:
220  return _penalty * _phi_slave[_j][_qp] * _test_slave[_i][_qp];
221 
222  case Moose::SlaveMaster:
223  return _penalty * -_phi_master[_j][_qp] * _test_slave[_i][_qp];
224 
225  case Moose::MasterSlave:
226  {
227  const double slave_jac = (*_jacobian)(_current_node->dof_number(0, _vars[_component], 0),
228  _connected_dof_indices[_j]);
229  return slave_jac * _test_master[_i][_qp];
230  }
231 
232  case Moose::MasterMaster:
233  return 0.0;
234  }
235 
236  return 0.0;
237 }
238 
239 Real
240 GluedContactConstraint::computeQpOffDiagJacobian(Moose::ConstraintJacobianType type,
241  unsigned int /*jvar*/)
242 {
243  Real retVal = 0.0;
244 
245  switch (type)
246  {
247  case Moose::SlaveSlave:
248  case Moose::SlaveMaster:
249  case Moose::MasterMaster:
250  break;
251 
252  case Moose::MasterSlave:
253  {
254  const double slave_jac = (*_jacobian)(_current_node->dof_number(0, _vars[_component], 0),
255  _connected_dof_indices[_j]);
256  retVal = slave_jac * _test_master[_i][_qp];
257  break;
258  }
259  }
260 
261  return retVal;
262 }
std::vector< unsigned int > _vars
const unsigned int _component
InputParameters validParams< SparsityBasedContactConstraint >()
ContactFormulation contactFormulation(const std::string &the_name)
virtual void updateContactSet(bool beginning_of_step=false)
virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType type, unsigned int jvar)
Compute off-diagonal Jacobian entries.
virtual Real computeQpResidual(Moose::ConstraintType type)
InputParameters validParams< GluedContactConstraint >()
NumericVector< Number > & _residual_copy
ContactModel contactModel(const std::string &the_name)
virtual Real computeQpJacobian(Moose::ConstraintJacobianType type)
GluedContactConstraint(const InputParameters &parameters)