www.mooseframework.org
NodeFaceConstraint.h
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 #ifndef NODEFACECONSTRAINT_H
16 #define NODEFACECONSTRAINT_H
17 
18 // MOOSE includes
19 #include "Constraint.h"
21 
22 // Forward Declarations
23 class NodeFaceConstraint;
24 
25 // libMesh forward declarations
26 namespace libMesh
27 {
28 template <typename T>
29 class SparseMatrix;
30 }
31 
32 template <>
34 
46 {
47 public:
48  NodeFaceConstraint(const InputParameters & parameters);
49  virtual ~NodeFaceConstraint();
50 
54  void computeSlaveValue(NumericVector<Number> & current_solution);
55 
59  virtual void computeResidual();
60 
64  virtual void computeJacobian();
65 
69  virtual void computeOffDiagJacobian(unsigned int jvar);
70 
74  virtual void getConnectedDofIndices(unsigned int var_num);
75 
81  virtual bool shouldApply() { return true; }
82 
89  virtual bool overwriteSlaveResidual();
90 
97  virtual bool overwriteSlaveJacobian() { return overwriteSlaveResidual(); };
98 
102  virtual MooseVariable & masterVariable() { return _master_var; }
103 
104  // TODO: Make this protected or add an accessor
105  // Do the same for all the other public members
106  SparseMatrix<Number> * _jacobian;
107 
108 protected:
112  virtual Real computeQpSlaveValue() = 0;
113 
118  virtual Real computeQpResidual(Moose::ConstraintType type) = 0;
119 
124  virtual Real computeQpJacobian(Moose::ConstraintJacobianType type) = 0;
125 
130  unsigned int /*jvar*/)
131  {
132  return 0;
133  }
134 
136  virtual const VariableValue & coupledSlaveValue(const std::string & var_name,
137  unsigned int comp = 0)
138  {
139  return coupledValue(var_name, comp);
140  }
141  virtual const VariableValue & coupledSlaveValueOld(const std::string & var_name,
142  unsigned int comp = 0)
143  {
144  return coupledValueOld(var_name, comp);
145  }
146  virtual const VariableValue & coupledSlaveValueOlder(const std::string & var_name,
147  unsigned int comp = 0)
148  {
149  return coupledValueOlder(var_name, comp);
150  }
151 
152  virtual const VariableGradient & coupledSlaveGradient(const std::string & var_name,
153  unsigned int comp = 0)
154  {
155  return coupledGradient(var_name, comp);
156  }
157  virtual const VariableGradient & coupledSlaveGradientOld(const std::string & var_name,
158  unsigned int comp = 0)
159  {
160  return coupledGradientOld(var_name, comp);
161  }
162  virtual const VariableGradient & coupledSlaveGradientOlder(const std::string & var_name,
163  unsigned int comp = 0)
164  {
165  return coupledGradientOlder(var_name, comp);
166  }
167 
168  virtual const VariableSecond & coupledSlaveSecond(const std::string & var_name,
169  unsigned int comp = 0)
170  {
171  return coupledSecond(var_name, comp);
172  }
173 
174  virtual const VariableValue & coupledMasterValue(const std::string & var_name,
175  unsigned int comp = 0)
176  {
177  return coupledNeighborValue(var_name, comp);
178  }
179  virtual const VariableValue & coupledMasterValueOld(const std::string & var_name,
180  unsigned int comp = 0)
181  {
182  return coupledNeighborValueOld(var_name, comp);
183  }
184  virtual const VariableValue & coupledMasterValueOlder(const std::string & var_name,
185  unsigned int comp = 0)
186  {
187  return coupledNeighborValueOlder(var_name, comp);
188  }
189 
190  virtual const VariableGradient & coupledMasterGradient(const std::string & var_name,
191  unsigned int comp = 0)
192  {
193  return coupledNeighborGradient(var_name, comp);
194  }
195  virtual const VariableGradient & coupledMasterGradientOld(const std::string & var_name,
196  unsigned int comp = 0)
197  {
198  return coupledNeighborGradientOld(var_name, comp);
199  }
200  virtual const VariableGradient & coupledMasterGradientOlder(const std::string & var_name,
201  unsigned int comp = 0)
202  {
203  return coupledNeighborGradientOlder(var_name, comp);
204  }
205 
206  virtual const VariableSecond & coupledMasterSecond(const std::string & var_name,
207  unsigned int comp = 0)
208  {
209  return coupledNeighborSecond(var_name, comp);
210  }
211 
213  unsigned int _slave;
215  unsigned int _master;
216 
218  QBase *& _master_qrule;
219 
220 public:
222 
223 protected:
225  const Node *& _current_node;
226  const Elem *& _current_master;
227 
234 
237 
239  unsigned int _master_var_num;
240 
245 
250 
255 
257  const DofMap & _dof_map;
258 
259  const std::map<dof_id_type, std::vector<dof_id_type>> & _node_to_elem_map;
260 
268 
269 public:
270  std::vector<dof_id_type> _connected_dof_indices;
271 
272  DenseMatrix<Number> _Kne;
273  DenseMatrix<Number> _Kee;
274 };
275 
276 #endif
const VariableGradient & _grad_u_master
Holds the current solution gradient at the current quadrature point.
virtual const VariableGradient & coupledSlaveGradientOlder(const std::string &var_name, unsigned int comp=0)
virtual const VariableGradient & coupledSlaveGradient(const std::string &var_name, unsigned int comp=0)
ConstraintType
Definition: MooseTypes.h:198
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 bool overwriteSlaveJacobian()
Whether or not the slave&#39;s Jacobian row should be overwritten.
unsigned int _slave
Boundary ID for the slave surface.
virtual const VariableValue & coupledSlaveValue(const std::string &var_name, unsigned int comp=0)
coupling interface:
virtual const VariableSecond & coupledSlaveSecond(const std::string &var_name, unsigned int comp=0)
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
virtual const VariableGradient & coupledMasterGradient(const std::string &var_name, unsigned int comp=0)
virtual const VariableGradient & coupledMasterGradientOlder(const std::string &var_name, unsigned int comp=0)
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
const DofMap & _dof_map
DOF map.
const VariableValue & _u_master
Holds the current solution at the current quadrature point.
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
MooseVariable & _master_var
Master side variable.
const VariableTestValue & _test_master
Side test function.
SparseMatrix< Number > * _jacobian
DenseMatrix< Number > _Kee
const VariableTestGradient & _grad_test_master
Gradient of side shape function.
VariableTestValue _test_slave
Shape function on the slave side. This will always only have one entry and that entry will always be ...
virtual const VariableGradient & coupledSlaveGradientOld(const std::string &var_name, unsigned int comp=0)
const VariablePhiGradient & _grad_phi_master
Gradient of side shape function.
InputParameters validParams< NodeFaceConstraint >()
virtual const VariableGradient & coupledMasterGradientOld(const std::string &var_name, unsigned int comp=0)
A NodeFaceConstraint is used when you need to create constraints between two surfaces in a mesh...
VariablePhiValue _phi_slave
Shape function on the slave side. This will always.
const VariableValue & _u_slave
Value of the unknown variable this BC is action on.
const Node *& _current_node
current node being processed
virtual const VariableValue & coupledMasterValue(const std::string &var_name, unsigned int comp=0)
unsigned int _master_var_num
Number for the master variable.
virtual MooseVariable & masterVariable()
The variable on the Master side of the domain.
const MooseArray< Point > & _master_q_point
virtual const VariableValue & coupledMasterValueOlder(const std::string &var_name, unsigned int comp=0)
MatType type
const VariablePhiValue & _phi_master
Side shape function.
virtual const VariableValue & coupledMasterValueOld(const std::string &var_name, unsigned int comp=0)
virtual const VariableValue & coupledSlaveValueOlder(const std::string &var_name, unsigned int comp=0)
virtual bool shouldApply()
Whether or not this constraint should be applied.
ConstraintJacobianType
Definition: MooseTypes.h:204
std::vector< dof_id_type > _connected_dof_indices
bool _overwrite_slave_residual
Whether or not the slave&#39;s residual should be overwritten.
unsigned int _master
Boundary ID for the master surface.
virtual const VariableSecond & coupledMasterSecond(const std::string &var_name, unsigned int comp=0)
PenetrationLocator & _penetration_locator
virtual const VariableValue & coupledSlaveValueOld(const std::string &var_name, unsigned int comp=0)
const Elem *& _current_master
virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType, unsigned int)
This is the virtual that derived classes should override for computing the off-diag Jacobian...