17 #include "libmesh/sparse_matrix.h" 23 MooseEnum formulationtype(
"penalty kinematic",
"penalty");
26 "Formulation used to calculate constraint - penalty or kinematic.");
27 params.
addParam<NonlinearVariableName>(
"variable_secondary",
28 "The name of the variable for the secondary nodes, if it " 29 "is different from the primary nodes' variable");
38 _var(_sys.getFieldVariable<
Real>(_tid, parameters.
get<NonlinearVariableName>(
"variable"))),
39 _var_secondary(_sys.getFieldVariable<
Real>(
41 isParamValid(
"variable_secondary")
42 ? parameters.
get<NonlinearVariableName>(
"variable_secondary")
43 : parameters.
get<NonlinearVariableName>(
"variable"))),
44 _u_secondary(_var_secondary.dofValuesNeighbor()),
45 _u_primary(_var.dofValues())
50 MooseEnum temp_formulation = getParam<MooseEnum>(
"formulation");
51 if (temp_formulation ==
"penalty")
53 else if (temp_formulation ==
"kinematic")
56 mooseError(
"Formulation must be either Penalty or Kinematic");
74 for (
_i = 0;
_i < secondarydof.size(); ++
_i)
76 for (
_j = 0;
_j < primarydof.size(); ++
_j)
86 Real res = residual(secondarydof[
_i]);
95 if (!primarydof.empty())
97 if (!secondarydof.empty())
119 for (
_i = 0;
_i < secondarydof.size(); ++
_i)
121 for (
_j = 0;
_j < primarydof.size(); ++
_j)
133 Kne(
_i,
_j) += -jacobian(secondarydof[
_i], primarydof[
_j]) / primarydof.size() +
144 for (
_i = 0;
_i < secondarydof.size(); ++
_i)
153 value = -jacobian(secondarydof[
_i], secondarydof[
_i]) / primarydof.size() +
virtual void updateConnectivity()
Built the connectivity for this constraint.
virtual void zero() override final
virtual void zero() override final
virtual Real computeQpResidual(Moose::ConstraintType type)=0
This is the virtual that derived classes should override for computing the residual on neighboring el...
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.
void addResiduals(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided incoming residuals corresponding to the provided dof indices.
virtual void computeResidual() override final
Computes the nodal residual.
static InputParameters validParams()
unsigned int _i
Counter for primary and secondary nodes.
std::vector< dof_id_type > _primary_node_vector
node IDs of the primary node
void addJacobian(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided residual derivatives into the Jacobian for the provided dof indices.
Enhances MooseVariableInterface interface provide values from neighbor elements.
const std::vector< dof_id_type > & dofIndices() const final
Get local DoF indices.
virtual Real computeQpJacobian(Moose::ConstraintJacobianType type)=0
This is the virtual that derived classes should override for computing the Jacobian on neighboring el...
std::vector< Real > _weights
When the secondary node is constrained to move as a linear combination of the primary nodes...
Moose::ConstraintFormulationType _formulation
Specifies formulation type used to apply constraints.
VarKindType
Framework-wide stuff.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
virtual const OutputTools< Real >::VariableValue & value()
The value of the variable this object is operating on.
Assembly & _assembly
Reference to this Kernel's assembly object.
virtual void computeJacobian() override final
Computes the jacobian for the current element.
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const std::vector< dof_id_type > & dofIndicesNeighbor() const final
Get neighbor DOF indices for currently selected element.
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
void addJacobianElement(Assembly &assembly, Real value, dof_id_type row_index, dof_id_type column_index, Real scaling_factor)
Add into a single Jacobian element.
NodalConstraint(const InputParameters ¶meters)
MooseVariable & _var_secondary
static InputParameters validParams()
void scalingFactor(const std::vector< Real > &factor)
Set the scaling factor for this variable.