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

A MechanicalContactConstraint forces the value of a variable to be the same on both sides of an interface. More...

#include <MechanicalContactConstraint.h>

Inheritance diagram for MechanicalContactConstraint:
[legend]

Public Member Functions

 MechanicalContactConstraint (const InputParameters &parameters)
 
virtual ~MechanicalContactConstraint ()
 
virtual void timestepSetup () override
 
virtual void jacobianSetup () override
 
virtual bool AugmentedLagrangianContactConverged ()
 
virtual void updateAugmentedLagrangianMultiplier (bool beginning_of_step=false)
 
virtual void updateContactSet (bool beginning_of_step=false)
 
virtual Real computeQpSlaveValue () override
 
virtual Real computeQpResidual (Moose::ConstraintType type) override
 
virtual void computeJacobian () override
 Computes the jacobian for the current element. More...
 
virtual void computeOffDiagJacobian (unsigned int jvar) override
 Compute off-diagonal Jacobian entries. More...
 
virtual Real computeQpJacobian (Moose::ConstraintJacobianType type) override
 
virtual Real computeQpOffDiagJacobian (Moose::ConstraintJacobianType type, unsigned int jvar) override
 Compute off-diagonal Jacobian entries. More...
 
virtual void getConnectedDofIndices (unsigned int var_num) override
 Get the dof indices of the nodes connected to the slave node for a specific variable. More...
 
bool getCoupledVarComponent (unsigned int var_num, unsigned int &component)
 Determine whether the coupled variable is one of the displacement variables, and find its component. More...
 
virtual bool addCouplingEntriesToJacobian () override
 
bool shouldApply () override
 
void computeContactForce (PenetrationInfo *pinfo)
 

Protected Member Functions

Real nodalArea (PenetrationInfo &pinfo)
 
Real getPenalty (PenetrationInfo &pinfo)
 
Real getTangentialPenalty (PenetrationInfo &pinfo)
 

Protected Attributes

MooseSharedPointer< DisplacedProblem > _displaced_problem
 
const unsigned int _component
 
ContactModel _model
 
const ContactFormulation _formulation
 
const bool _normalize_penalty
 
const Real _penalty
 
Real _penalty_tangential
 
const Real _friction_coefficient
 
const Real _tension_release
 
const Real _capture_tolerance
 
const unsigned int _stick_lock_iterations
 
const Real _stick_unlock_factor
 
bool _update_contact_set
 
NumericVector< Number > & _residual_copy
 
const unsigned int _mesh_dimension
 
std::vector< unsigned int > _vars
 
MooseVariable * _nodal_area_var
 
SystemBase & _aux_system
 
const NumericVector< Number > * _aux_solution
 
const bool _master_slave_jacobian
 Whether to include coupling between the master and slave nodes in the Jacobian. More...
 
const bool _connected_slave_nodes_jacobian
 Whether to include coupling terms with the nodes connected to the slave nodes in the Jacobian. More...
 
const bool _non_displacement_vars_jacobian
 Whether to include coupling terms with non-displacement variables in the Jacobian. More...
 
Real _al_penetration_tolerance
 The tolerance of the penetration for augmented Lagrangian method. More...
 
Real _al_incremental_slip_tolerance
 The tolerance of the incremental slip for augmented Lagrangian method. More...
 
Real _al_frictional_force_tolerance
 The tolerance of the frictional force for augmented Lagrangian method. More...
 

Detailed Description

A MechanicalContactConstraint forces the value of a variable to be the same on both sides of an interface.

Definition at line 25 of file MechanicalContactConstraint.h.

Constructor & Destructor Documentation

MechanicalContactConstraint::MechanicalContactConstraint ( const InputParameters &  parameters)

Definition at line 111 of file MechanicalContactConstraint.C.

112  : NodeFaceConstraint(parameters),
113  _displaced_problem(parameters.get<FEProblemBase *>("_fe_problem_base")->getDisplacedProblem()),
114  _component(getParam<unsigned int>("component")),
115  _model(ContactMaster::contactModel(getParam<std::string>("model"))),
116  _formulation(ContactMaster::contactFormulation(getParam<std::string>("formulation"))),
117  _normalize_penalty(getParam<bool>("normalize_penalty")),
118  _penalty(getParam<Real>("penalty")),
119  _friction_coefficient(getParam<Real>("friction_coefficient")),
120  _tension_release(getParam<Real>("tension_release")),
121  _capture_tolerance(getParam<Real>("capture_tolerance")),
122  _stick_lock_iterations(getParam<unsigned int>("stick_lock_iterations")),
123  _stick_unlock_factor(getParam<Real>("stick_unlock_factor")),
124  _update_contact_set(true),
125  _residual_copy(_sys.residualGhosted()),
126  _mesh_dimension(_mesh.dimension()),
127  _vars(3, libMesh::invalid_uint),
128  _nodal_area_var(getVar("nodal_area", 0)),
130  _aux_solution(_aux_system.currentSolution()),
131  _master_slave_jacobian(getParam<bool>("master_slave_jacobian")),
132  _connected_slave_nodes_jacobian(getParam<bool>("connected_slave_nodes_jacobian")),
133  _non_displacement_vars_jacobian(getParam<bool>("non_displacement_variables_jacobian"))
134 {
135  _overwrite_slave_residual = false;
136 
137  if (isParamValid("displacements"))
138  {
139  // modern parameter scheme for displacements
140  for (unsigned int i = 0; i < coupledComponents("displacements"); ++i)
141  _vars[i] = coupled("displacements", i);
142  }
143  else
144  {
145  // Legacy parameter scheme for displacements
146  if (isParamValid("disp_x"))
147  _vars[0] = coupled("disp_x");
148  if (isParamValid("disp_y"))
149  _vars[1] = coupled("disp_y");
150  if (isParamValid("disp_z"))
151  _vars[2] = coupled("disp_z");
152 
153  mooseDeprecated("use the `displacements` parameter rather than the `disp_*` parameters (those "
154  "will go away with the deprecation of the Solid Mechanics module).");
155  }
156 
157  if (parameters.isParamValid("tangential_tolerance"))
158  _penetration_locator.setTangentialTolerance(getParam<Real>("tangential_tolerance"));
159 
160  if (parameters.isParamValid("normal_smoothing_distance"))
161  _penetration_locator.setNormalSmoothingDistance(getParam<Real>("normal_smoothing_distance"));
162 
163  if (parameters.isParamValid("normal_smoothing_method"))
164  _penetration_locator.setNormalSmoothingMethod(
165  parameters.get<std::string>("normal_smoothing_method"));
166 
168  mooseError("The 'tangential_penalty' formulation can only be used with the 'coulomb' model");
169 
170  if (_model == CM_GLUED)
171  _penetration_locator.setUpdate(false);
172 
173  if (_friction_coefficient < 0)
174  mooseError("The friction coefficient must be nonnegative");
175 
176  // set _penalty_tangential to the value of _penalty for now
178 
180  {
181  if (_model == CM_GLUED)
182  mooseError("The Augmented Lagrangian contact formulation does not support GLUED case.");
183 
184  FEProblemBase * fe_problem = getParam<FEProblemBase *>("_fe_problem_base");
185  if (dynamic_cast<AugmentedLagrangianContactProblem *>(fe_problem) == NULL)
186  mooseError("The Augmented Lagrangian contact formulation must use "
187  "AugmentedLagrangianContactProblem.");
188 
189  if (!parameters.isParamValid("al_penetration_tolerance"))
190  mooseError("For Augmented Lagrangian contact, al_penetration_tolerance must be provided.");
191  else
192  _al_penetration_tolerance = parameters.get<Real>("al_penetration_tolerance");
193 
194  if (_model != CM_FRICTIONLESS)
195  {
196  if (!parameters.isParamValid("al_incremental_slip_tolerance") ||
197  !parameters.isParamValid("al_frictional_force_tolerance"))
198  {
199  mooseError("For the Augmented Lagrangian frictional contact formualton, "
200  "al_incremental_slip_tolerance and "
201  "al_frictional_force_tolerance must be provided.");
202  }
203  else
204  {
205  _al_incremental_slip_tolerance = parameters.get<Real>("al_incremental_slip_tolerance");
206  _al_frictional_force_tolerance = parameters.get<Real>("al_frictional_force_tolerance");
207  }
208  }
209  }
210 }
NumericVector< Number > & _residual_copy
const ContactFormulation _formulation
Real _al_penetration_tolerance
The tolerance of the penetration for augmented Lagrangian method.
const bool _non_displacement_vars_jacobian
Whether to include coupling terms with non-displacement variables in the Jacobian.
MooseSharedPointer< DisplacedProblem > _displaced_problem
const bool _connected_slave_nodes_jacobian
Whether to include coupling terms with the nodes connected to the slave nodes in the Jacobian...
static ContactFormulation contactFormulation(std::string name)
std::vector< unsigned int > _vars
Real _al_incremental_slip_tolerance
The tolerance of the incremental slip for augmented Lagrangian method.
static ContactModel contactModel(std::string name)
Real _al_frictional_force_tolerance
The tolerance of the frictional force for augmented Lagrangian method.
const NumericVector< Number > * _aux_solution
const bool _master_slave_jacobian
Whether to include coupling between the master and slave nodes in the Jacobian.
virtual MechanicalContactConstraint::~MechanicalContactConstraint ( )
inlinevirtual

