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

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

#include <GluedContactConstraint.h>

Inheritance diagram for GluedContactConstraint:
[legend]

Public Member Functions

 GluedContactConstraint (const InputParameters &parameters)
 
virtual ~GluedContactConstraint ()
 
virtual void timestepSetup ()
 
virtual void jacobianSetup ()
 
virtual void updateContactSet (bool beginning_of_step=false)
 
virtual Real computeQpSlaveValue ()
 
virtual Real computeQpResidual (Moose::ConstraintType type)
 
virtual Real computeQpJacobian (Moose::ConstraintJacobianType type)
 
virtual Real computeQpOffDiagJacobian (Moose::ConstraintJacobianType type, unsigned int jvar)
 Compute off-diagonal Jacobian entries. More...
 
bool shouldApply ()
 
virtual void getConnectedDofIndices ()
 Gets the indices for all dofs conected to the constraint Get indices for all the slave Jacobian columns, not only those based on the mesh connectivity. More...
 

Protected Attributes

const unsigned int _component
 
const ContactModel _model
 
const ContactFormulation _formulation
 
const Real _penalty
 
const Real _friction_coefficient
 
const Real _tension_release
 
bool _updateContactSet
 
NumericVector< Number > & _residual_copy
 
unsigned int _x_var
 
unsigned int _y_var
 
unsigned int _z_var
 
std::vector< unsigned int > _vars
 
MooseVariable * _nodal_area_var
 
SystemBase & _aux_system
 
const NumericVector< Number > * _aux_solution
 

Detailed Description

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

Definition at line 28 of file GluedContactConstraint.h.

Constructor & Destructor Documentation

◆ GluedContactConstraint()

GluedContactConstraint::GluedContactConstraint ( const InputParameters &  parameters)

Definition at line 71 of file GluedContactConstraint.C.

72  : SparsityBasedContactConstraint(parameters),
73  _component(getParam<unsigned int>("component")),
74  _model(ContactMaster::contactModel(getParam<std::string>("model"))),
75  _formulation(ContactMaster::contactFormulation(getParam<std::string>("formulation"))),
76  _penalty(getParam<Real>("penalty")),
77  _friction_coefficient(getParam<Real>("friction_coefficient")),
78  _tension_release(getParam<Real>("tension_release")),
79  _updateContactSet(true),
80  _residual_copy(_sys.residualGhosted()),
81  _vars(3, libMesh::invalid_uint),
82  _nodal_area_var(getVar("nodal_area", 0)),
84  _aux_solution(_aux_system.currentSolution())
85 {
86  if (parameters.isParamValid("tangential_tolerance"))
87  _penetration_locator.setTangentialTolerance(getParam<Real>("tangential_tolerance"));
88 
89  if (parameters.isParamValid("normal_smoothing_distance"))
90  _penetration_locator.setNormalSmoothingDistance(getParam<Real>("normal_smoothing_distance"));
91 
92  if (parameters.isParamValid("normal_smoothing_method"))
93  _penetration_locator.setNormalSmoothingMethod(
94  parameters.get<std::string>("normal_smoothing_method"));
95 
96  if (isParamValid("displacements"))
97  {
98  // modern parameter scheme for displacements
99  for (unsigned int i = 0; i < coupledComponents("displacements"); ++i)
100  _vars[i] = coupled("displacements", i);
101  }
102  else
103  {
104  // Legacy parameter scheme for displacements
105  if (isParamValid("disp_x"))
106  _vars[0] = coupled("disp_x");
107  if (isParamValid("disp_y"))
108  _vars[1] = coupled("disp_y");
109  if (isParamValid("disp_z"))
110  _vars[2] = coupled("disp_z");
111 
112  mooseDeprecated("use the `displacements` parameter rather than the `disp_*` parameters (those "
113  "will go away with the deprecation of the Solid Mechanics module).");
114  }
115 
116  _penetration_locator.setUpdate(false);
117 }
std::vector< unsigned int > _vars
const unsigned int _component
SparsityBasedContactConstraint(const InputParameters &parameters)
const ContactFormulation _formulation
NumericVector< Number > & _residual_copy
static ContactFormulation contactFormulation(std::string name)
static ContactModel contactModel(std::string name)
const NumericVector< Number > * _aux_solution

◆ ~GluedContactConstraint()

virtual GluedContactConstraint::~GluedContactConstraint ( )
inlinevirtual

Definition at line 32 of file GluedContactConstraint.h.

32 {}

Member Function Documentation

◆ computeQpJacobian()

Real GluedContactConstraint::computeQpJacobian ( Moose::ConstraintJacobianType  type)
virtual

Reimplemented from SparsityBasedContactConstraint.

Definition at line 219 of file GluedContactConstraint.C.

