20 #include "libmesh/string_to_enum.h" 25 MooseEnum orders(
"FIRST SECOND THIRD FOURTH",
"FIRST");
27 params.
addParam<BoundaryName>(
"secondary",
"The boundary ID associated with the secondary side");
28 params.
addParam<BoundaryName>(
"primary",
"The boundary ID associated with the primary side");
30 "Tangential distance to extend edges of contact surfaces");
32 "normal_smoothing_distance",
33 "Distance from edge in parametric coordinates over which to smooth contact normal");
34 params.
addParam<std::string>(
"normal_smoothing_method",
35 "Method to use to smooth normals (edge_based|nodal_normal_based)");
36 params.
addParam<
MooseEnum>(
"order", orders,
"The finite element order used for projections");
38 params.
addCoupledVar(
"primary_variable",
"The variable on the primary side of the domain");
50 _secondary(_mesh.
getBoundaryID(getParam<BoundaryName>(
"secondary"))),
51 _primary(_mesh.
getBoundaryID(getParam<BoundaryName>(
"primary"))),
52 _var(_sys.getFieldVariable<
Real>(_tid, parameters.
get<NonlinearVariableName>(
"variable"))),
54 _primary_q_point(_assembly.qPointsFace()),
55 _primary_qrule(_assembly.qRuleFace()),
58 getPenetrationLocator(getParam<BoundaryName>(
"primary"),
59 getParam<BoundaryName>(
"secondary"),
62 _current_node(_var.node()),
63 _current_primary(_var.neighbor()),
64 _u_secondary(_var.dofValues()),
68 _primary_var(*getVar(
"primary_variable", 0)),
69 _primary_var_num(_primary_var.number()),
71 _phi_primary(_assembly.phiFaceNeighbor(_primary_var)),
72 _grad_phi_primary(_assembly.gradPhiFaceNeighbor(_primary_var)),
74 _test_primary(_var.phiFaceNeighbor()),
75 _grad_test_primary(_var.gradPhiFaceNeighbor()),
77 _u_primary(_primary_var.slnNeighbor()),
78 _grad_u_primary(_primary_var.gradSlnNeighbor()),
80 _dof_map(_sys.dofMap()),
81 _node_to_elem_map(_mesh.nodeToElemMap()),
83 _overwrite_secondary_residual(true),
84 _primary_JxW(_assembly.JxWNeighbor())
130 "The secondary residual has not yet been computed, so the value will be garbage!");
238 for (
_j = 0;
_j < primary_jsize;
_j++)
250 for (
_j = 0;
_j < primary_jsize;
_j++)
261 std::set<dof_id_type> unique_dof_indices;
264 mooseAssert(node_to_elem_pair !=
_node_to_elem_map.end(),
"Missing entry in node to elem map");
265 const std::vector<dof_id_type> & elems = node_to_elem_pair->second;
268 for (
const auto & cur_elem : elems)
270 std::vector<dof_id_type> dof_indices;
274 for (
const auto & dof : dof_indices)
275 unique_dof_indices.insert(dof);
278 for (
const auto & dof : unique_dof_indices)
288 const std::set<BoundaryID> &
296 std::set<SubdomainID>
MooseMesh & _mesh
Reference to this Kernel's mesh object.
bool _secondary_residual_computed
Whether the secondary residual has been computed.
BoundaryID _secondary
Boundary ID for the secondary surface.
virtual void computeResidual() override
Computes the residual Nodal residual.
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
virtual void computeSecondaryValue(NumericVector< Number > ¤t_solution)
Compute the value the secondary node should have at the beginning of a timestep.
virtual Elem * elemPtr(const dof_id_type i)
void setNormalSmoothingDistance(Real normal_smoothing_distance)
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpResidual(Moose::ConstraintType type)=0
This is the virtual that derived classes should override for computing the residual on neighboring el...
MooseVariable & _primary_var
Primary side variable.
std::set< SubdomainID > getSecondaryConnectedBlocks() const
T * get(const std::unique_ptr< T > &u)
The MooseUtils::get() specializations are used to support making forwards-compatible code changes fro...
Base class for all Constraint types.
DenseMatrix< Number > _Kne
The Jacobian corresponding to the derivatives of the neighbor/primary residual with respect to the el...
const std::map< dof_id_type, std::vector< dof_id_type > > & _node_to_elem_map
Real secondaryResidual() const
void assignTaggedLocalResidual()
Local residual blocks will assigned as the current local kernel residual.
DenseMatrix< Number > _Ken
The Jacobian corresponding to the derivatives of the elemental/secondary residual with respect to the...
static InputParameters validParams()
virtual void getDofIndices(const Elem *, std::vector< dof_id_type > &) const
This class provides an interface for common operations on field variables of both FE and FV types wit...
const VariableTestValue & _test_primary
Side test function.
void setTangentialTolerance(Real tangential_tolerance)
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
const MooseVariableFieldBase & getVariable(unsigned int jvar_num) const
Retrieve the variable object from our system associated with jvar_num.
void residualSetup() override
Gets called just before the residual is computed and before this object is asked to do its job...
const Node *const & _current_node
current node being processed
DenseMatrix< Number > _Kee
The Jacobian corresponding to the derivatives of the elemental/secondary residual with respect to the...
SystemBase & _sys
Reference to the EquationSystem object.
void prepareMatrixTagNeighbor(Assembly &assembly, unsigned int ivar, unsigned int jvar, Moose::DGJacobianType type)
Prepare data for computing element jacobian according to the active tags for DG and interface kernels...
BoundaryID getBoundaryID(const BoundaryName &boundary_name, const MeshBase &mesh)
Gets the boundary ID associated with the given BoundaryName.
NodeFaceConstraint(const InputParameters ¶meters)
Enhances MooseVariableInterface interface provide values from neighbor elements.
virtual const std::vector< dof_id_type > & dofIndicesNeighbor() const =0
Get neighbor DOF indices for currently selected element.
unsigned int size() const
The number of elements that can currently be stored in the array.
virtual void computeJacobian() override
Computes the jacobian for the current element.
std::set< BoundaryID > _boundary_ids
the union of the secondary and primary boundary ids
T string_to_enum(const std::string &s)
VarKindType
Framework-wide stuff.
virtual void getConnectedDofIndices(unsigned int var_num)
Gets the indices for all dofs connected to the constraint.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
std::set< SubdomainID > getBoundaryConnectedBlocks(const BoundaryID bid) const
Get the list of subdomains associated with the given boundary.
Real _secondary_residual
The value of the secondary residual.
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
const std::set< BoundaryID > & buildBoundaryIDs()
Builds the _boundary_ids data member and returns it.
Assembly & _assembly
Reference to this Kernel's assembly object.
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar...
VariableTestValue _test_secondary
Shape function on the secondary side. This will always only have one entry and that entry will always...
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
void prepareVectorTagNeighbor(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual the according to active tags for DG and interface kernels...
virtual Real computeQpSecondaryValue()=0
Compute the value the secondary node should have at the beginning of a timestep.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const dof_id_type & nodalDofIndex() const override
void release()
Manually deallocates the data pointer.
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
void resize(const unsigned int new_m, const unsigned int new_n)
void resize(unsigned int size)
Change the number of elements the array can store.
std::vector< dof_id_type > _connected_dof_indices
const InputParameters & parameters() const
Get the parameters of the object.
const VariablePhiValue & _phi_primary
Side shape function.
BoundaryID _primary
Boundary ID for the primary surface.
bool _overwrite_secondary_residual
Whether or not the secondary's residual should be overwritten.
virtual void set(const numeric_index_type i, const Number value)=0
VariablePhiValue _phi_secondary
Shape function on the secondary side. This will always.
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
void setNormalSmoothingMethod(std::string nsmString)
MooseVariableFieldBase & getVariable(THREAD_ID tid, const std::string &var_name) const
Gets a reference to a variable of with specified name.
PenetrationLocator & _penetration_locator
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual according to active tags.
virtual bool overwriteSecondaryResidual()
Whether or not the secondary's residual should be overwritten.
static InputParameters validParams()
virtual Real computeQpJacobian(Moose::ConstraintJacobianType type)=0
This is the virtual that derived classes should override for computing the Jacobian on neighboring el...
virtual ~NodeFaceConstraint()
virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType, unsigned int)
This is the virtual that derived classes should override for computing the off-diag Jacobian...