www.mooseframework.org
NodeFaceConstraint.h
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 #pragma once
11 
12 // MOOSE includes
13 #include "Constraint.h"
15 
16 // Forward Declarations
17 class NodeFaceConstraint;
18 
19 // libMesh forward declarations
20 namespace libMesh
21 {
22 template <typename T>
23 class SparseMatrix;
24 }
25 
26 template <>
28 
41 {
42 public:
44  virtual ~NodeFaceConstraint();
45 
49  void computeSlaveValue(NumericVector<Number> & current_solution);
50 
54  virtual void computeResidual();
55 
59  virtual void computeJacobian();
60 
64  virtual void computeOffDiagJacobian(unsigned int jvar);
65 
69  virtual void getConnectedDofIndices(unsigned int var_num);
70 
76  virtual bool shouldApply() { return true; }
77 
84  virtual bool overwriteSlaveResidual();
85 
92  virtual bool overwriteSlaveJacobian() { return overwriteSlaveResidual(); };
93 
97  virtual MooseVariable & masterVariable() { return _master_var; }
98 
102  MooseVariable & variable() { return _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 
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 
220  const QBase * const & _master_qrule;
221 
222 public:
224 
225 protected:
227  const Node * const & _current_node;
228  const Elem * const & _current_master;
229 
236 
239 
241  unsigned int _master_var_num;
242 
247 
252 
257 
259  const DofMap & _dof_map;
260 
261  const std::map<dof_id_type, std::vector<dof_id_type>> & _node_to_elem_map;
262 
270 
271 public:
272  std::vector<dof_id_type> _connected_dof_indices;
273 
274  DenseMatrix<Number> _Kne;
275  DenseMatrix<Number> _Kee;
276 };
277 
MooseVariable & variable()
The variable number that this object operates on.
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 void computeOffDiagJacobian(unsigned int jvar)
Computes d-residual / d-jvar...
const Elem *const & _current_master
virtual const VariableGradient & coupledSlaveGradient(const std::string &var_name, unsigned int comp=0)
ConstraintType
Definition: MooseTypes.h:523
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
MooseVariable & _var
virtual bool overwriteSlaveJacobian()
Whether or not the slave&#39;s Jacobian row should be overwritten.
virtual Real computeQpResidual(Moose::ConstraintType type)=0
This is the virtual that derived classes should override for computing the residual on neighboring el...
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:39
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
virtual void computeResidual()
Computes the residual Nodal residual.
const Node *const & _current_node
current node being processed
DenseMatrix< Number > _Kee
virtual const VariableSecond & coupledNeighborSecond(const std::string &var_name, unsigned int i=0)
const VariableTestGradient & _grad_test_master
Gradient of side shape function.
const std::string & type() const
Get the type of this object.
Definition: MooseObject.h:53
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.
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.
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseObject.h:65
InputParameters validParams< NodeFaceConstraint >()
OutputTools< Real >::VariableTestValue VariableTestValue
Definition: MooseTypes.h:205
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...
virtual Real computeQpSlaveValue()=0
Compute the value the slave node should have at the beginning of a timestep.
virtual const VariableGradient & coupledGradientOlder(const std::string &var_name, unsigned int comp=0)
Returns an old gradient from two time steps previous of a coupled variable.
Definition: Coupleable.C:931
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 VariableValue & _u_slave
Value of the unknown variable this BC is action on.
virtual const VariableValue & coupledMasterValue(const std::string &var_name, unsigned int comp=0)
virtual void computeJacobian()
Computes the jacobian for the current element.
virtual const VariableGradient & coupledNeighborGradientOlder(const std::string &var_name, unsigned int comp=0)
virtual const VariableValue & coupledValue(const std::string &var_name, unsigned int comp=0)
Returns value of a coupled variable.
Definition: Coupleable.C:361
unsigned int _master_var_num
Number for the master variable.
virtual const VariableValue & coupledNeighborValueOld(const std::string &var_name, unsigned int comp=0)
virtual const VariableValue & coupledNeighborValueOlder(const std::string &var_name, unsigned int comp=0)
virtual const VariableGradient & coupledNeighborGradientOld(const std::string &var_name, unsigned int comp=0)
virtual MooseVariable & masterVariable()
The variable on the Master side of the domain.
virtual const VariableValue & coupledNeighborValue(const std::string &var_name, unsigned int comp=0)
const MooseArray< Point > & _master_q_point
virtual const VariableValue & coupledMasterValueOlder(const std::string &var_name, unsigned int comp=0)
OutputTools< Real >::VariableValue VariableValue
Definition: MooseTypes.h:197
OutputTools< Real >::VariablePhiGradient VariablePhiGradient
Definition: MooseTypes.h:202
virtual const VariableValue & coupledValueOld(const std::string &var_name, unsigned int comp=0)
Returns an old value from previous time step of a coupled variable.
Definition: Coupleable.C:457
forward declarations
Definition: MooseArray.h:16
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:543
const QBase *const & _master_qrule
virtual const VariableGradient & coupledGradient(const std::string &var_name, unsigned int comp=0)
Returns gradient of a coupled variable.
Definition: Coupleable.C:892
virtual const VariableSecond & coupledSecond(const std::string &var_name, unsigned int comp=0)
Returns second derivative of a coupled variable.
Definition: Coupleable.C:1145
std::vector< dof_id_type > _connected_dof_indices
OutputTools< Real >::VariableSecond VariableSecond
Definition: MooseTypes.h:199
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
OutputTools< Real >::VariableTestGradient VariableTestGradient
Definition: MooseTypes.h:206
virtual const VariableGradient & coupledNeighborGradient(const std::string &var_name, unsigned int comp=0)
virtual Real computeQpJacobian(Moose::ConstraintJacobianType type)=0
This is the virtual that derived classes should override for computing the Jacobian on neighboring el...
void computeSlaveValue(NumericVector< Number > &current_solution)
Compute the value the slave node should have at the beginning of a timestep.
virtual const VariableGradient & coupledGradientOld(const std::string &var_name, unsigned int comp=0)
Returns an old gradient from previous time step of a coupled variable.
Definition: Coupleable.C:911
virtual const VariableValue & coupledSlaveValueOld(const std::string &var_name, unsigned int comp=0)
virtual const VariableValue & coupledValueOlder(const std::string &var_name, unsigned int comp=0)
Returns an old value from two time steps previous of a coupled variable.
Definition: Coupleable.C:484
virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType, unsigned int)
This is the virtual that derived classes should override for computing the off-diag Jacobian...