www.mooseframework.org
NodeFaceConstraint.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "NodeFaceConstraint.h"
11 
12 // MOOSE includes
13 #include "Assembly.h"
14 #include "MooseEnum.h"
15 #include "MooseMesh.h"
16 #include "MooseVariableFE.h"
17 #include "PenetrationLocator.h"
18 #include "SystemBase.h"
19 
20 #include "libmesh/string_to_enum.h"
21 
22 template <>
25 {
26  MooseEnum orders("FIRST SECOND THIRD FOURTH", "FIRST");
28  params.addRequiredParam<BoundaryName>("slave", "The boundary ID associated with the slave side");
29  params.addRequiredParam<BoundaryName>("master",
30  "The boundary ID associated with the master side");
31  params.addParam<Real>("tangential_tolerance",
32  "Tangential distance to extend edges of contact surfaces");
33  params.addParam<Real>(
34  "normal_smoothing_distance",
35  "Distance from edge in parametric coordinates over which to smooth contact normal");
36  params.addParam<std::string>("normal_smoothing_method",
37  "Method to use to smooth normals (edge_based|nodal_normal_based)");
38  params.addParam<MooseEnum>("order", orders, "The finite element order used for projections");
39 
40  params.addRequiredCoupledVar("master_variable", "The variable on the master side of the domain");
41  params.addRequiredParam<NonlinearVariableName>(
42  "variable", "The name of the variable that this constraint is applied to.");
43 
44  return params;
45 }
46 
48  : Constraint(parameters),
49  // The slave side is at nodes (hence passing 'true'). The neighbor side is the master side and
50  // it is not at nodes (so passing false)
54  _slave(_mesh.getBoundaryID(getParam<BoundaryName>("slave"))),
55  _master(_mesh.getBoundaryID(getParam<BoundaryName>("master"))),
56  _var(_sys.getFieldVariable<Real>(_tid, parameters.get<NonlinearVariableName>("variable"))),
57 
58  _master_q_point(_assembly.qPointsFace()),
59  _master_qrule(_assembly.qRuleFace()),
60 
61  _penetration_locator(
62  getPenetrationLocator(getParam<BoundaryName>("master"),
63  getParam<BoundaryName>("slave"),
64  Utility::string_to_enum<Order>(getParam<MooseEnum>("order")))),
65 
66  _current_node(_var.node()),
67  _current_master(_var.neighbor()),
68  _u_slave(_var.dofValues()),
69  _phi_slave(1), // One entry
70  _test_slave(1), // One entry
71 
72  _master_var(*getVar("master_variable", 0)),
73  _master_var_num(_master_var.number()),
74 
75  _phi_master(_assembly.phiFaceNeighbor(_master_var)),
76  _grad_phi_master(_assembly.gradPhiFaceNeighbor(_master_var)),
77 
78  _test_master(_var.phiFaceNeighbor()),
79  _grad_test_master(_var.gradPhiFaceNeighbor()),
80 
81  _u_master(_master_var.slnNeighbor()),
82  _grad_u_master(_master_var.gradSlnNeighbor()),
83 
84  _dof_map(_sys.dofMap()),
85  _node_to_elem_map(_mesh.nodeToElemMap()),
86 
87  _overwrite_slave_residual(true)
88 {
90 
91  if (parameters.isParamValid("tangential_tolerance"))
92  {
93  _penetration_locator.setTangentialTolerance(getParam<Real>("tangential_tolerance"));
94  }
95  if (parameters.isParamValid("normal_smoothing_distance"))
96  {
97  _penetration_locator.setNormalSmoothingDistance(getParam<Real>("normal_smoothing_distance"));
98  }
99  if (parameters.isParamValid("normal_smoothing_method"))
100  {
102  parameters.get<std::string>("normal_smoothing_method"));
103  }
104  // Put a "1" into test_slave
105  // will always only have one entry that is 1
106  _test_slave[0].push_back(1);
107 }
108 
110 {
112  _test_slave.release();
113 }
114 
115 void
116 NodeFaceConstraint::computeSlaveValue(NumericVector<Number> & current_solution)
117 {
118  const dof_id_type & dof_idx = _var.nodalDofIndex();
119  _qp = 0;
120  current_solution.set(dof_idx, computeQpSlaveValue());
121 }
122 
123 void
125 {
126  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
127  DenseVector<Number> & neighbor_re = _assembly.residualBlockNeighbor(_master_var.number());
128 
129  _qp = 0;
130 
131  for (_i = 0; _i < _test_master.size(); _i++)
132  neighbor_re(_i) += computeQpResidual(Moose::Master);
133 
134  _i = 0;
136 }
137 
138 void
140 {
142 
143  // DenseMatrix<Number> & Kee = _assembly.jacobianBlock(_var.number(), _var.number());
144  DenseMatrix<Number> & Ken =
146 
147  // DenseMatrix<Number> & Kne = _assembly.jacobianBlockNeighbor(Moose::NeighborElement,
148  // _var.number(), _var.number());
149  DenseMatrix<Number> & Knn =
151 
152  _Kee.resize(_test_slave.size(), _connected_dof_indices.size());
153  _Kne.resize(_test_master.size(), _connected_dof_indices.size());
154 
156 
157  _qp = 0;
158 
159  // Fill up _phi_slave so that it is 1 when j corresponds to this dof and 0 for every other dof
160  // This corresponds to evaluating all of the connected shape functions at _this_ node
161  for (unsigned int j = 0; j < _connected_dof_indices.size(); j++)
162  {
163  _phi_slave[j].resize(1);
164 
166  _phi_slave[j][_qp] = 1.0;
167  else
168  _phi_slave[j][_qp] = 0.0;
169  }
170 
171  for (_i = 0; _i < _test_slave.size(); _i++)
172  // Loop over the connected dof indices so we can get all the jacobian contributions
173  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
175 
176  if (Ken.m() && Ken.n())
177  for (_i = 0; _i < _test_slave.size(); _i++)
178  for (_j = 0; _j < _phi_master.size(); _j++)
180 
181  for (_i = 0; _i < _test_master.size(); _i++)
182  // Loop over the connected dof indices so we can get all the jacobian contributions
183  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
185 
186  if (Knn.m() && Knn.n())
187  for (_i = 0; _i < _test_master.size(); _i++)
188  for (_j = 0; _j < _phi_master.size(); _j++)
190 }
191 
192 void
194 {
196 
197  _Kee.resize(_test_slave.size(), _connected_dof_indices.size());
198  _Kne.resize(_test_master.size(), _connected_dof_indices.size());
199 
200  DenseMatrix<Number> & Ken =
202  DenseMatrix<Number> & Knn =
204 
206 
207  _qp = 0;
208 
209  // Fill up _phi_slave so that it is 1 when j corresponds to this dof and 0 for every other dof
210  // This corresponds to evaluating all of the connected shape functions at _this_ node
211  for (unsigned int j = 0; j < _connected_dof_indices.size(); j++)
212  {
213  _phi_slave[j].resize(1);
214 
216  _phi_slave[j][_qp] = 1.0;
217  else
218  _phi_slave[j][_qp] = 0.0;
219  }
220 
221  for (_i = 0; _i < _test_slave.size(); _i++)
222  // Loop over the connected dof indices so we can get all the jacobian contributions
223  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
225 
226  for (_i = 0; _i < _test_slave.size(); _i++)
227  for (_j = 0; _j < _phi_master.size(); _j++)
229 
230  if (_Kne.m() && _Kne.n())
231  for (_i = 0; _i < _test_master.size(); _i++)
232  // Loop over the connected dof indices so we can get all the jacobian contributions
233  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
235 
236  for (_i = 0; _i < _test_master.size(); _i++)
237  for (_j = 0; _j < _phi_master.size(); _j++)
239 }
240 
241 void
243 {
244  MooseVariableFEBase & var = _sys.getVariable(0, var_num);
245 
246  _connected_dof_indices.clear();
247  std::set<dof_id_type> unique_dof_indices;
248 
249  auto node_to_elem_pair = _node_to_elem_map.find(_current_node->id());
250  mooseAssert(node_to_elem_pair != _node_to_elem_map.end(), "Missing entry in node to elem map");
251  const std::vector<dof_id_type> & elems = node_to_elem_pair->second;
252 
253  // Get the dof indices from each elem connected to the node
254  for (const auto & cur_elem : elems)
255  {
256  std::vector<dof_id_type> dof_indices;
257 
258  var.getDofIndices(_mesh.elemPtr(cur_elem), dof_indices);
259 
260  for (const auto & dof : dof_indices)
261  unique_dof_indices.insert(dof);
262  }
263 
264  for (const auto & dof : unique_dof_indices)
265  _connected_dof_indices.push_back(dof);
266 }
267 
268 bool
270 {
272 }
VarFieldType
Definition: MooseTypes.h:488
virtual void computeOffDiagJacobian(unsigned int jvar)
Computes d-residual / d-jvar...
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
DenseVector< Number > & residualBlockNeighbor(unsigned int var_num, TagID tag_id=0)
Definition: Assembly.h:725
virtual Elem * elemPtr(const dof_id_type i)
Definition: MooseMesh.C:2267
void setNormalSmoothingDistance(Real normal_smoothing_distance)
MooseVariable & _var
unsigned int number() const
Get variable number coming from libMesh.
void resize(unsigned int size)
Change the number of elements the array can store.
Definition: MooseArray.h:218
virtual Real computeQpResidual(Moose::ConstraintType type)=0
This is the virtual that derived classes should override for computing the residual on neighboring el...
InputParameters validParams< NodeFaceConstraint >()
Base class for all Constraint types.
Definition: Constraint.h:39
DenseMatrix< Number > _Kne
const std::map< dof_id_type, std::vector< dof_id_type > > & _node_to_elem_map
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int _i
Definition: Constraint.h:75
MooseVariable & _master_var
Master side variable.
void setTangentialTolerance(Real tangential_tolerance)
const VariableTestValue & _test_master
Side test function.
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
virtual void computeResidual()
Computes the residual Nodal residual.
void addMooseVariableDependency(MooseVariableFEBase *var)
Call this function to add the passed in MooseVariableFEBase as a variable that this object depends on...
const Node *const & _current_node
current node being processed
DenseMatrix< Number > _Kee
MooseVariableFEBase & getVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:105
unsigned int _j
Definition: Constraint.h:75
NodeFaceConstraint(const InputParameters &parameters)
Enhances MooseVariableInterface interface provide values from neighbor elements.
virtual bool overwriteSlaveResidual()
Whether or not the slave&#39;s residual should be overwritten.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:259
VariableTestValue _test_slave
Shape function on the slave side. This will always only have one entry and that entry will always be ...
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseObject.h:65
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
MooseMesh & _mesh
Definition: Constraint.h:73
virtual Real computeQpSlaveValue()=0
Compute the value the slave node should have at the beginning of a timestep.
VariablePhiValue _phi_slave
Shape function on the slave side. This will always.
virtual void getConnectedDofIndices(unsigned int var_num)
Gets the indices for all dofs connected to the constraint.
virtual void getDofIndices(const Elem *elem, std::vector< dof_id_type > &dof_indices) const =0
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
virtual void computeJacobian()
Computes the jacobian for the current element.
Assembly & _assembly
Definition: Constraint.h:72
void addRequiredCoupledVar(const std::string &name, const std::string &doc_string)
This method adds a coupled variable name pair.
const VariablePhiValue & _phi_master
Side shape function.
const dof_id_type & nodalDofIndex() const override
void release()
Manually deallocates the data pointer.
Definition: MooseArray.h:61
SystemBase & _sys
Definition: Constraint.h:68
std::vector< dof_id_type > _connected_dof_indices
DenseMatrix< Number > & jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar, TagID tag=0)
Definition: Assembly.C:2033
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
bool _overwrite_slave_residual
Whether or not the slave&#39;s residual should be overwritten.
Definition: Moose.h:112
void setNormalSmoothingMethod(std::string nsmString)
DenseVector< Number > & residualBlock(unsigned int var_num, TagID tag_id=0)
Definition: Assembly.h:720
PenetrationLocator & _penetration_locator
InputParameters validParams< Constraint >()
Definition: Constraint.C:16
virtual Real computeQpJacobian(Moose::ConstraintJacobianType type)=0
This is the virtual that derived classes should override for computing the Jacobian on neighboring el...
unsigned int _qp
Definition: Constraint.h:76
void computeSlaveValue(NumericVector< Number > &current_solution)
Compute the value the slave node should have at the beginning of a timestep.
virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType, unsigned int)
This is the virtual that derived classes should override for computing the off-diag Jacobian...
bool isParamValid(const std::string &name) const
This method returns parameters that have been initialized in one fashion or another, i.e.