Definition at line 29 of file MechanicalContactConstraint.h.

29 {}

Member Function Documentation

virtual bool MechanicalContactConstraint::addCouplingEntriesToJacobian ( )
inlineoverridevirtual

Definition at line 81 of file MechanicalContactConstraint.h.

81 { return _master_slave_jacobian; }
const bool _master_slave_jacobian
Whether to include coupling between the master and slave nodes in the Jacobian.
bool MechanicalContactConstraint::AugmentedLagrangianContactConverged ( )
virtual

Definition at line 313 of file MechanicalContactConstraint.C.

Referenced by ~MechanicalContactConstraint().

314 {
315  Real contactResidual = 0.0;
316  unsigned int converged = 0;
317 
318  for (auto & pinfo_pair : _penetration_locator._penetration_info)
319  {
320  const dof_id_type slave_node_num = pinfo_pair.first;
321  PenetrationInfo * pinfo = pinfo_pair.second;
322 
323  // Skip this pinfo if there are no DOFs on this node.
324  if (!pinfo || pinfo->_node->n_comp(_sys.number(), _vars[_component]) < 1)
325  continue;
326 
327  const Real distance = pinfo->_normal * (pinfo->_closest_point - _mesh.nodeRef(slave_node_num));
328 
329  if (pinfo->isCaptured())
330  {
331  if (contactResidual < std::abs(distance))
332  contactResidual = std::abs(distance);
333 
334  // penetration < tol
335  if (contactResidual > _al_penetration_tolerance)
336  {
337  converged = 1;
338  break;
339  }
340 
341  if (_model == CM_COULOMB)
342  {
343  RealVectorValue contact_force_normal((pinfo->_contact_force * pinfo->_normal) *
344  pinfo->_normal);
345  RealVectorValue contact_force_tangential(pinfo->_contact_force - contact_force_normal);
346 
347  RealVectorValue tangential_inc_slip =
348  pinfo->_incremental_slip - (pinfo->_incremental_slip * pinfo->_normal) * pinfo->_normal;
349 
350  const Real tan_mag(contact_force_tangential.norm());
351  const Real tangential_inc_slip_mag = tangential_inc_slip.norm();
352 
353  RealVectorValue distance_vec =
354  (pinfo->_normal * (_mesh.nodeRef(slave_node_num) - pinfo->_closest_point)) *
355  pinfo->_normal;
356 
357  Real penalty = getPenalty(*pinfo);
358  RealVectorValue pen_force_normal =
359  penalty * distance_vec + pinfo->_lagrange_multiplier * pinfo->_normal;
360 
361  // Frictional capacity
362  Real capacity(_friction_coefficient * (pen_force_normal * pinfo->_normal < 0
363  ? -pen_force_normal * pinfo->_normal
364  : 0.0));
365 
366  // incremental slip <= tol for all pinfo_pair such that tan_mag < capacity
367  if (MooseUtils::absoluteFuzzyLessThan(tan_mag, capacity) &&
368  pinfo->_mech_status == PenetrationInfo::MS_STICKING)
369  {
370  if (MooseUtils::absoluteFuzzyGreaterThan(tangential_inc_slip_mag,
372  {
373  converged = 2;
374  break;
375  }
376  }
377 
378  // for all pinfo_pair, tag_mag should be less than (1 + tol) * (capacity + tol)
379  if (tan_mag >
381  {
382  converged = 3;
383  break;
384  }
385  }
386  }
387  }
388 
389  _communicator.max(converged);
390 
391  if (converged == 1)
392  _console
393  << "The Augmented Lagrangian contact tangential sliding enforcement is NOT satisfied \n";
394  else if (converged == 2)
395  _console
396  << "The Augmented Lagrangian contact tangential sliding enforcement is NOT satisfied \n";
397  else if (converged == 3)
398  _console << "The Augmented Lagrangian contact frictional force enforcement is NOT satisfied \n";
399  else
400  return true;
401 
402  return false;
403 }
Real _al_penetration_tolerance
The tolerance of the penetration for augmented Lagrangian method.
Real getPenalty(PenetrationInfo &pinfo)
std::vector< unsigned int > _vars
Real _al_incremental_slip_tolerance
The tolerance of the incremental slip for augmented Lagrangian method.
Real _al_frictional_force_tolerance
The tolerance of the frictional force for augmented Lagrangian method.
void MechanicalContactConstraint::computeContactForce ( PenetrationInfo *  pinfo)

Definition at line 492 of file MechanicalContactConstraint.C.

Referenced by addCouplingEntriesToJacobian(), and shouldApply().

493 {
494  const Node * node = pinfo->_node;
495 
496  // Build up residual vector
497  RealVectorValue res_vec;
498  for (unsigned int i = 0; i < _mesh_dimension; ++i)
499  {
500  dof_id_type dof_number = node->dof_number(0, _vars[i], 0);
501  res_vec(i) = _residual_copy(dof_number);
502  }
503 
504  RealVectorValue distance_vec(_mesh.nodeRef(node->id()) - pinfo->_closest_point);
505  const Real penalty = getPenalty(*pinfo);
506  const Real penalty_slip = getTangentialPenalty(*pinfo);
507 
508  RealVectorValue pen_force(penalty * distance_vec);
509 
510  switch (_model)
511  {
512  case CM_FRICTIONLESS:
513  switch (_formulation)
514  {
515  case CF_KINEMATIC:
516  pinfo->_contact_force = -pinfo->_normal * (pinfo->_normal * res_vec);
517  break;
518 
519  case CF_PENALTY:
520  pinfo->_contact_force = pinfo->_normal * (pinfo->_normal * pen_force);
521  break;
522 
524  pinfo->_contact_force =
525  (pinfo->_normal *
526  (pinfo->_normal * (pen_force + pinfo->_lagrange_multiplier * pinfo->_normal)));
527  break;
528 
529  default:
530  mooseError("Invalid contact formulation");
531  break;
532  }
533  pinfo->_mech_status = PenetrationInfo::MS_SLIPPING;
534  break;
535  case CM_COULOMB:
536  switch (_formulation)
537  {
538  case CF_KINEMATIC:
539  {
540  // Frictional capacity
541  const Real capacity(_friction_coefficient *
542  (res_vec * pinfo->_normal > 0 ? res_vec * pinfo->_normal : 0));
543 
544  // Normal and tangential components of predictor force
545  pinfo->_contact_force = -res_vec;
546  RealVectorValue contact_force_normal((pinfo->_contact_force * pinfo->_normal) *
547  pinfo->_normal);
548  RealVectorValue contact_force_tangential(pinfo->_contact_force - contact_force_normal);
549 
550  RealVectorValue tangential_inc_slip =
551  pinfo->_incremental_slip -
552  (pinfo->_incremental_slip * pinfo->_normal) * pinfo->_normal;
553 
554  // Magnitude of tangential predictor force
555  const Real tan_mag(contact_force_tangential.size());
556  const Real tangential_inc_slip_mag = tangential_inc_slip.size();
557  const Real slip_tol = capacity / penalty;
558  pinfo->_slip_tol = slip_tol;
559 
560  if ((tangential_inc_slip_mag > slip_tol || tan_mag > capacity) &&
561  (pinfo->_stick_locked_this_step < _stick_lock_iterations ||
562  tan_mag > capacity * _stick_unlock_factor))
563  {
564  if (pinfo->_stick_locked_this_step >= _stick_lock_iterations)
565  pinfo->_stick_locked_this_step = 0;
566 
567  // Check for scenario where node has slipped too far in a previous iteration
568  // for special treatment (only done if the incremental slip is non-trivial)
569  bool slipped_too_far = false;
570  RealVectorValue slip_inc_direction;
571  if (tangential_inc_slip_mag > slip_tol)
572  {
573  slip_inc_direction = tangential_inc_slip / tangential_inc_slip_mag;
574  Real slip_dot_tang_force = slip_inc_direction * contact_force_tangential;
575  if (slip_dot_tang_force < capacity)
576  slipped_too_far = true;
577  }
578 
579  if (slipped_too_far) // slip back along slip increment
580  pinfo->_contact_force = contact_force_normal + capacity * slip_inc_direction;
581  else
582  {
583  if (tan_mag > 0) // slip along tangential force direction
584  pinfo->_contact_force =
585  contact_force_normal + capacity * contact_force_tangential / tan_mag;
586  else // treat as frictionless
587  pinfo->_contact_force = contact_force_normal;
588  }
589  if (capacity == 0)
590  pinfo->_mech_status = PenetrationInfo::MS_SLIPPING;
591  else
592  pinfo->_mech_status = PenetrationInfo::MS_SLIPPING_FRICTION;
593  }
594  else
595  {
596  if (pinfo->_mech_status != PenetrationInfo::MS_STICKING &&
597  pinfo->_mech_status != PenetrationInfo::MS_NO_CONTACT)
598  ++pinfo->_stick_locked_this_step;
599  pinfo->_mech_status = PenetrationInfo::MS_STICKING;
600  }
601  break;
602  }
603 
604  case CF_PENALTY:
605  {
606  distance_vec = pinfo->_incremental_slip +
607  (pinfo->_normal * (_mesh.nodeRef(node->id()) - pinfo->_closest_point)) *
608  pinfo->_normal;
609  pen_force = penalty * distance_vec;
610 
611  // Frictional capacity
612  // const Real capacity( _friction_coefficient * (pen_force * pinfo->_normal < 0 ?
613  // -pen_force * pinfo->_normal : 0) );
614  const Real capacity(_friction_coefficient *
615  (res_vec * pinfo->_normal > 0 ? res_vec * pinfo->_normal : 0));
616 
617  // Elastic predictor
618  pinfo->_contact_force =
619  pen_force + (pinfo->_contact_force_old -
620  pinfo->_normal * (pinfo->_normal * pinfo->_contact_force_old));
621  RealVectorValue contact_force_normal((pinfo->_contact_force * pinfo->_normal) *
622  pinfo->_normal);
623  RealVectorValue contact_force_tangential(pinfo->_contact_force - contact_force_normal);
624 
625  // Tangential magnitude of elastic predictor
626  const Real tan_mag(contact_force_tangential.norm());
627 
628  if (tan_mag > capacity)
629  {
630  pinfo->_contact_force =
631  contact_force_normal + capacity * contact_force_tangential / tan_mag;
632  if (MooseUtils::absoluteFuzzyEqual(capacity, 0))
633  pinfo->_mech_status = PenetrationInfo::MS_SLIPPING;
634  else
635  pinfo->_mech_status = PenetrationInfo::MS_SLIPPING_FRICTION;
636  }
637  else
638  pinfo->_mech_status = PenetrationInfo::MS_STICKING;
639  break;
640  }
641 
643  {
644  distance_vec = (pinfo->_normal * (_mesh.nodeRef(node->id()) - pinfo->_closest_point)) *
645  pinfo->_normal;
646 
647  RealVectorValue contact_force_normal =
648  penalty * distance_vec + pinfo->_lagrange_multiplier * pinfo->_normal;
649 
650  RealVectorValue tangential_inc_slip =
651  pinfo->_incremental_slip -
652  (pinfo->_incremental_slip * pinfo->_normal) * pinfo->_normal;
653 
654  RealVectorValue contact_force_tangential =
655  pinfo->_lagrange_multiplier_slip +
656  (pinfo->_contact_force_old -
657  pinfo->_normal * (pinfo->_normal * pinfo->_contact_force_old));
658 
659  RealVectorValue inc_pen_force_tangential = penalty_slip * tangential_inc_slip;
660 
661  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
662  pinfo->_contact_force =
663  contact_force_normal + contact_force_tangential + inc_pen_force_tangential;
664  else
665  pinfo->_contact_force = contact_force_normal + contact_force_tangential;
666 
667  break;
668  }
669 
671  {
672  // Frictional capacity (kinematic formulation)
673  const Real capacity = _friction_coefficient * std::max(res_vec * pinfo->_normal, 0.0);
674 
675  // Normal component of contact force (kinematic formulation)
676  RealVectorValue contact_force_normal((-res_vec * pinfo->_normal) * pinfo->_normal);
677 
678  // Predictor tangential component of contact force (penalty formulation)
679  RealVectorValue inc_pen_force_tangential = penalty * pinfo->_incremental_slip;
680  RealVectorValue contact_force_tangential =
681  inc_pen_force_tangential +
682  (pinfo->_contact_force_old -
683  pinfo->_normal * (pinfo->_normal * pinfo->_contact_force_old));
684 
685  // Magnitude of tangential predictor
686  const Real tan_mag(contact_force_tangential.norm());
687 
688  if (tan_mag > capacity)
689  {
690  pinfo->_contact_force =
691  contact_force_normal + capacity * contact_force_tangential / tan_mag;
692  if (MooseUtils::absoluteFuzzyEqual(capacity, 0))
693  pinfo->_mech_status = PenetrationInfo::MS_SLIPPING;
694  else
695  pinfo->_mech_status = PenetrationInfo::MS_SLIPPING_FRICTION;
696  }
697  else
698  {
699  pinfo->_contact_force = contact_force_normal + contact_force_tangential;
700  pinfo->_mech_status = PenetrationInfo::MS_STICKING;
701  }
702  break;
703  }
704 
705  default:
706  mooseError("Invalid contact formulation");
707  break;
708  }
709  break;
710 
711  case CM_GLUED:
712  switch (_formulation)
713  {
714  case CF_KINEMATIC:
715  pinfo->_contact_force = -res_vec;
716  break;
717 
718  case CF_PENALTY:
719  pinfo->_contact_force = pen_force;
720  break;
721 
723  pinfo->_contact_force =
724  pen_force + pinfo->_lagrange_multiplier * distance_vec / distance_vec.norm();
725  break;
726 
727  default:
728  mooseError("Invalid contact formulation");
729  break;
730  }
731  pinfo->_mech_status = PenetrationInfo::MS_STICKING;
732  break;
733 
734  default:
735  mooseError("Invalid or unavailable contact model");
736  break;
737  }
738 }
NumericVector< Number > & _residual_copy
const ContactFormulation _formulation
Real getPenalty(PenetrationInfo &pinfo)
Real getTangentialPenalty(PenetrationInfo &pinfo)
std::vector< unsigned int > _vars
void MechanicalContactConstraint::computeJacobian ( )
overridevirtual

Computes the jacobian for the current element.

Definition at line 1632 of file MechanicalContactConstraint.C.

Referenced by ~MechanicalContactConstraint().

1633 {
1634  getConnectedDofIndices(_var.number());
1635 
1636  DenseMatrix<Number> & Knn =
1637  _assembly.jacobianBlockNeighbor(Moose::NeighborNeighbor, _master_var.number(), _var.number());
1638 
1639  _Kee.resize(_test_slave.size(), _connected_dof_indices.size());
1640 
1641  for (_i = 0; _i < _test_slave.size(); _i++)
1642  // Loop over the connected dof indices so we can get all the jacobian contributions
1643  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
1644  _Kee(_i, _j) += computeQpJacobian(Moose::SlaveSlave);
1645 
1647  {
1648  DenseMatrix<Number> & Ken =
1649  _assembly.jacobianBlockNeighbor(Moose::ElementNeighbor, _var.number(), _var.number());
1650  if (Ken.m() && Ken.n())
1651  for (_i = 0; _i < _test_slave.size(); _i++)
1652  for (_j = 0; _j < _phi_master.size(); _j++)
1653  Ken(_i, _j) += computeQpJacobian(Moose::SlaveMaster);
1654 
1655  _Kne.resize(_test_master.size(), _connected_dof_indices.size());
1656  for (_i = 0; _i < _test_master.size(); _i++)
1657  // Loop over the connected dof indices so we can get all the jacobian contributions
1658  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
1659  _Kne(_i, _j) += computeQpJacobian(Moose::MasterSlave);
1660  }
1661 
1662  if (Knn.m() && Knn.n())
1663  for (_i = 0; _i < _test_master.size(); _i++)
1664  for (_j = 0; _j < _phi_master.size(); _j++)
1665  Knn(_i, _j) += computeQpJacobian(Moose::MasterMaster);
1666 }
virtual void getConnectedDofIndices(unsigned int var_num) override
Get the dof indices of the nodes connected to the slave node for a specific variable.
virtual Real computeQpJacobian(Moose::ConstraintJacobianType type) override
const bool _master_slave_jacobian
Whether to include coupling between the master and slave nodes in the Jacobian.
void MechanicalContactConstraint::computeOffDiagJacobian ( unsigned int  jvar)
overridevirtual

Compute off-diagonal Jacobian entries.

Parameters
jvarThe index of the coupled variable

Definition at line 1669 of file MechanicalContactConstraint.C.

Referenced by ~MechanicalContactConstraint().

1670 {
1671  getConnectedDofIndices(jvar);
1672 
1673  _Kee.resize(_test_slave.size(), _connected_dof_indices.size());
1674 
1675  DenseMatrix<Number> & Knn =
1676  _assembly.jacobianBlockNeighbor(Moose::NeighborNeighbor, _master_var.number(), jvar);
1677 
1678  for (_i = 0; _i < _test_slave.size(); _i++)
1679  // Loop over the connected dof indices so we can get all the jacobian contributions
1680  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
1681  _Kee(_i, _j) += computeQpOffDiagJacobian(Moose::SlaveSlave, jvar);
1682 
1684  {
1685  DenseMatrix<Number> & Ken =
1686  _assembly.jacobianBlockNeighbor(Moose::ElementNeighbor, _var.number(), jvar);
1687  for (_i = 0; _i < _test_slave.size(); _i++)
1688  for (_j = 0; _j < _phi_master.size(); _j++)
1689  Ken(_i, _j) += computeQpOffDiagJacobian(Moose::SlaveMaster, jvar);
1690 
1691  _Kne.resize(_test_master.size(), _connected_dof_indices.size());
1692  if (_Kne.m() && _Kne.n())
1693  for (_i = 0; _i < _test_master.size(); _i++)
1694  // Loop over the connected dof indices so we can get all the jacobian contributions
1695  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
1696  _Kne(_i, _j) += computeQpOffDiagJacobian(Moose::MasterSlave, jvar);
1697  }
1698 
1699  for (_i = 0; _i < _test_master.size(); _i++)
1700  for (_j = 0; _j < _phi_master.size(); _j++)
1701  Knn(_i, _j) += computeQpOffDiagJacobian(Moose::MasterMaster, jvar);
1702 }
virtual void getConnectedDofIndices(unsigned int var_num) override
Get the dof indices of the nodes connected to the slave node for a specific variable.
virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType type, unsigned int jvar) override
Compute off-diagonal Jacobian entries.
const bool _master_slave_jacobian
Whether to include coupling between the master and slave nodes in the Jacobian.
Real MechanicalContactConstraint::computeQpJacobian ( Moose::ConstraintJacobianType  type)
overridevirtual

Definition at line 793 of file MechanicalContactConstraint.C.

Referenced by computeJacobian(), and ~MechanicalContactConstraint().

794 {
795  PenetrationInfo * pinfo = _penetration_locator._penetration_info[_current_node->id()];
796 
797  const Real penalty = getPenalty(*pinfo);
798  const Real penalty_slip = getTangentialPenalty(*pinfo);
799 
800  switch (type)
801  {
802  case Moose::SlaveSlave:
803  switch (_model)
804  {
805  case CM_FRICTIONLESS:
806  switch (_formulation)
807  {
808  case CF_KINEMATIC:
809  {
810  RealVectorValue jac_vec;
811  for (unsigned int i = 0; i < _mesh_dimension; ++i)
812  {
813  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
814  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
815  }
816  return -pinfo->_normal(_component) * (pinfo->_normal * jac_vec) +
817  (_phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp]) *
818  pinfo->_normal(_component) * pinfo->_normal(_component);
819  }
820 
821  case CF_PENALTY:
823  return _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp] *
824  pinfo->_normal(_component) * pinfo->_normal(_component);
825 
826  default:
827  mooseError("Invalid contact formulation");
828  }
829 
830  case CM_COULOMB:
831  switch (_formulation)
832  {
833  case CF_KINEMATIC:
834  {
835  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
836  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
837  {
838  RealVectorValue jac_vec;
839  for (unsigned int i = 0; i < _mesh_dimension; ++i)
840  {
841  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
842  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
843  }
844  return -pinfo->_normal(_component) * (pinfo->_normal * jac_vec) +
845  (_phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp]) *
846  pinfo->_normal(_component) * pinfo->_normal(_component);
847  }
848  else
849  {
850  const Real curr_jac = (*_jacobian)(
851  _current_node->dof_number(0, _vars[_component], 0), _connected_dof_indices[_j]);
852  return -curr_jac + _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp];
853  }
854  }
855 
856  case CF_PENALTY:
857  {
858  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
859  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
860  return _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp] *
861  pinfo->_normal(_component) * pinfo->_normal(_component);
862  else
863  return _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp];
864  }
866  {
867  Real normal_comp = _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp] *
868  pinfo->_normal(_component) * pinfo->_normal(_component);
869 
870  Real tang_comp = 0.0;
871  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
872  tang_comp = _phi_slave[_j][_qp] * penalty_slip * _test_slave[_i][_qp] *
873  (1.0 - pinfo->_normal(_component) * pinfo->_normal(_component));
874  return normal_comp + tang_comp;
875  }
876 
878  {
879  RealVectorValue jac_vec;
880  for (unsigned int i = 0; i < _mesh_dimension; ++i)
881  {
882  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
883  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
884  }
885  Real normal_comp = -pinfo->_normal(_component) * (pinfo->_normal * jac_vec) +
886  (_phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp]) *
887  pinfo->_normal(_component) * pinfo->_normal(_component);
888 
889  Real tang_comp = 0.0;
890  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
891  tang_comp = _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp] *
892  (1.0 - pinfo->_normal(_component) * pinfo->_normal(_component));
893 
894  return normal_comp + tang_comp;
895  }
896 
897  default:
898  mooseError("Invalid contact formulation");
899  }
900 
901  case CM_GLUED:
902  switch (_formulation)
903  {
904  case CF_KINEMATIC:
905  {
906  const Real curr_jac = (*_jacobian)(_current_node->dof_number(0, _vars[_component], 0),
907  _connected_dof_indices[_j]);
908  return -curr_jac + _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp];
909  }
910 
911  case CF_PENALTY:
913  return _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp];
914 
915  default:
916  mooseError("Invalid contact formulation");
917  }
918  default:
919  mooseError("Invalid or unavailable contact model");
920  }
921 
922  case Moose::SlaveMaster:
923  switch (_model)
924  {
925  case CM_FRICTIONLESS:
926  switch (_formulation)
927  {
928  case CF_KINEMATIC:
929  {
930  Node * curr_master_node = _current_master->get_node(_j);
931 
932  RealVectorValue jac_vec;
933  for (unsigned int i = 0; i < _mesh_dimension; ++i)
934  {
935  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
936  jac_vec(i) =
937  (*_jacobian)(dof_number, curr_master_node->dof_number(0, _vars[_component], 0));
938  }
939  return -pinfo->_normal(_component) * (pinfo->_normal * jac_vec) -
940  (_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp]) *
941  pinfo->_normal(_component) * pinfo->_normal(_component);
942  }
943 
944  case CF_PENALTY:
946  return -_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp] *
947  pinfo->_normal(_component) * pinfo->_normal(_component);
948 
949  default:
950  mooseError("Invalid contact formulation");
951  }
952 
953  case CM_COULOMB:
954  switch (_formulation)
955  {
956  case CF_KINEMATIC:
957  {
958  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
959  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
960  {
961  Node * curr_master_node = _current_master->get_node(_j);
962 
963  RealVectorValue jac_vec;
964  for (unsigned int i = 0; i < _mesh_dimension; ++i)
965  {
966  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
967  jac_vec(i) = (*_jacobian)(dof_number,
968  curr_master_node->dof_number(0, _vars[_component], 0));
969  }
970  return -pinfo->_normal(_component) * (pinfo->_normal * jac_vec) -
971  (_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp]) *
972  pinfo->_normal(_component) * pinfo->_normal(_component);
973  }
974  else
975  {
976  Node * curr_master_node = _current_master->get_node(_j);
977  const Real curr_jac =
978  (*_jacobian)(_current_node->dof_number(0, _vars[_component], 0),
979  curr_master_node->dof_number(0, _vars[_component], 0));
980  return -curr_jac - _phi_master[_j][_qp] * penalty * _test_slave[_i][_qp];
981  }
982  }
983 
984  case CF_PENALTY:
985  {
986  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
987  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
988  return -_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp] *
989  pinfo->_normal(_component) * pinfo->_normal(_component);
990  else
991  return -_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp];
992  }
994  {
995  Real normal_comp = -_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp] *
996  pinfo->_normal(_component) * pinfo->_normal(_component);
997 
998  Real tang_comp = 0.0;
999  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
1000  tang_comp = -_phi_master[_j][_qp] * penalty_slip * _test_slave[_i][_qp] *
1001  (1.0 - pinfo->_normal(_component) * pinfo->_normal(_component));
1002  return normal_comp + tang_comp;
1003  }
1004 
1005  case CF_TANGENTIAL_PENALTY:
1006  {
1007  Node * curr_master_node = _current_master->get_node(_j);
1008 
1009  RealVectorValue jac_vec;
1010  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1011  {
1012  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1013  jac_vec(i) =
1014  (*_jacobian)(dof_number, curr_master_node->dof_number(0, _vars[_component], 0));
1015  }
1016  Real normal_comp = -pinfo->_normal(_component) * (pinfo->_normal * jac_vec) -
1017  (_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp]) *
1018  pinfo->_normal(_component) * pinfo->_normal(_component);
1019 
1020  Real tang_comp = 0.0;
1021  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
1022  tang_comp = -_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp] *
1023  (1.0 - pinfo->_normal(_component) * pinfo->_normal(_component));
1024 
1025  return normal_comp + tang_comp;
1026  }
1027 
1028  default:
1029  mooseError("Invalid contact formulation");
1030  }
1031  case CM_GLUED:
1032  switch (_formulation)
1033  {
1034  case CF_KINEMATIC:
1035  {
1036  Node * curr_master_node = _current_master->get_node(_j);
1037  const Real curr_jac =
1038  (*_jacobian)(_current_node->dof_number(0, _vars[_component], 0),
1039  curr_master_node->dof_number(0, _vars[_component], 0));
1040  return -curr_jac - _phi_master[_j][_qp] * penalty * _test_slave[_i][_qp];
1041  }
1042 
1043  case CF_PENALTY:
1044  case CF_AUGMENTED_LAGRANGE:
1045  return -_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp];
1046 
1047  default:
1048  mooseError("Invalid contact formulation");
1049  }
1050 
1051  default:
1052  mooseError("Invalid or unavailable contact model");
1053  }
1054 
1055  case Moose::MasterSlave:
1056  switch (_model)
1057  {
1058  case CM_FRICTIONLESS:
1059  switch (_formulation)
1060  {
1061  case CF_KINEMATIC:
1062  {
1063  RealVectorValue jac_vec;
1064  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1065  {
1066  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1067  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
1068  }
1069  return pinfo->_normal(_component) * (pinfo->_normal * jac_vec) *
1070  _test_master[_i][_qp];
1071  }
1072 
1073  case CF_PENALTY:
1074  case CF_AUGMENTED_LAGRANGE:
1075  return -_test_master[_i][_qp] * penalty * _phi_slave[_j][_qp] *
1076  pinfo->_normal(_component) * pinfo->_normal(_component);
1077 
1078  default:
1079  mooseError("Invalid contact formulation");
1080  }
1081  case CM_COULOMB:
1082  switch (_formulation)
1083  {
1084  case CF_KINEMATIC:
1085  {
1086  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1087  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
1088  {
1089  RealVectorValue jac_vec;
1090  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1091  {
1092  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1093  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
1094  }
1095  return pinfo->_normal(_component) * (pinfo->_normal * jac_vec) *
1096  _test_master[_i][_qp];
1097  }
1098  else
1099  {
1100  const Real slave_jac = (*_jacobian)(
1101  _current_node->dof_number(0, _vars[_component], 0), _connected_dof_indices[_j]);
1102  return slave_jac * _test_master[_i][_qp];
1103  }
1104  }
1105 
1106  case CF_PENALTY:
1107  {
1108  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1109  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
1110  return -_test_master[_i][_qp] * penalty * _phi_slave[_j][_qp] *
1111  pinfo->_normal(_component) * pinfo->_normal(_component);
1112  else
1113  return -_test_master[_i][_qp] * penalty * _phi_slave[_j][_qp];
1114  }
1115  case CF_AUGMENTED_LAGRANGE:
1116  {
1117  Real normal_comp = -_phi_slave[_j][_qp] * penalty * _test_master[_i][_qp] *
1118  pinfo->_normal(_component) * pinfo->_normal(_component);
1119 
1120  Real tang_comp = 0.0;
1121  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
1122  tang_comp = -_phi_slave[_j][_qp] * penalty_slip * _test_master[_i][_qp] *
1123  (1.0 - pinfo->_normal(_component) * pinfo->_normal(_component));
1124  return normal_comp + tang_comp;
1125  }
1126 
1127  case CF_TANGENTIAL_PENALTY:
1128  {
1129  RealVectorValue jac_vec;
1130  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1131  {
1132  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1133  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
1134  }
1135  Real normal_comp =
1136  pinfo->_normal(_component) * (pinfo->_normal * jac_vec) * _test_master[_i][_qp];
1137 
1138  Real tang_comp = 0.0;
1139  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
1140  tang_comp = -_test_master[_i][_qp] * penalty * _phi_slave[_j][_qp] *
1141  (1.0 - pinfo->_normal(_component) * pinfo->_normal(_component));
1142 
1143  return normal_comp + tang_comp;
1144  }
1145 
1146  default:
1147  mooseError("Invalid contact formulation");
1148  }
1149 
1150  case CM_GLUED:
1151  switch (_formulation)
1152  {
1153  case CF_KINEMATIC:
1154  {
1155  const Real slave_jac = (*_jacobian)(
1156  _current_node->dof_number(0, _vars[_component], 0), _connected_dof_indices[_j]);
1157  return slave_jac * _test_master[_i][_qp];
1158  }
1159 
1160  case CF_PENALTY:
1161  case CF_AUGMENTED_LAGRANGE:
1162  return -_test_master[_i][_qp] * penalty * _phi_slave[_j][_qp];
1163 
1164  default:
1165  mooseError("Invalid contact formulation");
1166  }
1167 
1168  default:
1169  mooseError("Invalid or unavailable contact model");
1170  }
1171 
1172  case Moose::MasterMaster:
1173  switch (_model)
1174  {
1175  case CM_FRICTIONLESS:
1176  switch (_formulation)
1177  {
1178  case CF_KINEMATIC:
1179  return 0.0;
1180 
1181  case CF_PENALTY:
1182  case CF_AUGMENTED_LAGRANGE:
1183  return _test_master[_i][_qp] * penalty * _phi_master[_j][_qp] *
1184  pinfo->_normal(_component) * pinfo->_normal(_component);
1185 
1186  default:
1187  mooseError("Invalid contact formulation");
1188  }
1189 
1190  case CM_COULOMB:
1191  case CM_GLUED:
1192  switch (_formulation)
1193  {
1194  case CF_KINEMATIC:
1195  return 0.0;
1196 
1197  case CF_PENALTY:
1198  {
1199  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1200  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
1201  return _test_master[_i][_qp] * penalty * _phi_master[_j][_qp] *
1202  pinfo->_normal(_component) * pinfo->_normal(_component);
1203  else
1204  return _test_master[_i][_qp] * penalty * _phi_master[_j][_qp];
1205  }
1206 
1207  case CF_TANGENTIAL_PENALTY:
1208  {
1209  Real tang_comp = 0.0;
1210  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
1211  tang_comp = _test_master[_i][_qp] * penalty * _phi_master[_j][_qp] *
1212  (1.0 - pinfo->_normal(_component) * pinfo->_normal(_component));
1213  return tang_comp; // normal component is zero
1214  }
1215 
1216  case CF_AUGMENTED_LAGRANGE:
1217  {
1218  Real normal_comp = _phi_master[_j][_qp] * penalty * _test_master[_i][_qp] *
1219  pinfo->_normal(_component) * pinfo->_normal(_component);
1220 
1221  Real tang_comp = 0.0;
1222  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
1223  tang_comp = _phi_master[_j][_qp] * penalty_slip * _test_master[_i][_qp] *
1224  (1.0 - pinfo->_normal(_component) * pinfo->_normal(_component));
1225  return normal_comp + tang_comp;
1226  }
1227 
1228  default:
1229  mooseError("Invalid contact formulation");
1230  }
1231 
1232  default:
1233  mooseError("Invalid or unavailable contact model");
1234  }
1235  }
1236 
1237  return 0.0;
1238 }
const ContactFormulation _formulation
Real getPenalty(PenetrationInfo &pinfo)
Real getTangentialPenalty(PenetrationInfo &pinfo)
std::vector< unsigned int > _vars
Real MechanicalContactConstraint::computeQpOffDiagJacobian ( Moose::ConstraintJacobianType  type,
unsigned int  jvar 
)
overridevirtual

