www.mooseframework.org
LinearNodalConstraint.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 // MOOSE includes
16 #include "LinearNodalConstraint.h"
17 #include "MooseMesh.h"
18 
19 template <>
22 {
24  params.addRequiredParam<std::vector<unsigned int>>("master", "The master node IDs.");
25  params.addParam<std::vector<unsigned int>>("slave_node_ids", "The list of slave node ids");
26  params.addParam<BoundaryName>(
27  "slave_node_set", "NaN", "The boundary ID associated with the slave side");
28  params.addRequiredParam<Real>("penalty", "The penalty used for the boundary term");
29  params.addRequiredParam<std::vector<Real>>(
30  "weights",
31  "The weights associated with the master node ids. Must be of the same size as master nodes");
32  return params;
33 }
34 
36  : NodalConstraint(parameters),
37  _master_node_ids(getParam<std::vector<unsigned int>>("master")),
38  _slave_node_ids(getParam<std::vector<unsigned int>>("slave_node_ids")),
39  _slave_node_set_id(getParam<BoundaryName>("slave_node_set")),
40  _penalty(getParam<Real>("penalty"))
41 {
42  _weights = getParam<std::vector<Real>>("weights");
43 
44  if (_master_node_ids.size() != _weights.size())
45  mooseError("master and weights should be of equal size.");
46 
47  if ((_slave_node_ids.size() == 0) && (_slave_node_set_id == "NaN"))
48  mooseError("Please specify slave_node_ids or slave_node_set.");
49  else if ((_slave_node_ids.size() == 0) && (_slave_node_set_id != "NaN"))
50  {
51  std::vector<dof_id_type> nodelist = _mesh.getNodeList(_mesh.getBoundaryID(_slave_node_set_id));
52  std::vector<dof_id_type>::iterator in;
53 
54  for (in = nodelist.begin(); in != nodelist.end(); ++in)
55  {
56  if (_mesh.nodeRef(*in).processor_id() == _subproblem.processor_id())
57  _connected_nodes.push_back(*in); // defining slave nodes in the base class
58  }
59  }
60  else if ((_slave_node_ids.size() > 0) && (_slave_node_set_id == "NaN"))
61  {
62  for (const auto & dof : _slave_node_ids)
63  if (_mesh.nodeRef(dof).processor_id() == _subproblem.processor_id())
64  _connected_nodes.push_back(dof);
65  }
66 
67  const auto & node_to_elem_map = _mesh.nodeToElemMap();
68 
69  // Add elements connected to master node to Ghosted Elements
70  for (const auto & dof : _master_node_ids)
71  {
72  // defining master nodes in base class
73  _master_node_vector.push_back(dof);
74 
75  auto node_to_elem_pair = node_to_elem_map.find(dof);
76  mooseAssert(node_to_elem_pair != node_to_elem_map.end(), "Missing entry in node to elem map");
77  const std::vector<dof_id_type> & elems = node_to_elem_pair->second;
78 
79  for (const auto & elem_id : elems)
80  _subproblem.addGhostedElem(elem_id);
81  }
82 }
83 
84 Real
86 {
96  unsigned int master_size = _master_node_ids.size();
97 
98  switch (type)
99  {
100  case Moose::Master:
101  return (_u_master[_j] * _weights[_j] - _u_slave[_i] / master_size) * _penalty;
102  case Moose::Slave:
103  return (_u_slave[_i] / master_size - _u_master[_j] * _weights[_j]) * _penalty;
104  }
105  return 0.;
106 }
107 
108 Real
110 {
111  unsigned int master_size = _master_node_ids.size();
112 
113  switch (type)
114  {
115  case Moose::MasterMaster:
116  return _penalty * _weights[_j];
117  case Moose::MasterSlave:
118  return -_penalty / master_size;
119  case Moose::SlaveSlave:
120  return _penalty / master_size;
121  case Moose::SlaveMaster:
122  return -_penalty * _weights[_j];
123  }
124  return 0.;
125 }
virtual Real computeQpResidual(Moose::ConstraintType type) override
Computes the residual for the current slave node.
const VariableValue & _u_master
Holds the current solution at the current quadrature point.
ConstraintType
Definition: MooseTypes.h:198
std::vector< unsigned int > _master_node_ids
unsigned int _j
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int _i
Counter for master and slave nodes.
const VariableValue & _u_slave
Value of the unknown variable this BC is action on.
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...
InputParameters validParams< NodalConstraint >()
BoundaryID getBoundaryID(const BoundaryName &boundary_name) const
Get the associated BoundaryID for the boundary name.
Definition: MooseMesh.C:929
virtual Real computeQpJacobian(Moose::ConstraintJacobianType type) override
Computes the jacobian for the constraint.
std::vector< dof_id_type > _connected_nodes
node IDs connected to the master node (slave nodes)
std::vector< Real > _weights
When the slave node is constrained to move as a linear combination of the master nodes, the coefficients associated with each master node is stored in _weights.
MooseMesh & _mesh
Definition: Constraint.h:81
std::vector< unsigned int > _slave_node_ids
SubProblem & _subproblem
Definition: Constraint.h:74
LinearNodalConstraint(const InputParameters &parameters)
MatType type
std::vector< dof_id_type > _master_node_vector
node IDs of the master node
InputParameters validParams< LinearNodalConstraint >()
virtual void addGhostedElem(dof_id_type elem_id)=0
Will make sure that all dofs connected to elem_id are ghosted to this processor.
ConstraintJacobianType
Definition: MooseTypes.h:204
virtual const Node & nodeRef(const dof_id_type i) const
Definition: MooseMesh.C:398
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...
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
const std::vector< dof_id_type > & getNodeList(boundary_id_type nodeset_id) const
Return a writable reference to a vector of node IDs that belong to nodeset_id.
Definition: MooseMesh.C:2382
const std::map< dof_id_type, std::vector< dof_id_type > > & nodeToElemMap()
If not already created, creates a map from every node to all elements to which they are connected...
Definition: MooseMesh.C:642