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 26 of file GluedContactConstraint.h.

Constructor & Destructor Documentation

GluedContactConstraint::GluedContactConstraint ( const InputParameters &  parameters)

Definition at line 67 of file GluedContactConstraint.C.

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)),
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 }
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
virtual GluedContactConstraint::~GluedContactConstraint ( )
inlinevirtual

Definition at line 30 of file GluedContactConstraint.h.

30 {}

Member Function Documentation

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

Reimplemented from SparsityBasedContactConstraint.

Definition at line 215 of file GluedContactConstraint.C.

Referenced by ~GluedContactConstraint().

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 }
std::vector< unsigned int > _vars
const unsigned int _component
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 240 of file GluedContactConstraint.C.

Referenced by ~GluedContactConstraint().

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
Real GluedContactConstraint::computeQpResidual ( Moose::ConstraintType  type)
virtual

Reimplemented from SparsityBasedContactConstraint.

Definition at line 180 of file GluedContactConstraint.C.

Referenced by ~GluedContactConstraint().

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 }
std::vector< unsigned int > _vars
const unsigned int _component
NumericVector< Number > & _residual_copy
Real GluedContactConstraint::computeQpSlaveValue ( )
virtual

Reimplemented from SparsityBasedContactConstraint.

Definition at line 174 of file GluedContactConstraint.C.

Referenced by ~GluedContactConstraint().

175 {
176  return _u_slave[_qp];
177 }
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 25 of file SparsityBasedContactConstraint.C.

Referenced by SparsityBasedContactConstraint::computeQpJacobian().

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

Definition at line 126 of file GluedContactConstraint.C.

Referenced by ~GluedContactConstraint().

127 {
128  if (_component == 0)
129  {
130  if (_updateContactSet)
132 
133  _updateContactSet = true;
134  }
135 }
const unsigned int _component
virtual void updateContactSet(bool beginning_of_step=false)
bool GluedContactConstraint::shouldApply ( )

Definition at line 167 of file GluedContactConstraint.C.

Referenced by ~GluedContactConstraint().

168 {
169  auto hpit = _penetration_locator._has_penetrated.find(_current_node->id());
170  return hpit != _penetration_locator._has_penetrated.end();
171 }
void GluedContactConstraint::timestepSetup ( )
virtual

Definition at line 116 of file GluedContactConstraint.C.

Referenced by ~GluedContactConstraint().

117 {
118  if (_component == 0)
119  {
120  updateContactSet(true);
121  _updateContactSet = false;
122  }
123 }
const unsigned int _component
virtual void updateContactSet(bool beginning_of_step=false)
void GluedContactConstraint::updateContactSet ( bool  beginning_of_step = false)
virtual

Definition at line 138 of file GluedContactConstraint.C.

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

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 }
std::vector< unsigned int > _vars
const unsigned int _component

Member Data Documentation

const NumericVector<Number>* GluedContactConstraint::_aux_solution
protected

Definition at line 72 of file GluedContactConstraint.h.

SystemBase& GluedContactConstraint::_aux_system
protected

Definition at line 71 of file GluedContactConstraint.h.

const unsigned int GluedContactConstraint::_component
protected
const ContactFormulation GluedContactConstraint::_formulation
protected

Definition at line 55 of file GluedContactConstraint.h.

const Real GluedContactConstraint::_friction_coefficient
protected

Definition at line 58 of file GluedContactConstraint.h.

const ContactModel GluedContactConstraint::_model
protected

Definition at line 54 of file GluedContactConstraint.h.

MooseVariable* GluedContactConstraint::_nodal_area_var
protected

Definition at line 70 of file GluedContactConstraint.h.

const Real GluedContactConstraint::_penalty
protected

Definition at line 57 of file GluedContactConstraint.h.

Referenced by computeQpJacobian(), and computeQpResidual().

NumericVector<Number>& GluedContactConstraint::_residual_copy
protected

Definition at line 62 of file GluedContactConstraint.h.

Referenced by computeQpResidual().

const Real GluedContactConstraint::_tension_release
protected

Definition at line 59 of file GluedContactConstraint.h.

bool GluedContactConstraint::_updateContactSet
protected

Definition at line 60 of file GluedContactConstraint.h.

Referenced by jacobianSetup(), and timestepSetup().

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

Definition at line 64 of file GluedContactConstraint.h.

unsigned int GluedContactConstraint::_y_var
protected

Definition at line 65 of file GluedContactConstraint.h.

unsigned int GluedContactConstraint::_z_var
protected

Definition at line 66 of file GluedContactConstraint.h.


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