www.mooseframework.org
NodeFaceConstraint.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #include "NodeFaceConstraint.h"
16 
17 // MOOSE includes
18 #include "Assembly.h"
19 #include "MooseEnum.h"
20 #include "MooseMesh.h"
21 #include "MooseVariable.h"
22 #include "PenetrationLocator.h"
23 #include "SystemBase.h"
24 
25 #include "libmesh/string_to_enum.h"
26 
27 template <>
30 {
31  MooseEnum orders("FIRST SECOND THIRD FOURTH", "FIRST");
33  params.addRequiredParam<BoundaryName>("slave", "The boundary ID associated with the slave side");
34  params.addRequiredParam<BoundaryName>("master",
35  "The boundary ID associated with the master side");
36  params.addParam<Real>("tangential_tolerance",
37  "Tangential distance to extend edges of contact surfaces");
38  params.addParam<Real>(
39  "normal_smoothing_distance",
40  "Distance from edge in parametric coordinates over which to smooth contact normal");
41  params.addParam<std::string>("normal_smoothing_method",
42  "Method to use to smooth normals (edge_based|nodal_normal_based)");
43  params.addParam<MooseEnum>("order", orders, "The finite element order used for projections");
44 
45  params.addRequiredCoupledVar("master_variable", "The variable on the master side of the domain");
46 
47  return params;
48 }
49 
51  : Constraint(parameters),
52  // The slave side is at nodes (hence passing 'true'). The neighbor side is the master side and
53  // it is not at nodes (so passing false)
55  _slave(_mesh.getBoundaryID(getParam<BoundaryName>("slave"))),
56  _master(_mesh.getBoundaryID(getParam<BoundaryName>("master"))),
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.nodalSln()),
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()),
76  _grad_phi_master(_assembly.gradPhiFaceNeighbor()),
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 {
89  if (parameters.isParamValid("tangential_tolerance"))
90  {
91  _penetration_locator.setTangentialTolerance(getParam<Real>("tangential_tolerance"));
92  }
93  if (parameters.isParamValid("normal_smoothing_distance"))
94  {
95  _penetration_locator.setNormalSmoothingDistance(getParam<Real>("normal_smoothing_distance"));
96  }
97  if (parameters.isParamValid("normal_smoothing_method"))
98  {
100  parameters.get<std::string>("normal_smoothing_method"));
101  }
102  // Put a "1" into test_slave
103  // will always only have one entry that is 1
104  _test_slave[0].push_back(1);
105 }
106 
108 {
111 }
112 
113 void
114 NodeFaceConstraint::computeSlaveValue(NumericVector<Number> & current_solution)
115 {
116  dof_id_type & dof_idx = _var.nodalDofIndex();
117  _qp = 0;
118  current_solution.set(dof_idx, computeQpSlaveValue());
119 }
120 
121 void
123 {
124  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
125  DenseVector<Number> & neighbor_re = _assembly.residualBlockNeighbor(_master_var.number());
126 
127  _qp = 0;
128 
129  for (_i = 0; _i < _test_master.size(); _i++)
130  neighbor_re(_i) += computeQpResidual(Moose::Master);
131 
132  _i = 0;
134 }
135 
136 void
138 {
140 
141  // DenseMatrix<Number> & Kee = _assembly.jacobianBlock(_var.number(), _var.number());
142  DenseMatrix<Number> & Ken =
144 
145  // DenseMatrix<Number> & Kne = _assembly.jacobianBlockNeighbor(Moose::NeighborElement,
146  // _var.number(), _var.number());
147  DenseMatrix<Number> & Knn =
149 
150  _Kee.resize(_test_slave.size(), _connected_dof_indices.size());
151  _Kne.resize(_test_master.size(), _connected_dof_indices.size());
152 
154 
155  _qp = 0;
156 
157  // Fill up _phi_slave so that it is 1 when j corresponds to this dof and 0 for every other dof
158  // This corresponds to evaluating all of the connected shape functions at _this_ node
159  for (unsigned int j = 0; j < _connected_dof_indices.size(); j++)
160  {
161  _phi_slave[j].resize(1);
162 
164  _phi_slave[j][_qp] = 1.0;
165  else
166  _phi_slave[j][_qp] = 0.0;
167  }
168 
169  for (_i = 0; _i < _test_slave.size(); _i++)
170  // Loop over the connected dof indices so we can get all the jacobian contributions
171  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
173 
174  if (Ken.m() && Ken.n())
175  for (_i = 0; _i < _test_slave.size(); _i++)
176  for (_j = 0; _j < _phi_master.size(); _j++)
178 
179  for (_i = 0; _i < _test_master.size(); _i++)
180  // Loop over the connected dof indices so we can get all the jacobian contributions
181  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
183 
184  if (Knn.m() && Knn.n())
185  for (_i = 0; _i < _test_master.size(); _i++)
186  for (_j = 0; _j < _phi_master.size(); _j++)
188 }
189 
190 void
192 {
194 
195  _Kee.resize(_test_slave.size(), _connected_dof_indices.size());
196  _Kne.resize(_test_master.size(), _connected_dof_indices.size());
197 
198  DenseMatrix<Number> & Ken =
200  DenseMatrix<Number> & Knn =
202 
204 
205  _qp = 0;
206 
207  // Fill up _phi_slave so that it is 1 when j corresponds to this dof and 0 for every other dof
208  // This corresponds to evaluating all of the connected shape functions at _this_ node
209  for (unsigned int j = 0; j < _connected_dof_indices.size(); j++)
210  {
211  _phi_slave[j].resize(1);
212 
214  _phi_slave[j][_qp] = 1.0;
215  else
216  _phi_slave[j][_qp] = 0.0;
217  }
218 
219  for (_i = 0; _i < _test_slave.size(); _i++)
220  // Loop over the connected dof indices so we can get all the jacobian contributions
221  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
223 
224  for (_i = 0; _i < _test_slave.size(); _i++)
225  for (_j = 0; _j < _phi_master.size(); _j++)
227 
228  if (_Kne.m() && _Kne.n())
229  for (_i = 0; _i < _test_master.size(); _i++)
230  // Loop over the connected dof indices so we can get all the jacobian contributions
231  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
233 
234  for (_i = 0; _i < _test_master.size(); _i++)
235  for (_j = 0; _j < _phi_master.size(); _j++)
237 }
238 
239 void
241 {
242  MooseVariable & var = _sys.getVariable(0, var_num);
243 
244  _connected_dof_indices.clear();
245  std::set<dof_id_type> unique_dof_indices;
246 
247  auto node_to_elem_pair = _node_to_elem_map.find(_current_node->id());
248  mooseAssert(node_to_elem_pair != _node_to_elem_map.end(), "Missing entry in node to elem map");
249  const std::vector<dof_id_type> & elems = node_to_elem_pair->second;
250 
251  // Get the dof indices from each elem connected to the node
252  for (const auto & cur_elem : elems)
253  {
254  std::vector<dof_id_type> dof_indices;
255 
256  var.getDofIndices(_mesh.elemPtr(cur_elem), dof_indices);
257 
258  for (const auto & dof : dof_indices)
259  unique_dof_indices.insert(dof);
260  }
261 
262  for (const auto & dof : unique_dof_indices)
263  _connected_dof_indices.push_back(dof);
264 }
265 
266 bool
268 {
270 }
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...
Class for stuff related to variables.
Definition: MooseVariable.h:43
virtual Elem * elemPtr(const dof_id_type i)
Definition: MooseMesh.C:2115
void setNormalSmoothingDistance(Real normal_smoothing_distance)
void resize(const unsigned int size)
Change the number of elements the array can store.
Definition: MooseArray.h:205
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:42
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:83
MooseVariable & _master_var
Master side variable.
void setTangentialTolerance(Real tangential_tolerance)
const VariableTestValue & _test_master
Side test function.
DenseVector< Number > & residualBlock(unsigned int var_num, Moose::KernelType type=Moose::KT_NONTIME)
Definition: Assembly.h:504
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.
DenseMatrix< Number > _Kee
dof_id_type & nodalDofIndex()
unsigned int _j
Definition: Constraint.h:83
NodeFaceConstraint(const InputParameters &parameters)
virtual bool overwriteSlaveResidual()
Whether or not the slave&#39;s residual should be overwritten.
bool isParamValid(const std::string &name) const
This method returns parameters that have been initialized in one fashion or another, i.e.
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:103
VariableTestValue _test_slave
Shape function on the slave side. This will always only have one entry and that entry will always be ...
MooseMesh & _mesh
Definition: Constraint.h:81
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.
const Node *& _current_node
current node being processed
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:37
virtual void computeJacobian()
Computes the jacobian for the current element.
Assembly & _assembly
Definition: Constraint.h:79
void addRequiredCoupledVar(const std::string &name, const std::string &doc_string)
This method adds a coupled variable name pair.
unsigned int number() const
Get variable number coming from libMesh.
const VariablePhiValue & _phi_master
Side shape function.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
void release()
Manually deallocates the data pointer.
Definition: MooseArray.h:56
SystemBase & _sys
Definition: Constraint.h:75
std::vector< dof_id_type > _connected_dof_indices
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.
void getDofIndices(const Elem *elem, std::vector< dof_id_type > &dof_indices)
Get dof indices for the variable.
void setNormalSmoothingMethod(std::string nsmString)
PenetrationLocator & _penetration_locator
InputParameters validParams< Constraint >()
Definition: Constraint.C:20
DenseVector< Number > & residualBlockNeighbor(unsigned int var_num, Moose::KernelType type=Moose::KT_NONTIME)
Definition: Assembly.h:509
DenseMatrix< Number > & jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:901
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:84
void computeSlaveValue(NumericVector< Number > &current_solution)
Compute the value the slave node should have at the beginning of a timestep.
MooseVariable & _var
Definition: Constraint.h:80
virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType, unsigned int)
This is the virtual that derived classes should override for computing the off-diag Jacobian...