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 // libMesh forward declarations
17 namespace libMesh
18 {
19 template <typename T>
20 class SparseMatrix;
21 }
22 
31 {
32 public:
34 
36  virtual ~NodeElemConstraint();
37 
39  void computeSecondaryValue(NumericVector<Number> & current_solution);
40 
42  virtual void computeResidual() override;
43 
45  virtual void computeJacobian() override;
46 
48  virtual void computeOffDiagJacobian(unsigned int jvar) override;
49 
51  virtual void getConnectedDofIndices(unsigned int var_num);
52 
57  virtual bool shouldApply() { return true; }
58 
64  virtual bool overwriteSecondaryResidual();
65 
72 
77  virtual MooseVariable & primaryVariable() { return _primary_var; }
78 
82  const MooseVariable & variable() const override { return _var; }
83 
84 protected:
86  virtual void prepareSecondaryToPrimaryMap() = 0;
87 
89  virtual Real computeQpSecondaryValue() = 0;
90 
93 
96 
99  unsigned int /*jvar*/)
100  {
101  return 0;
102  }
103 
104  // coupling interface:
105  virtual const VariableValue & coupledSecondaryValue(const std::string & var_name,
106  unsigned int comp = 0)
107  {
108  return coupledValue(var_name, comp);
109  }
110  virtual const VariableValue & coupledSecondaryValueOld(const std::string & var_name,
111  unsigned int comp = 0)
112  {
113  return coupledValueOld(var_name, comp);
114  }
115  virtual const VariableValue & coupledSecondaryValueOlder(const std::string & var_name,
116  unsigned int comp = 0)
117  {
118  return coupledValueOlder(var_name, comp);
119  }
120 
121  virtual const VariableGradient & coupledSecondaryGradient(const std::string & var_name,
122  unsigned int comp = 0)
123  {
124  return coupledGradient(var_name, comp);
125  }
126  virtual const VariableGradient & coupledSecondaryGradientOld(const std::string & var_name,
127  unsigned int comp = 0)
128  {
129  return coupledGradientOld(var_name, comp);
130  }
131  virtual const VariableGradient & coupledSecondaryGradientOlder(const std::string & var_name,
132  unsigned int comp = 0)
133  {
134  return coupledGradientOlder(var_name, comp);
135  }
136 
137  virtual const VariableSecond & coupledSecondarySecond(const std::string & var_name,
138  unsigned int comp = 0)
139  {
140  return coupledSecond(var_name, comp);
141  }
142 
143  virtual const VariableValue & coupledPrimaryValue(const std::string & var_name,
144  unsigned int comp = 0)
145  {
146  return coupledNeighborValue(var_name, comp);
147  }
148  virtual const VariableValue & coupledPrimaryValueOld(const std::string & var_name,
149  unsigned int comp = 0)
150  {
151  return coupledNeighborValueOld(var_name, comp);
152  }
153  virtual const VariableValue & coupledPrimaryValueOlder(const std::string & var_name,
154  unsigned int comp = 0)
155  {
156  return coupledNeighborValueOlder(var_name, comp);
157  }
158 
159  virtual const VariableGradient & coupledPrimaryGradient(const std::string & var_name,
160  unsigned int comp = 0)
161  {
162  return coupledNeighborGradient(var_name, comp);
163  }
164  virtual const VariableGradient & coupledPrimaryGradientOld(const std::string & var_name,
165  unsigned int comp = 0)
166  {
167  return coupledNeighborGradientOld(var_name, comp);
168  }
169  virtual const VariableGradient & coupledPrimaryGradientOlder(const std::string & var_name,
170  unsigned int comp = 0)
171  {
172  return coupledNeighborGradientOlder(var_name, comp);
173  }
174 
175  virtual const VariableSecond & coupledPrimarySecond(const std::string & var_name,
176  unsigned int comp = 0)
177  {
178  return coupledNeighborSecond(var_name, comp);
179  }
180 
182  unsigned short _secondary;
184  unsigned short _primary;
185 
187 
189  const QBase * const & _primary_qrule;
190 
192  const Node * const & _current_node;
193  const Elem * const & _current_elem;
194 
203 
206 
208  unsigned int _primary_var_num;
209 
214 
219 
226 
228  const DofMap & _dof_map;
229 
230  const std::map<dof_id_type, std::vector<dof_id_type>> & _node_to_elem_map;
231 
233  std::map<dof_id_type, dof_id_type> _secondary_to_primary_map;
234 
242 
243 public:
246  std::vector<dof_id_type> _connected_dof_indices;
251 };
OutputTools< Real >::VariableGradient VariableGradient
Definition: MooseTypes.h:303
unsigned short _secondary
secondary block id
virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType, unsigned int)
This is the virtual that derived classes should override for computing the off-diag Jacobian...
static InputParameters validParams()
NodeElemConstraint(const InputParameters &parameters)
SparseMatrix< Number > * _jacobian
ConstraintType
Definition: MooseTypes.h:670
const VariableTestValue & _test_primary
Side test function.
virtual const VariableSecond & coupledSecond(const std::string &var_name, unsigned int comp=0) const
Returns second spatial derivatives of a coupled variable.
Definition: Coupleable.C:1722
Class for stuff related to variables.
Definition: Adaptivity.h:31
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
const QBase *const & _primary_qrule
const DofMap & _dof_map
DOF map.
bool _overwrite_secondary_residual
Whether or not the secondary&#39;s residual should be overwritten.
DenseMatrix< Number > _Kee
stiffness matrix holding secondary-secondary jacobian
unsigned short _primary
primary block id
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
virtual const VariableGradient & coupledNeighborGradient(const std::string &var_name, unsigned int comp=0) const
virtual MooseVariable & primaryVariable()
The variable on the primary elem.
Base class for all Constraint types.
Definition: Constraint.h:19
virtual const VariableSecond & coupledSecondarySecond(const std::string &var_name, unsigned int comp=0)
virtual const VariableValue & coupledNeighborValueOlder(const std::string &var_name, unsigned int comp=0) const
virtual const VariableGradient & coupledNeighborGradientOld(const std::string &var_name, unsigned int comp=0) const
virtual void computeJacobian() override
Computes the jacobian for the current element.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual const VariableGradient & coupledSecondaryGradientOld(const std::string &var_name, unsigned int comp=0)
virtual const VariableGradient & coupledSecondaryGradient(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.
virtual const VariableValue & coupledValueOld(const std::string &var_name, unsigned int comp=0) const
Returns an old value from previous time step of a coupled variable.
Definition: Coupleable.C:981
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
const VariableValue & _u_secondary
Value of the unknown variable on the secondary node.
virtual const VariableValue & coupledSecondaryValueOld(const std::string &var_name, unsigned int comp=0)
virtual bool overwriteSecondaryResidual()
Whether or not the secondary&#39;s residual should be overwritten.
virtual const VariableGradient & coupledGradient(const std::string &var_name, unsigned int comp=0) const
Returns gradient of a coupled variable.
Definition: Coupleable.C:1477
virtual const VariableValue & coupledValueOlder(const std::string &var_name, unsigned int comp=0) const
Returns an old value from two time steps previous of a coupled variable.
Definition: Coupleable.C:1003
virtual Real computeQpSecondaryValue()=0
Compute the value the secondary node should have at the beginning of a timestep.
MooseVariable & _var
virtual const VariableValue & coupledValue(const std::string &var_name, unsigned int comp=0) const
Returns value of a coupled variable.
Definition: Coupleable.C:478
virtual void getConnectedDofIndices(unsigned int var_num)
Gets the indices for all dofs connected to the constraint.
virtual bool overwriteSecondaryJacobian()
Whether or not the secondary&#39;s Jacobian row should be overwritten.
Enhances MooseVariableInterface interface provide values from neighbor elements.
virtual const VariableSecond & coupledPrimarySecond(const std::string &var_name, unsigned int comp=0)
virtual const VariableGradient & coupledGradientOld(const std::string &var_name, unsigned int comp=0) const
Returns an old gradient from previous time step of a coupled variable.
Definition: Coupleable.C:1493
virtual const VariableValue & coupledNeighborValue(const std::string &var_name, unsigned int comp=0) const
virtual const VariableValue & coupledPrimaryValue(const std::string &var_name, unsigned int comp=0)
const VariableValue & _u_primary
Holds the current solution at the current quadrature point.
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar...
OutputTools< Real >::VariableTestValue VariableTestValue
Definition: MooseTypes.h:312
std::map< dof_id_type, dof_id_type > _secondary_to_primary_map
maps secondary node ids to primary element ids
const VariableTestGradient & _grad_test_primary
Gradient of side shape function.
const std::string & type() const
Get the type of this class.
Definition: MooseBase.h:51
virtual const VariableGradient & coupledNeighborGradientOlder(const std::string &var_name, unsigned int comp=0) const
const VariableValue & _u_primary_old
Holds the old solution at the current quadrature point.
virtual const VariableValue & coupledSecondaryValue(const std::string &var_name, unsigned int comp=0)
MooseVariable & _primary_var
Primary side variable.
const Elem *const & _current_elem
OutputTools< Real >::VariablePhiGradient VariablePhiGradient
Definition: MooseTypes.h:308
const MooseVariable & variable() const override
The variable number that this object operates on.
forward declarations
Definition: MooseArray.h:17
const VariablePhiValue & _phi_primary
Side shape function.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const MooseArray< Point > & _primary_q_point
unsigned int _primary_var_num
Number for the primary variable.
virtual const VariableValue & coupledPrimaryValueOld(const std::string &var_name, unsigned int comp=0)
void computeSecondaryValue(NumericVector< Number > &current_solution)
Compute the value the secondary node should have at the beginning of a timestep.
const VariableGradient & _grad_u_primary
Holds the current solution gradient at the current quadrature point.
virtual void computeResidual() override
Computes the residual Nodal residual.
const Node *const & _current_node
current node being processed
ConstraintJacobianType
Definition: MooseTypes.h:709
virtual bool shouldApply()
Whether or not this constraint should be applied.
std::vector< dof_id_type > _connected_dof_indices
dofs connected to the secondary node
virtual const VariableGradient & coupledPrimaryGradientOlder(const std::string &var_name, unsigned int comp=0)
const InputParameters & parameters() const
Get the parameters of the object.
OutputTools< Real >::VariableSecond VariableSecond
Definition: MooseTypes.h:304
virtual const VariableValue & coupledNeighborValueOld(const std::string &var_name, unsigned int comp=0) const
virtual const VariableValue & coupledSecondaryValueOlder(const std::string &var_name, unsigned int comp=0)
virtual const VariableValue & coupledPrimaryValueOlder(const std::string &var_name, unsigned int comp=0)
OutputTools< Real >::VariableTestGradient VariableTestGradient
Definition: MooseTypes.h:313
A NodeElemConstraint is used when you need to create constraints between a secondary node and a prima...
VariableTestValue _test_secondary
Shape function on the secondary side. This will always only have one entry and that entry will always...
const VariableValue & _u_secondary_old
old solution
DenseMatrix< Number > _Kne
stiffness matrix holding primary-secondary jacobian
virtual void prepareSecondaryToPrimaryMap()=0
prepare the _secondary_to_primary_map
VariablePhiValue _phi_secondary
Shape function on the secondary side.
virtual const VariableGradient & coupledPrimaryGradient(const std::string &var_name, unsigned int comp=0)
virtual const VariableGradient & coupledPrimaryGradientOld(const std::string &var_name, unsigned int comp=0)
virtual const VariableGradient & coupledGradientOlder(const std::string &var_name, unsigned int comp=0) const
Returns an old gradient from two time steps previous of a coupled variable.
Definition: Coupleable.C:1509
virtual const VariableGradient & coupledSecondaryGradientOlder(const std::string &var_name, unsigned int comp=0)
const VariablePhiGradient & _grad_phi_primary
Gradient of side shape function.
virtual const VariableSecond & coupledNeighborSecond(const std::string &var_name, unsigned int i=0) const