www.mooseframework.org
NodeElemConstraint.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 NodeElemConstraint;
18 
19 // libMesh forward declarations
20 namespace libMesh
21 {
22 template <typename T>
23 class SparseMatrix;
24 }
25 
26 template <>
28 
37 {
38 public:
40  virtual ~NodeElemConstraint();
41 
43  void computeSlaveValue(NumericVector<Number> & current_solution);
44 
46  virtual void computeResidual();
47 
49  virtual void computeJacobian();
50 
52  virtual void computeOffDiagJacobian(unsigned int jvar);
53 
55  virtual void getConnectedDofIndices(unsigned int var_num);
56 
61  virtual bool shouldApply() { return true; }
62 
68  virtual bool overwriteSlaveResidual();
69 
75  virtual bool overwriteSlaveJacobian() { return overwriteSlaveResidual(); };
76 
81  virtual MooseVariable & masterVariable() { return _master_var; }
82 
86  MooseVariable & variable() { return _var; }
87 
88 protected:
90  virtual void prepareSlaveToMasterMap() = 0;
91 
93  virtual Real computeQpSlaveValue() = 0;
94 
97 
100 
103  unsigned int /*jvar*/)
104  {
105  return 0;
106  }
107 
108  // coupling interface:
109  virtual const VariableValue & coupledSlaveValue(const std::string & var_name,
110  unsigned int comp = 0)
111  {
112  return coupledValue(var_name, comp);
113  }
114  virtual const VariableValue & coupledSlaveValueOld(const std::string & var_name,
115  unsigned int comp = 0)
116  {
117  return coupledValueOld(var_name, comp);
118  }
119  virtual const VariableValue & coupledSlaveValueOlder(const std::string & var_name,
120  unsigned int comp = 0)
121  {
122  return coupledValueOlder(var_name, comp);
123  }
124 
125  virtual const VariableGradient & coupledSlaveGradient(const std::string & var_name,
126  unsigned int comp = 0)
127  {
128  return coupledGradient(var_name, comp);
129  }
130  virtual const VariableGradient & coupledSlaveGradientOld(const std::string & var_name,
131  unsigned int comp = 0)
132  {
133  return coupledGradientOld(var_name, comp);
134  }
135  virtual const VariableGradient & coupledSlaveGradientOlder(const std::string & var_name,
136  unsigned int comp = 0)
137  {
138  return coupledGradientOlder(var_name, comp);
139  }
140 
141  virtual const VariableSecond & coupledSlaveSecond(const std::string & var_name,
142  unsigned int comp = 0)
143  {
144  return coupledSecond(var_name, comp);
145  }
146 
147  virtual const VariableValue & coupledMasterValue(const std::string & var_name,
148  unsigned int comp = 0)
149  {
150  return coupledNeighborValue(var_name, comp);
151  }
152  virtual const VariableValue & coupledMasterValueOld(const std::string & var_name,
153  unsigned int comp = 0)
154  {
155  return coupledNeighborValueOld(var_name, comp);
156  }
157  virtual const VariableValue & coupledMasterValueOlder(const std::string & var_name,
158  unsigned int comp = 0)
159  {
160  return coupledNeighborValueOlder(var_name, comp);
161  }
162 
163  virtual const VariableGradient & coupledMasterGradient(const std::string & var_name,
164  unsigned int comp = 0)
165  {
166  return coupledNeighborGradient(var_name, comp);
167  }
168  virtual const VariableGradient & coupledMasterGradientOld(const std::string & var_name,
169  unsigned int comp = 0)
170  {
171  return coupledNeighborGradientOld(var_name, comp);
172  }
173  virtual const VariableGradient & coupledMasterGradientOlder(const std::string & var_name,
174  unsigned int comp = 0)
175  {
176  return coupledNeighborGradientOlder(var_name, comp);
177  }
178 
179  virtual const VariableSecond & coupledMasterSecond(const std::string & var_name,
180  unsigned int comp = 0)
181  {
182  return coupledNeighborSecond(var_name, comp);
183  }
184 
186  unsigned short _slave;
188  unsigned short _master;
189 
191 
193  const QBase * const & _master_qrule;
194 
196  const Node * const & _current_node;
197  const Elem * const & _current_elem;
198 
207 
210 
212  unsigned int _master_var_num;
213 
218 
223 
230 
232  const DofMap & _dof_map;
233 
234  const std::map<dof_id_type, std::vector<dof_id_type>> & _node_to_elem_map;
235 
237  std::map<dof_id_type, dof_id_type> _slave_to_master_map;
238 
246 
247 public:
248  SparseMatrix<Number> * _jacobian;
250  std::vector<dof_id_type> _connected_dof_indices;
252  DenseMatrix<Number> _Kne;
254  DenseMatrix<Number> _Kee;
255 };
256 
OutputTools< Real >::VariableGradient VariableGradient
Definition: MooseTypes.h:198
virtual const VariableGradient & coupledSlaveGradientOlder(const std::string &var_name, unsigned int comp=0)
std::map< dof_id_type, dof_id_type > _slave_to_master_map
maps slave node ids to master element ids
virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType, unsigned int)
This is the virtual that derived classes should override for computing the off-diag Jacobian...
const QBase *const & _master_qrule
virtual const VariableValue & coupledMasterValue(const std::string &var_name, unsigned int comp=0)
NodeElemConstraint(const InputParameters &parameters)
SparseMatrix< Number > * _jacobian
ConstraintType
Definition: MooseTypes.h:523
const VariableTestGradient & _grad_test_master
Gradient of side shape function.
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
const DofMap & _dof_map
DOF map.
DenseMatrix< Number > _Kee
stiffness matrix holding slave-slave jacobian
const VariableGradient & _grad_u_master
Holds the current solution gradient at the current quadrature point.
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 const VariableGradient & coupledMasterGradientOld(const std::string &var_name, unsigned int comp=0)
virtual MooseVariable & masterVariable()
The variable on the master elem.
virtual const VariableSecond & coupledMasterSecond(const std::string &var_name, unsigned int comp=0)
virtual void computeResidual()
Computes the residual Nodal residual.
const VariableValue & _u_master
Holds the current solution at the current quadrature point.
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.
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
const VariableValue & _u_master_old
Holds the old solution at the current quadrature point.
virtual Real computeQpSlaveValue()=0
Compute the value the slave node should have at the beginning of a timestep.
const VariableTestValue & _test_master
Side test function.
virtual void computeOffDiagJacobian(unsigned int jvar)
Computes d-residual / d-jvar...
virtual const VariableGradient & coupledMasterGradient(const std::string &var_name, unsigned int comp=0)
MooseVariable & variable()
The variable number that this object operates on.
virtual const VariableSecond & coupledNeighborSecond(const std::string &var_name, unsigned int i=0)
const std::string & type() const
Get the type of this object.
Definition: MooseObject.h:53
MooseVariable & _var
virtual void getConnectedDofIndices(unsigned int var_num)
Gets the indices for all dofs connected to the constraint.
unsigned short _master
master block id
Enhances MooseVariableInterface interface provide values from neighbor elements.
virtual const VariableGradient & coupledSlaveGradientOld(const std::string &var_name, unsigned int comp=0)
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseObject.h:65
OutputTools< Real >::VariableTestValue VariableTestValue
Definition: MooseTypes.h:205
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
virtual void prepareSlaveToMasterMap()=0
prepare the _slave_to_master_map
const VariablePhiValue & _phi_master
Side shape function.
virtual const VariableGradient & coupledSlaveGradient(const std::string &var_name, unsigned int comp=0)
VariablePhiValue _phi_slave
Shape function on the slave side.
unsigned short _slave
slave block id
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
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 const VariableValue & coupledNeighborValue(const std::string &var_name, unsigned int comp=0)
const VariableValue & _u_slave
Value of the unknown variable on the slave node.
const VariableValue & _u_slave_old
old solution
VariableTestValue _test_slave
Shape function on the slave side. This will always only have one entry and that entry will always be ...
const Elem *const & _current_elem
virtual const VariableValue & coupledMasterValueOld(const std::string &var_name, unsigned int comp=0)
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
virtual const VariableValue & coupledSlaveValue(const std::string &var_name, unsigned int comp=0)
unsigned int _master_var_num
Number for the master variable.
const Node *const & _current_node
current node being processed
virtual const VariableValue & coupledSlaveValueOld(const std::string &var_name, unsigned int comp=0)
ConstraintJacobianType
Definition: MooseTypes.h:543
virtual const VariableGradient & coupledGradient(const std::string &var_name, unsigned int comp=0)
Returns gradient of a coupled variable.
Definition: Coupleable.C:892
const VariablePhiGradient & _grad_phi_master
Gradient of side shape function.
virtual bool shouldApply()
Whether or not this constraint should be applied.
std::vector< dof_id_type > _connected_dof_indices
dofs connected to the slave node
virtual const VariableSecond & coupledSecond(const std::string &var_name, unsigned int comp=0)
Returns second derivative of a coupled variable.
Definition: Coupleable.C:1145
virtual const VariableValue & coupledMasterValueOlder(const std::string &var_name, unsigned int comp=0)
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.
const MooseArray< Point > & _master_q_point
OutputTools< Real >::VariableSecond VariableSecond
Definition: MooseTypes.h:199
virtual const VariableValue & coupledSlaveValueOlder(const std::string &var_name, unsigned int comp=0)
virtual void computeJacobian()
Computes the jacobian for the current element.
OutputTools< Real >::VariableTestGradient VariableTestGradient
Definition: MooseTypes.h:206
A NodeElemConstraint is used when you need to create constraints between a slave node and a master el...
virtual const VariableGradient & coupledNeighborGradient(const std::string &var_name, unsigned int comp=0)
virtual bool overwriteSlaveJacobian()
Whether or not the slave&#39;s Jacobian row should be overwritten.
DenseMatrix< Number > _Kne
stiffness matrix holding master-slave jacobian
InputParameters validParams< NodeElemConstraint >()
virtual const VariableSecond & coupledSlaveSecond(const std::string &var_name, unsigned int comp=0)
virtual const VariableGradient & coupledMasterGradientOlder(const std::string &var_name, unsigned int comp=0)
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 & 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
bool _overwrite_slave_residual
Whether or not the slave&#39;s residual should be overwritten.