Compute off-diagonal Jacobian entries.

Parameters
typeThe type of coupling
jvarThe index of the coupled variable

Definition at line 1241 of file MechanicalContactConstraint.C.

Referenced by computeOffDiagJacobian(), and ~MechanicalContactConstraint().

1243 {
1244  PenetrationInfo * pinfo = _penetration_locator._penetration_info[_current_node->id()];
1245 
1246  const Real penalty = getPenalty(*pinfo);
1247  const Real penalty_slip = getTangentialPenalty(*pinfo);
1248 
1249  unsigned int coupled_component;
1250  Real normal_component_in_coupled_var_dir = 1.0;
1251  if (getCoupledVarComponent(jvar, coupled_component))
1252  normal_component_in_coupled_var_dir = pinfo->_normal(coupled_component);
1253 
1254  switch (type)
1255  {
1256  case Moose::SlaveSlave:
1257  switch (_model)
1258  {
1259  case CM_FRICTIONLESS:
1260  switch (_formulation)
1261  {
1262  case CF_KINEMATIC:
1263  {
1264  RealVectorValue jac_vec;
1265  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1266  {
1267  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1268  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
1269  }
1270  return -pinfo->_normal(_component) * (pinfo->_normal * jac_vec) +
1271  (_phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp]) *
1272  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1273  }
1274 
1275  case CF_PENALTY:
1276  case CF_AUGMENTED_LAGRANGE:
1277  return _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp] *
1278  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1279 
1280  default:
1281  mooseError("Invalid contact formulation");
1282  }
1283 
1284  case CM_COULOMB:
1285  {
1287  (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1288  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION))
1289  {
1290  RealVectorValue jac_vec;
1291  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1292  {
1293  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1294  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
1295  }
1296  return -pinfo->_normal(_component) * (pinfo->_normal * jac_vec) +
1297  (_phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp]) *
1298  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1299  }
1300  else if ((_formulation == CF_PENALTY) &&
1301  (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1302  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION))
1303  return _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp] *
1304  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1305 
1306  else if (_formulation == CF_AUGMENTED_LAGRANGE)
1307  {
1308  Real normal_comp = _phi_slave[_j][_qp] * penalty * _test_slave[_i][_qp] *
1309  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1310 
1311  Real tang_comp = 0.0;
1312  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
1313  tang_comp = _phi_slave[_j][_qp] * penalty_slip * _test_slave[_i][_qp] *
1314  (-pinfo->_normal(_component) * normal_component_in_coupled_var_dir);
1315  return normal_comp + tang_comp;
1316  }
1317  else
1318  {
1319  const Real curr_jac = (*_jacobian)(_current_node->dof_number(0, _vars[_component], 0),
1320  _connected_dof_indices[_j]);
1321  return -curr_jac;
1322  }
1323  }
1324 
1325  case CM_GLUED:
1326  {
1327  const Real curr_jac = (*_jacobian)(_current_node->dof_number(0, _vars[_component], 0),
1328  _connected_dof_indices[_j]);
1329  return -curr_jac;
1330  }
1331  default:
1332  mooseError("Invalid or unavailable contact model");
1333  }
1334 
1335  case Moose::SlaveMaster:
1336  switch (_model)
1337  {
1338  case CM_FRICTIONLESS:
1339  switch (_formulation)
1340  {
1341  case CF_KINEMATIC:
1342  {
1343  Node * curr_master_node = _current_master->get_node(_j);
1344 
1345  RealVectorValue jac_vec;
1346  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1347  {
1348  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1349  jac_vec(i) =
1350  (*_jacobian)(dof_number, curr_master_node->dof_number(0, _vars[_component], 0));
1351  }
1352  return -pinfo->_normal(_component) * (pinfo->_normal * jac_vec) -
1353  (_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp]) *
1354  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1355  }
1356 
1357  case CF_PENALTY:
1358  case CF_AUGMENTED_LAGRANGE:
1359  return -_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp] *
1360  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1361 
1362  default:
1363  mooseError("Invalid contact formulation");
1364  }
1365 
1366  case CM_COULOMB:
1368  (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1369  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION))
1370  {
1371  Node * curr_master_node = _current_master->get_node(_j);
1372 
1373  RealVectorValue jac_vec;
1374  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1375  {
1376  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1377  jac_vec(i) =
1378  (*_jacobian)(dof_number, curr_master_node->dof_number(0, _vars[_component], 0));
1379  }
1380  return -pinfo->_normal(_component) * (pinfo->_normal * jac_vec) -
1381  (_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp]) *
1382  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1383  }
1384  else if ((_formulation == CF_PENALTY) &&
1385  (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1386  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION))
1387  return -_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp] *
1388  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1389 
1390  else if (_formulation == CF_AUGMENTED_LAGRANGE)
1391  {
1392  Real normal_comp = -_phi_master[_j][_qp] * penalty * _test_slave[_i][_qp] *
1393  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1394 
1395  Real tang_comp = 0.0;
1396  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
1397  tang_comp = -_phi_master[_j][_qp] * penalty_slip * _test_slave[_i][_qp] *
1398  (-pinfo->_normal(_component) * normal_component_in_coupled_var_dir);
1399  return normal_comp + tang_comp;
1400  }
1401  else
1402  return 0.0;
1403 
1404  case CM_GLUED:
1405  return 0;
1406 
1407  default:
1408  mooseError("Invalid or unavailable contact model");
1409  }
1410 
1411  case Moose::MasterSlave:
1412  switch (_model)
1413  {
1414  case CM_FRICTIONLESS:
1415  switch (_formulation)
1416  {
1417  case CF_KINEMATIC:
1418  {
1419  RealVectorValue jac_vec;
1420  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1421  {
1422  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1423  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
1424  }
1425  return pinfo->_normal(_component) * (pinfo->_normal * jac_vec) *
1426  _test_master[_i][_qp];
1427  }
1428 
1429  case CF_PENALTY:
1430  case CF_AUGMENTED_LAGRANGE:
1431  return -_test_master[_i][_qp] * penalty * _phi_slave[_j][_qp] *
1432  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1433 
1434  default:
1435  mooseError("Invalid contact formulation");
1436  }
1437  case CM_COULOMB:
1438  switch (_formulation)
1439  {
1440  case CF_KINEMATIC:
1441  {
1442  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1443  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
1444  {
1445  RealVectorValue jac_vec;
1446  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1447  {
1448  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1449  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
1450  }
1451  return pinfo->_normal(_component) * (pinfo->_normal * jac_vec) *
1452  _test_master[_i][_qp];
1453  }
1454  else
1455  {
1456  const Real slave_jac = (*_jacobian)(
1457  _current_node->dof_number(0, _vars[_component], 0), _connected_dof_indices[_j]);
1458  return slave_jac * _test_master[_i][_qp];
1459  }
1460  }
1461 
1462  case CF_PENALTY:
1463  {
1464  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1465  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
1466  return -_test_master[_i][_qp] * penalty * _phi_slave[_j][_qp] *
1467  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1468  else
1469  return 0.0;
1470  }
1471  case CF_AUGMENTED_LAGRANGE:
1472  {
1473  Real normal_comp = -_phi_slave[_j][_qp] * penalty * _test_master[_i][_qp] *
1474  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1475 
1476  Real tang_comp = 0.0;
1477  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
1478  tang_comp = -_phi_slave[_j][_qp] * penalty_slip * _test_master[_i][_qp] *
1479  (-pinfo->_normal(_component) * normal_component_in_coupled_var_dir);
1480  return normal_comp + tang_comp;
1481  }
1482  case CF_TANGENTIAL_PENALTY:
1483  {
1484  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1485  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
1486  {
1487  RealVectorValue jac_vec;
1488  for (unsigned int i = 0; i < _mesh_dimension; ++i)
1489  {
1490  dof_id_type dof_number = _current_node->dof_number(0, _vars[i], 0);
1491  jac_vec(i) = (*_jacobian)(dof_number, _connected_dof_indices[_j]);
1492  }
1493  return pinfo->_normal(_component) * (pinfo->_normal * jac_vec) *
1494  _test_master[_i][_qp];
1495  }
1496  else
1497  return 0.0;
1498  }
1499 
1500  default:
1501  mooseError("Invalid contact formulation");
1502  }
1503 
1504  case CM_GLUED:
1505  switch (_formulation)
1506  {
1507  case CF_KINEMATIC:
1508  {
1509  const Real slave_jac = (*_jacobian)(
1510  _current_node->dof_number(0, _vars[_component], 0), _connected_dof_indices[_j]);
1511  return slave_jac * _test_master[_i][_qp];
1512  }
1513 
1514  case CF_PENALTY:
1515  case CF_AUGMENTED_LAGRANGE:
1516  return 0.0;
1517 
1518  default:
1519  mooseError("Invalid contact formulation");
1520  }
1521 
1522  default:
1523  mooseError("Invalid or unavailable contact model");
1524  }
1525 
1526  case Moose::MasterMaster:
1527  switch (_model)
1528  {
1529  case CM_FRICTIONLESS:
1530  switch (_formulation)
1531  {
1532  case CF_KINEMATIC:
1533  return 0.0;
1534 
1535  case CF_PENALTY:
1536  case CF_AUGMENTED_LAGRANGE:
1537  return _test_master[_i][_qp] * penalty * _phi_master[_j][_qp] *
1538  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1539 
1540  default:
1541  mooseError("Invalid contact formulation");
1542  }
1543 
1544  case CM_COULOMB:
1545  {
1547  {
1548 
1549  Real normal_comp = _phi_master[_j][_qp] * penalty * _test_master[_i][_qp] *
1550  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1551 
1552  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1553  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
1554  return normal_comp;
1555  }
1556  else if (_formulation == CF_PENALTY &&
1557  (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1558  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION))
1559  return _test_master[_i][_qp] * penalty * _phi_master[_j][_qp] *
1560  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1561  else
1562  return 0.0;
1563  }
1564 
1565  case CM_GLUED:
1566  if (_formulation == CF_PENALTY &&
1567  (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
1568  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION))
1569  return _test_master[_i][_qp] * penalty * _phi_master[_j][_qp] *
1570  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1571  else if (_formulation == CF_AUGMENTED_LAGRANGE)
1572  {
1573  Real normal_comp = _phi_master[_j][_qp] * penalty * _test_master[_i][_qp] *
1574  pinfo->_normal(_component) * normal_component_in_coupled_var_dir;
1575 
1576  Real tang_comp = 0.0;
1577  if (pinfo->_mech_status == PenetrationInfo::MS_STICKING)
1578  tang_comp = _phi_master[_j][_qp] * penalty_slip * _test_master[_i][_qp] *
1579  (-pinfo->_normal(_component) * normal_component_in_coupled_var_dir);
1580  return normal_comp + tang_comp;
1581  }
1582  else
1583  return 0.0;
1584 
1585  default:
1586  mooseError("Invalid or unavailable contact model");
1587  }
1588  }
1589 
1590  return 0.0;
1591 }
bool getCoupledVarComponent(unsigned int var_num, unsigned int &component)
Determine whether the coupled variable is one of the displacement variables, and find its component...
const ContactFormulation _formulation
Real getPenalty(PenetrationInfo &pinfo)
Real getTangentialPenalty(PenetrationInfo &pinfo)
std::vector< unsigned int > _vars
Real MechanicalContactConstraint::computeQpResidual ( Moose::ConstraintType  type)
overridevirtual

