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 
19 template <>
22 {
24  MooseEnum formulationtype("penalty kinematic", "penalty");
25  params.addParam<MooseEnum>("formulation",
26  formulationtype,
27  "Formulation used to calculate constraint - penalty or kinematic.");
28  params.addRequiredParam<NonlinearVariableName>(
29  "variable", "The name of the variable that this constraint is applied to.");
30  return params;
31 }
32 
34  : Constraint(parameters),
38  _var(_sys.getFieldVariable<Real>(_tid, parameters.get<NonlinearVariableName>("variable"))),
39  _u_slave(_var.dofValuesNeighbor()),
40  _u_master(_var.dofValues())
41 {
43 
44  MooseEnum temp_formulation = getParam<MooseEnum>("formulation");
45  if (temp_formulation == "penalty")
47  else if (temp_formulation == "kinematic")
49  else
50  mooseError("Formulation must be either Penalty or Kinematic");
51 }
52 
53 void
54 NodalConstraint::computeResidual(NumericVector<Number> & residual)
55 {
56  if ((_weights.size() == 0) && (_master_node_vector.size() == 1))
57  _weights.push_back(1.0);
58 
59  std::vector<dof_id_type> masterdof = _var.dofIndices();
60  std::vector<dof_id_type> slavedof = _var.dofIndicesNeighbor();
61  DenseVector<Number> re(masterdof.size());
62  DenseVector<Number> neighbor_re(slavedof.size());
63 
64  re.zero();
65  neighbor_re.zero();
66 
67  for (_i = 0; _i < slavedof.size(); ++_i)
68  {
69  for (_j = 0; _j < masterdof.size(); ++_j)
70  {
71  switch (_formulation)
72  {
73  case Moose::Penalty:
75  neighbor_re(_i) += computeQpResidual(Moose::Slave) * _var.scalingFactor();
76  break;
77  case Moose::Kinematic:
78  // Transfer the current residual of the slave node to the master nodes
79  Real res = residual(slavedof[_i]);
80  re(_j) += res * _weights[_j];
81  neighbor_re(_i) += -res / _master_node_vector.size() + computeQpResidual(Moose::Slave);
82  break;
83  }
84  }
85  }
86  _assembly.cacheResidualNodes(re, masterdof);
87  _assembly.cacheResidualNodes(neighbor_re, slavedof);
88 }
89 
90 void
91 NodalConstraint::computeJacobian(SparseMatrix<Number> & jacobian)
92 {
93  if ((_weights.size() == 0) && (_master_node_vector.size() == 1))
94  _weights.push_back(1.0);
95 
96  // Calculate Jacobian enteries and cache those entries along with the row and column indices
97  std::vector<dof_id_type> slavedof = _var.dofIndicesNeighbor();
98  std::vector<dof_id_type> masterdof = _var.dofIndices();
99 
100  DenseMatrix<Number> Kee(masterdof.size(), masterdof.size());
101  DenseMatrix<Number> Ken(masterdof.size(), slavedof.size());
102  DenseMatrix<Number> Kne(slavedof.size(), masterdof.size());
103  DenseMatrix<Number> Knn(slavedof.size(), slavedof.size());
104 
105  Kee.zero();
106  Ken.zero();
107  Kne.zero();
108  Knn.zero();
109 
110  for (_i = 0; _i < slavedof.size(); ++_i)
111  {
112  for (_j = 0; _j < masterdof.size(); ++_j)
113  {
114  switch (_formulation)
115  {
116  case Moose::Penalty:
121  break;
122  case Moose::Kinematic:
123  Kee(_j, _j) = 0.;
124  Ken(_j, _i) += jacobian(slavedof[_i], masterdof[_j]) * _weights[_j];
125  Kne(_i, _j) += -jacobian(slavedof[_i], masterdof[_j]) / masterdof.size() +
127  Knn(_i, _i) += -jacobian(slavedof[_i], slavedof[_i]) / masterdof.size() +
129  break;
130  }
131  }
132  }
133  _assembly.cacheJacobianBlock(Kee, masterdof, masterdof, _var.scalingFactor());
134  _assembly.cacheJacobianBlock(Ken, masterdof, slavedof, _var.scalingFactor());
135  _assembly.cacheJacobianBlock(Kne, slavedof, masterdof, _var.scalingFactor());
136  _assembly.cacheJacobianBlock(Knn, slavedof, slavedof, _var.scalingFactor());
137 }
138 
139 void
141 {
142 }
VarFieldType
Definition: MooseTypes.h:488
virtual void updateConnectivity()
Built the connectivity for this constraint.
void cacheResidualNodes(const DenseVector< Number > &res, const std::vector< dof_id_type > &dof_index, TagID tag=0)
Lets an external class cache residual at a set of nodes.
Definition: Assembly.C:2875
virtual void computeJacobian(SparseMatrix< Number > &jacobian)
Computes the jacobian for the current element.
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
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...
Base class for all Constraint types.
Definition: Constraint.h:39
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
unsigned int _i
Counter for master and slave nodes.
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
void addMooseVariableDependency(MooseVariableFEBase *var)
Call this function to add the passed in MooseVariableFEBase as a variable that this object depends on...
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 slave node is constrained to move as a linear combination of the master nodes, the coefficients associated with each master node is stored in _weights.
Moose::ConstraintFormulationType _formulation
Specifies formulation type used to apply constraints.
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
virtual void computeResidual(NumericVector< Number > &residual)
Computes the nodal residual.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
InputParameters validParams< NodalConstraint >()
Assembly & _assembly
Definition: Constraint.h:72
void cacheJacobianBlock(DenseMatrix< Number > &jac_block, const std::vector< dof_id_type > &idof_indices, const std::vector< dof_id_type > &jdof_indices, Real scaling_factor, TagID tag=0)
Definition: Assembly.C:3003
std::vector< dof_id_type > _master_node_vector
node IDs of the master node
const std::vector< dof_id_type > & dofIndicesNeighbor() const final
Get neighbor DOF indices for currently selected element.
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...
Definition: Moose.h:112
InputParameters validParams< Constraint >()
Definition: Constraint.C:16
NodalConstraint(const InputParameters &parameters)
MooseVariable & _var
void scalingFactor(Real factor)
Set the scaling factor for this variable.