www.mooseframework.org
NodalConstraint.C
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 #include "NodalConstraint.h"
11 
12 // MOOSE includes
13 #include "Assembly.h"
14 #include "MooseVariableFE.h"
15 #include "SystemBase.h"
16 
17 #include "libmesh/sparse_matrix.h"
18 
21 {
23  MooseEnum formulationtype("penalty kinematic", "penalty");
24  params.addParam<MooseEnum>("formulation",
25  formulationtype,
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");
30  return params;
31 }
32 
34  : Constraint(parameters),
38  _var(_sys.getFieldVariable<Real>(_tid, parameters.get<NonlinearVariableName>("variable"))),
39  _var_secondary(_sys.getFieldVariable<Real>(
40  _tid,
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())
46 {
49 
50  MooseEnum temp_formulation = getParam<MooseEnum>("formulation");
51  if (temp_formulation == "penalty")
53  else if (temp_formulation == "kinematic")
55  else
56  mooseError("Formulation must be either Penalty or Kinematic");
57 }
58 
59 void
61 {
62  if ((_weights.size() == 0) && (_primary_node_vector.size() == 1))
63  _weights.push_back(1.0);
64 
65  std::vector<dof_id_type> primarydof = _var.dofIndices();
66  std::vector<dof_id_type> secondarydof = _var_secondary.dofIndicesNeighbor();
67 
68  DenseVector<Number> re(primarydof.size());
69  DenseVector<Number> neighbor_re(secondarydof.size());
70 
71  re.zero();
72  neighbor_re.zero();
73 
74  for (_i = 0; _i < secondarydof.size(); ++_i)
75  {
76  for (_j = 0; _j < primarydof.size(); ++_j)
77  {
78  switch (_formulation)
79  {
80  case Moose::Penalty:
83  break;
84  case Moose::Kinematic:
85  // Transfer the current residual of the secondary node to the primary nodes
86  Real res = residual(secondarydof[_i]);
87  re(_j) += res * _weights[_j];
88  neighbor_re(_i) +=
90  break;
91  }
92  }
93  }
94  // We've already applied scaling
95  if (!primarydof.empty())
96  addResiduals(_assembly, re, primarydof, /*scaling_factor=*/1);
97  if (!secondarydof.empty())
98  addResiduals(_assembly, neighbor_re, secondarydof, /*scaling_factor=*/1);
99 }
100 
101 void
103 {
104  if ((_weights.size() == 0) && (_primary_node_vector.size() == 1))
105  _weights.push_back(1.0);
106 
107  // Calculate the dense-block Jacobian entries
108  std::vector<dof_id_type> secondarydof = _var_secondary.dofIndicesNeighbor();
109  std::vector<dof_id_type> primarydof = _var.dofIndices();
110 
111  DenseMatrix<Number> Kee(primarydof.size(), primarydof.size());
112  DenseMatrix<Number> Ken(primarydof.size(), secondarydof.size());
113  DenseMatrix<Number> Kne(secondarydof.size(), primarydof.size());
114 
115  Kee.zero();
116  Ken.zero();
117  Kne.zero();
118 
119  for (_i = 0; _i < secondarydof.size(); ++_i)
120  {
121  for (_j = 0; _j < primarydof.size(); ++_j)
122  {
123  switch (_formulation)
124  {
125  case Moose::Penalty:
129  break;
130  case Moose::Kinematic:
131  Kee(_j, _j) = 0.;
132  Ken(_j, _i) += jacobian(secondarydof[_i], primarydof[_j]) * _weights[_j];
133  Kne(_i, _j) += -jacobian(secondarydof[_i], primarydof[_j]) / primarydof.size() +
135  break;
136  }
137  }
138  }
139  addJacobian(_assembly, Kee, primarydof, primarydof, _var.scalingFactor());
140  addJacobian(_assembly, Ken, primarydof, secondarydof, _var.scalingFactor());
141  addJacobian(_assembly, Kne, secondarydof, primarydof, _var_secondary.scalingFactor());
142 
143  // Calculate and cache the diagonal secondary-secondary entries
144  for (_i = 0; _i < secondarydof.size(); ++_i)
145  {
146  Number value = 0.0;
147  switch (_formulation)
148  {
149  case Moose::Penalty:
151  break;
152  case Moose::Kinematic:
153  value = -jacobian(secondarydof[_i], secondarydof[_i]) / primarydof.size() +
155  break;
156  }
158  _assembly, value, secondarydof[_i], secondarydof[_i], _var_secondary.scalingFactor());
159  }
160 }
161 
162 void
164 {
165 }
VarFieldType
Definition: MooseTypes.h:634
virtual void updateConnectivity()
Built the connectivity for this constraint.
virtual void zero() override final
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
virtual void zero() override final
unsigned int _j
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...
Definition: MooseUtils.h:1147
Base class for all Constraint types.
Definition: Constraint.h:19
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.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
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.
Definition: MooseTypes.h:627
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
virtual const OutputTools< Real >::VariableValue & value()
The value of the variable this object is operating on.
Assembly & _assembly
Reference to this Kernel&#39;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.
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
Real Number
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 &parameters)
MooseVariable & _var_secondary
static InputParameters validParams()
Definition: Constraint.C:15
void scalingFactor(const std::vector< Real > &factor)
Set the scaling factor for this variable.
MooseVariable & _var