www.mooseframework.org
NodeElemConstraint.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 "NodeElemConstraint.h"
11 
12 // MOOSE includes
13 #include "Assembly.h"
14 #include "MooseEnum.h"
15 #include "MooseMesh.h"
16 #include "MooseVariableFE.h"
17 #include "SystemBase.h"
18 
19 #include "libmesh/string_to_enum.h"
20 
21 template <>
24 {
26  params.addRequiredParam<SubdomainName>("slave", "slave block id");
27  params.addRequiredParam<SubdomainName>("master", "master block id");
28  params.addRequiredCoupledVar("master_variable", "The variable on the master side of the domain");
29  params.addRequiredParam<NonlinearVariableName>(
30  "variable", "The name of the variable that this constraint is applied to.");
31 
32  return params;
33 }
34 
36  : Constraint(parameters),
37  // The slave side is at nodes (hence passing 'true'). The neighbor side is the master side and
38  // it is not at nodes (so passing false)
42 
43  _slave(_mesh.getSubdomainID(getParam<SubdomainName>("slave"))),
44  _master(_mesh.getSubdomainID(getParam<SubdomainName>("master"))),
45  _var(_sys.getFieldVariable<Real>(_tid, parameters.get<NonlinearVariableName>("variable"))),
46 
47  _master_q_point(_assembly.qPoints()),
48  _master_qrule(_assembly.qRule()),
49 
50  _current_node(_var.node()),
51  _current_elem(_var.neighbor()),
52 
53  _u_slave(_var.dofValues()),
54  _u_slave_old(_var.dofValuesOld()),
55  _phi_slave(1),
56  _test_slave(1), // One entry
57 
58  _master_var(*getVar("master_variable", 0)),
59  _master_var_num(_master_var.number()),
60 
61  _phi_master(_assembly.phiNeighbor(_master_var)),
62  _grad_phi_master(_assembly.gradPhiNeighbor(_master_var)),
63 
64  _test_master(_var.phiNeighbor()),
65  _grad_test_master(_var.gradPhiNeighbor()),
66 
67  _u_master(_master_var.slnNeighbor()),
68  _u_master_old(_master_var.slnOldNeighbor()),
69  _grad_u_master(_master_var.gradSlnNeighbor()),
70 
71  _dof_map(_sys.dofMap()),
72  _node_to_elem_map(_mesh.nodeToElemMap()),
73 
74  _overwrite_slave_residual(false)
75 {
76  _mesh.errorIfDistributedMesh("NodeElemConstraint");
77 
79  // Put a "1" into test_slave
80  // will always only have one entry that is 1
81  _test_slave[0].push_back(1);
82 }
83 
85 {
87  _test_slave.release();
88 }
89 
90 void
91 NodeElemConstraint::computeSlaveValue(NumericVector<Number> & current_solution)
92 {
93  const dof_id_type & dof_idx = _var.nodalDofIndex();
94  _qp = 0;
95  current_solution.set(dof_idx, computeQpSlaveValue());
96 }
97 
98 void
100 {
101  DenseVector<Number> & slave_re = _assembly.residualBlock(_var.number());
102  DenseVector<Number> & master_re = _assembly.residualBlockNeighbor(_master_var.number());
103 
104  _qp = 0;
105 
106  for (_i = 0; _i < _test_master.size(); _i++)
107  master_re(_i) += computeQpResidual(Moose::Master);
108 
109  for (_i = 0; _i < _test_slave.size(); _i++)
110  slave_re(_i) += computeQpResidual(Moose::Slave);
111 }
112 
113 void
115 {
117 
118  DenseMatrix<Number> & Ken =
120 
121  DenseMatrix<Number> & Knn =
123 
124  _Kee.resize(_test_slave.size(), _connected_dof_indices.size());
125  _Kne.resize(_test_master.size(), _connected_dof_indices.size());
126 
128 
129  _qp = 0;
130 
131  // Fill up _phi_slave so that it is 1 when j corresponds to this dof and 0 for every other dof
132  // This corresponds to evaluating all of the connected shape functions at _this_ node
133  for (unsigned int j = 0; j < _connected_dof_indices.size(); j++)
134  {
135  _phi_slave[j].resize(1);
136 
138  _phi_slave[j][_qp] = 1.0;
139  else
140  _phi_slave[j][_qp] = 0.0;
141  }
142 
143  for (_i = 0; _i < _test_slave.size(); _i++)
144  // Loop over the connected dof indices so we can get all the jacobian contributions
145  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
147 
148  if (Ken.m() && Ken.n())
149  for (_i = 0; _i < _test_slave.size(); _i++)
150  for (_j = 0; _j < _phi_master.size(); _j++)
152 
153  for (_i = 0; _i < _test_master.size(); _i++)
154  // Loop over the connected dof indices so we can get all the jacobian contributions
155  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
157 
158  if (Knn.m() && Knn.n())
159  for (_i = 0; _i < _test_master.size(); _i++)
160  for (_j = 0; _j < _phi_master.size(); _j++)
162 }
163 
164 void
166 {
168 
169  _Kee.resize(_test_slave.size(), _connected_dof_indices.size());
170  _Kne.resize(_test_master.size(), _connected_dof_indices.size());
171 
172  DenseMatrix<Number> & Ken =
174  DenseMatrix<Number> & Knn =
176 
178 
179  _qp = 0;
180 
181  // Fill up _phi_slave so that it is 1 when j corresponds to this dof and 0 for every other dof
182  // This corresponds to evaluating all of the connected shape functions at _this_ node
183  for (unsigned int j = 0; j < _connected_dof_indices.size(); j++)
184  {
185  _phi_slave[j].resize(1);
186 
188  _phi_slave[j][_qp] = 1.0;
189  else
190  _phi_slave[j][_qp] = 0.0;
191  }
192 
193  for (_i = 0; _i < _test_slave.size(); _i++)
194  // Loop over the connected dof indices so we can get all the jacobian contributions
195  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
197 
198  for (_i = 0; _i < _test_slave.size(); _i++)
199  for (_j = 0; _j < _phi_master.size(); _j++)
201 
202  if (_Kne.m() && _Kne.n())
203  for (_i = 0; _i < _test_master.size(); _i++)
204  // Loop over the connected dof indices so we can get all the jacobian contributions
205  for (_j = 0; _j < _connected_dof_indices.size(); _j++)
207 
208  for (_i = 0; _i < _test_master.size(); _i++)
209  for (_j = 0; _j < _phi_master.size(); _j++)
211 }
212 
213 void
215 {
216  MooseVariableFEBase & var = _sys.getVariable(0, var_num);
217 
218  _connected_dof_indices.clear();
219  std::set<dof_id_type> unique_dof_indices;
220 
221  auto node_to_elem_pair = _node_to_elem_map.find(_current_node->id());
222  mooseAssert(node_to_elem_pair != _node_to_elem_map.end(), "Missing entry in node to elem map");
223  const std::vector<dof_id_type> & elems = node_to_elem_pair->second;
224 
225  // Get the dof indices from each elem connected to the node
226  for (const auto & cur_elem : elems)
227  {
228  std::vector<dof_id_type> dof_indices;
229 
230  var.getDofIndices(_mesh.elemPtr(cur_elem), dof_indices);
231 
232  for (const auto & dof : dof_indices)
233  unique_dof_indices.insert(dof);
234  }
235 
236  for (const auto & dof : unique_dof_indices)
237  _connected_dof_indices.push_back(dof);
238 }
239 
240 bool
242 {
244 }
VarFieldType
Definition: MooseTypes.h:488
virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType, unsigned int)
This is the virtual that derived classes should override for computing the off-diag Jacobian...
NodeElemConstraint(const InputParameters &parameters)
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
DenseMatrix< Number > _Kee
stiffness matrix holding slave-slave jacobian
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.
const std::map< dof_id_type, std::vector< dof_id_type > > & _node_to_elem_map
Base class for all Constraint types.
Definition: Constraint.h:39
virtual void computeResidual()
Computes the residual Nodal residual.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual Real computeQpJacobian(Moose::ConstraintJacobianType type)=0
This is the virtual that derived classes should override for computing the Jacobian.
unsigned int _i
Definition: Constraint.h:75
virtual Real computeQpSlaveValue()=0
Compute the value the slave node should have at the beginning of a timestep.
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...
const VariableTestValue & _test_master
Side test function.
void addMooseVariableDependency(MooseVariableFEBase *var)
Call this function to add the passed in MooseVariableFEBase as a variable that this object depends on...
virtual void computeOffDiagJacobian(unsigned int jvar)
Computes d-residual / d-jvar...
MooseVariableFEBase & getVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:105
MooseVariable & _var
InputParameters validParams< NodeElemConstraint >()
virtual void getConnectedDofIndices(unsigned int var_num)
Gets the indices for all dofs connected to the constraint.
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
Definition: MooseMesh.C:2685
unsigned int _j
Definition: Constraint.h:75
Enhances MooseVariableInterface interface provide values from neighbor elements.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:259
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
MooseMesh & _mesh
Definition: Constraint.h:73
const VariablePhiValue & _phi_master
Side shape function.
virtual void getDofIndices(const Elem *elem, std::vector< dof_id_type > &dof_indices) const =0
VariablePhiValue _phi_slave
Shape function on the slave side.
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.
VariableTestValue _test_slave
Shape function on the slave side. This will always only have one entry and that entry will always be ...
const dof_id_type & nodalDofIndex() const override
void release()
Manually deallocates the data pointer.
Definition: MooseArray.h:61
const Node *const & _current_node
current node being processed
SystemBase & _sys
Definition: Constraint.h:68
std::vector< dof_id_type > _connected_dof_indices
dofs connected to the slave node
virtual bool overwriteSlaveResidual()
Whether or not the slave&#39;s residual should be overwritten.
void computeSlaveValue(NumericVector< Number > &current_solution)
Compute the value the slave node should have at the beginning of a timestep.
MooseVariable & _master_var
Master side variable.
DenseMatrix< Number > & jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar, TagID tag=0)
Definition: Assembly.C:2033
Definition: Moose.h:112
DenseVector< Number > & residualBlock(unsigned int var_num, TagID tag_id=0)
Definition: Assembly.h:720
virtual void computeJacobian()
Computes the jacobian for the current element.
InputParameters validParams< Constraint >()
Definition: Constraint.C:16
DenseMatrix< Number > _Kne
stiffness matrix holding master-slave jacobian
unsigned int _qp
Definition: Constraint.h:76
bool _overwrite_slave_residual
Whether or not the slave&#39;s residual should be overwritten.