220 {
221  switch (type)
222  {
223  case Moose::SlaveSlave:
224  return _penalty * _phi_slave[_j][_qp] * _test_slave[_i][_qp];
225 
226  case Moose::SlaveMaster:
227  return _penalty * -_phi_master[_j][_qp] * _test_slave[_i][_qp];
228 
229  case Moose::MasterSlave:
230  {
231  const double slave_jac = (*_jacobian)(_current_node->dof_number(0, _vars[_component], 0),
232  _connected_dof_indices[_j]);
233  return slave_jac * _test_master[_i][_qp];
234  }
235 
236  case Moose::MasterMaster:
237  return 0.0;
238  }
239 
240  return 0.0;
241 }
std::vector< unsigned int > _vars
const unsigned int _component

◆ computeQpOffDiagJacobian()

Real GluedContactConstraint::computeQpOffDiagJacobian ( Moose::ConstraintJacobianType  type,
unsigned int  jvar 
)
virtual

Compute off-diagonal Jacobian entries.

Parameters
typeThe type of coupling
jvarThe index of the coupled variable

Definition at line 244 of file GluedContactConstraint.C.

246 {
247  Real retVal = 0.0;
248 
249  switch (type)
250  {
251  case Moose::SlaveSlave:
252  case Moose::SlaveMaster:
253  case Moose::MasterMaster:
254  break;
255 
256  case Moose::MasterSlave:
257  {
258  const double slave_jac = (*_jacobian)(_current_node->dof_number(0, _vars[_component], 0),
259  _connected_dof_indices[_j]);
260  retVal = slave_jac * _test_master[_i][_qp];
261  break;
262  }
263  }
264 
265  return retVal;
266 }
std::vector< unsigned int > _vars
const unsigned int _component

◆ computeQpResidual()

Real GluedContactConstraint::computeQpResidual ( Moose::ConstraintType  type)
virtual

Reimplemented from SparsityBasedContactConstraint.

Definition at line 184 of file GluedContactConstraint.C.

185 {
186  switch (type)
187  {
188  case Moose::Slave:
189  {
190  PenetrationInfo * pinfo = _penetration_locator._penetration_info[_current_node->id()];
191  Real distance = (*_current_node)(_component)-pinfo->_closest_point(_component);
192  Real pen_force = _penalty * distance;
193 
194  Real resid = pen_force;
195  pinfo->_contact_force(_component) = resid;
196  pinfo->_mech_status = PenetrationInfo::MS_STICKING;
197 
198  return _test_slave[_i][_qp] * resid;
199  }
200 
201  case Moose::Master:
202  {
203  PenetrationInfo * pinfo = _penetration_locator._penetration_info[_current_node->id()];
204 
205  dof_id_type dof_number = _current_node->dof_number(0, _vars[_component], 0);
206  Real resid = _residual_copy(dof_number);
207 
208  pinfo->_contact_force(_component) = -resid;
209  pinfo->_mech_status = PenetrationInfo::MS_STICKING;
210 
211  return _test_master[_i][_qp] * resid;
212  }
213  }
214 
215  return 0.0;
216 }
std::vector< unsigned int > _vars
const unsigned int _component
NumericVector< Number > & _residual_copy

◆ computeQpSlaveValue()

Real GluedContactConstraint::computeQpSlaveValue ( )
virtual

Reimplemented from SparsityBasedContactConstraint.

Definition at line 178 of file GluedContactConstraint.C.

179 {
180  return _u_slave[_qp];
181 }

◆ getConnectedDofIndices()

void SparsityBasedContactConstraint::getConnectedDofIndices ( )
virtualinherited

Gets the indices for all dofs conected to the constraint Get indices for all the slave Jacobian columns, not only those based on the mesh connectivity.

Definition at line 29 of file SparsityBasedContactConstraint.C.

30 {
31 #if defined(LIBMESH_HAVE_PETSC) && !PETSC_VERSION_LESS_THAN(3, 3, 0)
32  _connected_dof_indices.clear();
33 
34  // An ugly hack: have to extract the row and look at it's sparsity structure, since otherwise I
35  // won't get the dofs connected to this row by virtue of intervariable coupling
36  // This is easier than sifting through all coupled variables, selecting those active on the
37  // current subdomain, dealing with the scalar variables, etc.
38  // Also, importantly, this will miss coupling to variables that might have introduced by prior
39  // constraints similar to this one!
40  PetscMatrix<Number> * petsc_jacobian = dynamic_cast<PetscMatrix<Number> *>(_jacobian);
41  mooseAssert(petsc_jacobian, "Expected a PETSc matrix");
42  Mat jac = petsc_jacobian->mat();
43  PetscErrorCode ierr;
44  PetscInt ncols;
45  const PetscInt * cols;
46  ierr = MatGetRow(jac, static_cast<PetscInt>(_var.nodalDofIndex()), &ncols, &cols, PETSC_NULL);
47  CHKERRABORT(_communicator.get(), ierr);
48  bool debug = false;
49  if (debug)
50  {
51  libMesh::out << "_connected_dof_indices: adding " << ncols << " dofs from Jacobian row["
52  << _var.nodalDofIndex() << "] = [";
53  }
54  for (PetscInt i = 0; i < ncols; ++i)
55  {
56  if (debug)
57  {
58  libMesh::out << cols[i] << " ";
59  }
60  _connected_dof_indices.push_back(cols[i]);
61  }
62  if (debug)
63  {
64  libMesh::out << "]\n";
65  }
66  ierr = MatRestoreRow(jac, static_cast<PetscInt>(_var.nodalDofIndex()), &ncols, &cols, PETSC_NULL);
67  CHKERRABORT(_communicator.get(), ierr);
68 #else
69  NodeFaceConstraint::getConnectedDofIndices();
70 #endif
71 }