Definition at line 747 of file MechanicalContactConstraint.C.

Referenced by ~MechanicalContactConstraint().

748 {
749  PenetrationInfo * pinfo = _penetration_locator._penetration_info[_current_node->id()];
750  Real resid = pinfo->_contact_force(_component);
751  switch (type)
752  {
753  case Moose::Slave:
754  if (_formulation == CF_KINEMATIC)
755  {
756  RealVectorValue distance_vec(*_current_node - pinfo->_closest_point);
757  const Real penalty = getPenalty(*pinfo);
758  RealVectorValue pen_force(penalty * distance_vec);
759 
760  if (_model == CM_FRICTIONLESS)
761  resid += pinfo->_normal(_component) * pinfo->_normal * pen_force;
762 
763  else if (_model == CM_COULOMB)
764  {
765  distance_vec = distance_vec - pinfo->_incremental_slip;
766  pen_force = penalty * distance_vec;
767  if (pinfo->_mech_status == PenetrationInfo::MS_SLIPPING ||
768  pinfo->_mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
769  resid += pinfo->_normal(_component) * pinfo->_normal * pen_force;
770  else
771  resid += pen_force(_component);
772  }
773  else if (_model == CM_GLUED)
774  resid += pen_force(_component);
775  }
777  {
778  RealVectorValue distance_vec(*_current_node - pinfo->_closest_point);
779  const Real penalty = getPenalty(*pinfo);
780  RealVectorValue pen_force(penalty * distance_vec);
781  resid += pinfo->_normal(_component) * pinfo->_normal * pen_force;
782  }
783  return _test_slave[_i][_qp] * resid;
784 
785  case Moose::Master:
786  return _test_master[_i][_qp] * -resid;
787  }
788 
789  return 0.0;
790 }
const ContactFormulation _formulation
Real getPenalty(PenetrationInfo &pinfo)
Real MechanicalContactConstraint::computeQpSlaveValue ( )
overridevirtual

Definition at line 741 of file MechanicalContactConstraint.C.

Referenced by ~MechanicalContactConstraint().

742 {
743  return _u_slave[_qp];
744 }
void MechanicalContactConstraint::getConnectedDofIndices ( unsigned int  var_num)
overridevirtual

Get the dof indices of the nodes connected to the slave node for a specific variable.

Parameters
var_numThe number of the variable for which dof indices are gathered
Returns
bool indicating whether the coupled variable is one of the displacement variables

Definition at line 1705 of file MechanicalContactConstraint.C.

Referenced by computeJacobian(), computeOffDiagJacobian(), and ~MechanicalContactConstraint().

1706 {
1707  unsigned int component;
1708  if (getCoupledVarComponent(var_num, component) || _non_displacement_vars_jacobian)
1709  {
1711  NodeFaceConstraint::getConnectedDofIndices(var_num);
1712  else
1713  {
1714  _connected_dof_indices.clear();
1715  MooseVariable & var = _sys.getVariable(0, var_num);
1716  _connected_dof_indices.push_back(var.nodalDofIndex());
1717  }
1718  }
1719 
1720  _phi_slave.resize(_connected_dof_indices.size());
1721 
1722  dof_id_type current_node_var_dof_index = _sys.getVariable(0, var_num).nodalDofIndex();
1723  _qp = 0;
1724 
1725  // Fill up _phi_slave so that it is 1 when j corresponds to the dof associated with this node
1726  // and 0 for every other dof
1727  // This corresponds to evaluating all of the connected shape functions at _this_ node
1728  for (unsigned int j = 0; j < _connected_dof_indices.size(); j++)
1729  {
1730  _phi_slave[j].resize(1);
1731 
1732  if (_connected_dof_indices[j] == current_node_var_dof_index)
1733  _phi_slave[j][_qp] = 1.0;
1734  else
1735  _phi_slave[j][_qp] = 0.0;
1736  }
1737 }
bool getCoupledVarComponent(unsigned int var_num, unsigned int &component)
Determine whether the coupled variable is one of the displacement variables, and find its component...
Real component(const SymmTensor &symm_tensor, unsigned int index)
const bool _non_displacement_vars_jacobian
Whether to include coupling terms with non-displacement variables in the Jacobian.
const bool _connected_slave_nodes_jacobian
Whether to include coupling terms with the nodes connected to the slave nodes in the Jacobian...
const bool _master_slave_jacobian
Whether to include coupling between the master and slave nodes in the Jacobian.
bool MechanicalContactConstraint::getCoupledVarComponent ( unsigned int  var_num,
unsigned int &  component 
)

Determine whether the coupled variable is one of the displacement variables, and find its component.

Parameters
var_numThe number of the variable to be checked
componentThe component index computed in this routine
Returns
bool indicating whether the coupled variable is one of the displacement variables

Definition at line 1740 of file MechanicalContactConstraint.C.

Referenced by computeQpOffDiagJacobian(), getConnectedDofIndices(), and ~MechanicalContactConstraint().

1741 {
1742  component = std::numeric_limits<unsigned int>::max();
1743  bool coupled_var_is_disp_var = false;
1744  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
1745  {
1746  if (var_num == _vars[i])
1747  {
1748  coupled_var_is_disp_var = true;
1749  component = i;
1750  break;
1751  }
1752  }
1753  return coupled_var_is_disp_var;
1754 }
Real component(const SymmTensor &symm_tensor, unsigned int index)
std::vector< unsigned int > _vars
Real MechanicalContactConstraint::getPenalty ( PenetrationInfo &  pinfo)
protected
Real MechanicalContactConstraint::getTangentialPenalty ( PenetrationInfo &  pinfo)
protected

Definition at line 1622 of file MechanicalContactConstraint.C.

Referenced by computeContactForce(), computeQpJacobian(), computeQpOffDiagJacobian(), and updateAugmentedLagrangianMultiplier().

1623 {
1624  Real penalty = _penalty_tangential;
1625  if (_normalize_penalty)
1626  penalty *= nodalArea(pinfo);
1627 
1628  return penalty;
1629 }
Real nodalArea(PenetrationInfo &pinfo)
void MechanicalContactConstraint::jacobianSetup ( )
overridevirtual

Definition at line 226 of file MechanicalContactConstraint.C.

Referenced by ~MechanicalContactConstraint().

227 {
228  if (_component == 0)
229  {
232  _update_contact_set = true;
233  }
234 }
virtual void updateContactSet(bool beginning_of_step=false)
Real MechanicalContactConstraint::nodalArea ( PenetrationInfo &  pinfo)
protected

Definition at line 1594 of file MechanicalContactConstraint.C.

Referenced by getPenalty(), getTangentialPenalty(), and updateContactSet().

1595 {
1596  const Node * node = pinfo._node;
1597 
1598  dof_id_type dof = node->dof_number(_aux_system.number(), _nodal_area_var->number(), 0);
1599 
1600  Real area = (*_aux_solution)(dof);
1601  if (area == 0.0)
1602  {
1603  if (_t_step > 1)
1604  mooseError("Zero nodal area found");
1605  else
1606  area = 1.0; // Avoid divide by zero during initialization
1607  }
1608  return area;
1609 }
bool MechanicalContactConstraint::shouldApply ( )
override

Definition at line 468 of file MechanicalContactConstraint.C.

Referenced by addCouplingEntriesToJacobian().

469 {
470  bool in_contact = false;
471 
472  std::map<dof_id_type, PenetrationInfo *>::iterator found =
473  _penetration_locator._penetration_info.find(_current_node->id());
474  if (found != _penetration_locator._penetration_info.end())
475  {
476  PenetrationInfo * pinfo = found->second;
477  if (pinfo != NULL && pinfo->isCaptured())
478  {
479  in_contact = true;
480 
481  // This computes the contact force once per constraint, rather than once per quad point
482  // and for both master and slave cases.
483  if (_component == 0)
484  computeContactForce(pinfo);
485  }
486  }
487 
488  return in_contact;
489 }
void computeContactForce(PenetrationInfo *pinfo)
void MechanicalContactConstraint::timestepSetup ( )
overridevirtual

Definition at line 213 of file MechanicalContactConstraint.C.

Referenced by ~MechanicalContactConstraint().

214 {
215  if (_component == 0)
216  {
217  updateContactSet(true);
220 
221  _update_contact_set = false;
222  }
223 }
const ContactFormulation _formulation
virtual void updateAugmentedLagrangianMultiplier(bool beginning_of_step=false)
virtual void updateContactSet(bool beginning_of_step=false)
void MechanicalContactConstraint::updateAugmentedLagrangianMultiplier ( bool  beginning_of_step = false)
virtual

Definition at line 237 of file MechanicalContactConstraint.C.

Referenced by timestepSetup(), and ~MechanicalContactConstraint().

238 {
239  for (auto & pinfo_pair : _penetration_locator._penetration_info)
240  {
241  const dof_id_type slave_node_num = pinfo_pair.first;
242  PenetrationInfo * pinfo = pinfo_pair.second;
243 
244  if (!pinfo || pinfo->_node->n_comp(_sys.number(), _vars[_component]) < 1)
245  continue;
246 
247  const Real distance = pinfo->_normal * (pinfo->_closest_point - _mesh.nodeRef(slave_node_num));
248 
249  if (beginning_of_step && _model == CM_COULOMB)
250  {
251  pinfo->_lagrange_multiplier_slip.zero();
252  if (pinfo->isCaptured())
253  pinfo->_mech_status = PenetrationInfo::MS_STICKING;
254  }
255 
256  if (pinfo->isCaptured())
257  {
258  if (_model == CM_FRICTIONLESS)
259  pinfo->_lagrange_multiplier -= getPenalty(*pinfo) * distance;
260 
261  if (_model == CM_COULOMB)
262  {
263  if (!beginning_of_step)
264  {
265  Real penalty = getPenalty(*pinfo);
266  RealVectorValue pen_force_normal =
267  penalty * (-distance) * pinfo->_normal + pinfo->_lagrange_multiplier * pinfo->_normal;
268 
269  // update normal lagrangian multiplier
270  pinfo->_lagrange_multiplier += penalty * (-distance);
271 
272  // Frictional capacity
273  const Real capacity(_friction_coefficient * (pen_force_normal * pinfo->_normal < 0
274  ? -pen_force_normal * pinfo->_normal
275  : 0));
276 
277  RealVectorValue tangential_inc_slip =
278  pinfo->_incremental_slip -
279  (pinfo->_incremental_slip * pinfo->_normal) * pinfo->_normal;
280 
281  Real penalty_slip = getTangentialPenalty(*pinfo);
282 
283  RealVectorValue inc_pen_force_tangential =
284  pinfo->_lagrange_multiplier_slip + penalty_slip * tangential_inc_slip;
285 
286  RealVectorValue tau_old = pinfo->_contact_force_old -
287  pinfo->_normal * (pinfo->_normal * pinfo->_contact_force_old);
288 
289  RealVectorValue contact_force_tangential = inc_pen_force_tangential + tau_old;
290  const Real tan_mag(contact_force_tangential.norm());
291 
292  if (tan_mag > capacity * (_al_frictional_force_tolerance + 1.0))
293  {
294  pinfo->_lagrange_multiplier_slip =
295  -tau_old + capacity * contact_force_tangential / tan_mag;
296  if (MooseUtils::absoluteFuzzyEqual(capacity, 0.0))
297  pinfo->_mech_status = PenetrationInfo::MS_SLIPPING;
298  else
299  pinfo->_mech_status = PenetrationInfo::MS_SLIPPING_FRICTION;
300  }
301  else
302  {
303  pinfo->_mech_status = PenetrationInfo::MS_STICKING;
304  pinfo->_lagrange_multiplier_slip += penalty_slip * tangential_inc_slip;
305  }
306  }
307  }
308  }
309  }
310 }
Real getPenalty(PenetrationInfo &pinfo)
Real getTangentialPenalty(PenetrationInfo &pinfo)
std::vector< unsigned int > _vars
Real _al_frictional_force_tolerance
The tolerance of the frictional force for augmented Lagrangian method.
void MechanicalContactConstraint::updateContactSet ( bool  beginning_of_step = false)
virtual

Definition at line 406 of file MechanicalContactConstraint.C.

Referenced by jacobianSetup(), timestepSetup(), and ~MechanicalContactConstraint().

407 {
408  for (auto & pinfo_pair : _penetration_locator._penetration_info)
409  {
410  const dof_id_type slave_node_num = pinfo_pair.first;
411  PenetrationInfo * pinfo = pinfo_pair.second;
412 
413  // Skip this pinfo if there are no DOFs on this node.
414  if (!pinfo || pinfo->_node->n_comp(_sys.number(), _vars[_component]) < 1)
415  continue;
416 
417  if (beginning_of_step)
418  {
419  if (_app.getExecutioner()->lastSolveConverged())
420  {
421  pinfo->_contact_force_old = pinfo->_contact_force;
422  pinfo->_accumulated_slip_old = pinfo->_accumulated_slip;
423  pinfo->_frictional_energy_old = pinfo->_frictional_energy;
424  pinfo->_mech_status_old = pinfo->_mech_status;
425  }
426  else if (pinfo->_mech_status_old == PenetrationInfo::MS_NO_CONTACT &&
427  pinfo->_mech_status != PenetrationInfo::MS_NO_CONTACT)
428  {
429  // The penetration info object could be based on a bad state so delete it
430  delete pinfo_pair.second;
431  pinfo_pair.second = NULL;
432  continue;
433  }
434 
435  pinfo->_locked_this_step = 0;
436  pinfo->_stick_locked_this_step = 0;
437  pinfo->_starting_elem = pinfo->_elem;
438  pinfo->_starting_side_num = pinfo->_side_num;
439  pinfo->_starting_closest_point_ref = pinfo->_closest_point_ref;
440  }
441  pinfo->_incremental_slip_prev_iter = pinfo->_incremental_slip;
442 
443  const Real contact_pressure = -(pinfo->_normal * pinfo->_contact_force) / nodalArea(*pinfo);
444  const Real distance = pinfo->_normal * (pinfo->_closest_point - _mesh.nodeRef(slave_node_num));
445 
446  // Capture
447  if (!pinfo->isCaptured() &&
448  MooseUtils::absoluteFuzzyGreaterEqual(distance, 0.0, _capture_tolerance))
449  {
450  pinfo->capture();
451 
452  // Increment the lock count every time the node comes back into contact from not being in
453  // contact.
455  ++pinfo->_locked_this_step;
456  }
457  // Release
458  else if (_model != CM_GLUED && pinfo->isCaptured() && _tension_release >= 0.0 &&
459  -contact_pressure >= _tension_release && pinfo->_locked_this_step < 2)
460  {
461  pinfo->release();
462  pinfo->_contact_force.zero();
463  }
464  }
465 }
Real nodalArea(PenetrationInfo &pinfo)
const ContactFormulation _formulation
std::vector< unsigned int > _vars

Member Data Documentation

Real MechanicalContactConstraint::_al_frictional_force_tolerance
protected

The tolerance of the frictional force for augmented Lagrangian method.

Definition at line 129 of file MechanicalContactConstraint.h.

Referenced by AugmentedLagrangianContactConverged(), MechanicalContactConstraint(), and updateAugmentedLagrangianMultiplier().

Real MechanicalContactConstraint::_al_incremental_slip_tolerance
protected

The tolerance of the incremental slip for augmented Lagrangian method.

Definition at line 127 of file MechanicalContactConstraint.h.

Referenced by AugmentedLagrangianContactConverged(), and MechanicalContactConstraint().

Real MechanicalContactConstraint::_al_penetration_tolerance
protected

The tolerance of the penetration for augmented Lagrangian method.

Definition at line 125 of file MechanicalContactConstraint.h.

Referenced by AugmentedLagrangianContactConverged(), and MechanicalContactConstraint().

const NumericVector<Number>* MechanicalContactConstraint::_aux_solution
protected

Definition at line 115 of file MechanicalContactConstraint.h.

SystemBase& MechanicalContactConstraint::_aux_system
protected

Definition at line 114 of file MechanicalContactConstraint.h.

Referenced by nodalArea().

const Real MechanicalContactConstraint::_capture_tolerance
protected

Definition at line 101 of file MechanicalContactConstraint.h.

Referenced by updateContactSet().

const unsigned int MechanicalContactConstraint::_component
protected
const bool MechanicalContactConstraint::_connected_slave_nodes_jacobian
protected

Whether to include coupling terms with the nodes connected to the slave nodes in the Jacobian.

Definition at line 120 of file MechanicalContactConstraint.h.

Referenced by getConnectedDofIndices().

MooseSharedPointer<DisplacedProblem> MechanicalContactConstraint::_displaced_problem
protected

Definition at line 87 of file MechanicalContactConstraint.h.

const ContactFormulation MechanicalContactConstraint::_formulation
protected
const Real MechanicalContactConstraint::_friction_coefficient
protected
const bool MechanicalContactConstraint::_master_slave_jacobian
protected

Whether to include coupling between the master and slave nodes in the Jacobian.

Definition at line 118 of file MechanicalContactConstraint.h.

Referenced by addCouplingEntriesToJacobian(), computeJacobian(), computeOffDiagJacobian(), and getConnectedDofIndices().

const unsigned int MechanicalContactConstraint::_mesh_dimension
protected
ContactModel MechanicalContactConstraint::_model
protected
MooseVariable* MechanicalContactConstraint::_nodal_area_var
protected

Definition at line 113 of file MechanicalContactConstraint.h.

Referenced by nodalArea().

const bool MechanicalContactConstraint::_non_displacement_vars_jacobian
protected

Whether to include coupling terms with non-displacement variables in the Jacobian.

Definition at line 122 of file MechanicalContactConstraint.h.

Referenced by getConnectedDofIndices().

const bool MechanicalContactConstraint::_normalize_penalty
protected

Definition at line 95 of file MechanicalContactConstraint.h.

Referenced by getPenalty(), and getTangentialPenalty().

const Real MechanicalContactConstraint::_penalty
protected

Definition at line 97 of file MechanicalContactConstraint.h.

Referenced by getPenalty(), and MechanicalContactConstraint().

Real MechanicalContactConstraint::_penalty_tangential
protected
NumericVector<Number>& MechanicalContactConstraint::_residual_copy
protected

Definition at line 106 of file MechanicalContactConstraint.h.

Referenced by computeContactForce().

const unsigned int MechanicalContactConstraint::_stick_lock_iterations
protected

Definition at line 102 of file MechanicalContactConstraint.h.

Referenced by computeContactForce().

const Real MechanicalContactConstraint::_stick_unlock_factor
protected

Definition at line 103 of file MechanicalContactConstraint.h.

Referenced by computeContactForce().

const Real MechanicalContactConstraint::_tension_release
protected

Definition at line 100 of file MechanicalContactConstraint.h.

Referenced by updateContactSet().

bool MechanicalContactConstraint::_update_contact_set
protected

Definition at line 104 of file MechanicalContactConstraint.h.

Referenced by jacobianSetup(), and timestepSetup().

std::vector<unsigned int> MechanicalContactConstraint::_vars
protected

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