◆ jacobianSetup()

void GluedContactConstraint::jacobianSetup ( )
virtual

Definition at line 130 of file GluedContactConstraint.C.

131 {
132  if (_component == 0)
133  {
134  if (_updateContactSet)
136 
137  _updateContactSet = true;
138  }
139 }
const unsigned int _component
virtual void updateContactSet(bool beginning_of_step=false)

◆ shouldApply()

bool GluedContactConstraint::shouldApply ( )

Definition at line 171 of file GluedContactConstraint.C.

172 {
173  auto hpit = _penetration_locator._has_penetrated.find(_current_node->id());
174  return hpit != _penetration_locator._has_penetrated.end();
175 }

◆ timestepSetup()

void GluedContactConstraint::timestepSetup ( )
virtual

Definition at line 120 of file GluedContactConstraint.C.

121 {
122  if (_component == 0)
123  {
124  updateContactSet(true);
125  _updateContactSet = false;
126  }
127 }
const unsigned int _component
virtual void updateContactSet(bool beginning_of_step=false)

◆ updateContactSet()

void GluedContactConstraint::updateContactSet ( bool  beginning_of_step = false)
virtual

Definition at line 142 of file GluedContactConstraint.C.

Referenced by jacobianSetup(), and timestepSetup().

143 {
144  auto & has_penetrated = _penetration_locator._has_penetrated;
145 
146  for (auto & pinfo_pair : _penetration_locator._penetration_info)
147  {
148  PenetrationInfo * pinfo = pinfo_pair.second;
149 
150  // Skip this pinfo if there are no DOFs on this node.
151  if (!pinfo || pinfo->_node->n_comp(_sys.number(), _vars[_component]) < 1)
152  continue;
153 
154  const dof_id_type slave_node_num = pinfo_pair.first;
155  auto hpit = has_penetrated.find(slave_node_num);
156 
157  if (beginning_of_step)
158  {
159  pinfo->_starting_elem = pinfo->_elem;
160  pinfo->_starting_side_num = pinfo->_side_num;
161  pinfo->_starting_closest_point_ref = pinfo->_closest_point_ref;
162  }
163 
164  if (pinfo->_distance >= 0)
165  if (hpit == has_penetrated.end())
166  has_penetrated.insert(slave_node_num);
167  }
168 }
std::vector< unsigned int > _vars
const unsigned int _component

Member Data Documentation

◆ _aux_solution

const NumericVector<Number>* GluedContactConstraint::_aux_solution
protected

Definition at line 74 of file GluedContactConstraint.h.

◆ _aux_system

SystemBase& GluedContactConstraint::_aux_system
protected

Definition at line 73 of file GluedContactConstraint.h.

◆ _component

const unsigned int GluedContactConstraint::_component
protected

◆ _formulation

const ContactFormulation GluedContactConstraint::_formulation
protected

Definition at line 57 of file GluedContactConstraint.h.

◆ _friction_coefficient

const Real GluedContactConstraint::_friction_coefficient
protected

Definition at line 60 of file GluedContactConstraint.h.

◆ _model

const ContactModel GluedContactConstraint::_model
protected

Definition at line 56 of file GluedContactConstraint.h.

◆ _nodal_area_var

MooseVariable* GluedContactConstraint::_nodal_area_var
protected

Definition at line 72 of file GluedContactConstraint.h.

◆ _penalty

const Real GluedContactConstraint::_penalty
protected

Definition at line 59 of file GluedContactConstraint.h.

Referenced by computeQpJacobian(), and computeQpResidual().

◆ _residual_copy

NumericVector<Number>& GluedContactConstraint::_residual_copy
protected

Definition at line 64 of file GluedContactConstraint.h.

Referenced by computeQpResidual().

◆ _tension_release

const Real GluedContactConstraint::_tension_release
protected

Definition at line 61 of file GluedContactConstraint.h.

◆ _updateContactSet

bool GluedContactConstraint::_updateContactSet
protected

Definition at line 62 of file GluedContactConstraint.h.

Referenced by jacobianSetup(), and timestepSetup().

◆ _vars

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

◆ _x_var

unsigned int GluedContactConstraint::_x_var
protected

Definition at line 66 of file GluedContactConstraint.h.

◆ _y_var

unsigned int GluedContactConstraint::_y_var
protected

Definition at line 67 of file GluedContactConstraint.h.

◆ _z_var

unsigned int GluedContactConstraint::_z_var
protected

Definition at line 68 of file GluedContactConstraint.h.


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