www.mooseframework.org
NodalConstraint.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #include "NodalConstraint.h"
16 
17 // MOOSE includes
18 #include "Assembly.h"
19 #include "MooseVariable.h"
20 #include "SystemBase.h"
21 
22 #include "libmesh/sparse_matrix.h"
23 
24 template <>
27 {
29  MooseEnum formulationtype("penalty kinematic", "penalty");
30  params.addParam<MooseEnum>("formulation",
31  formulationtype,
32  "Formulation used to calculate constraint - penalty or kinematic.");
33  return params;
34 }
35 
37  : Constraint(parameters),
39  _u_slave(_var.nodalSlnNeighbor()),
40  _u_master(_var.nodalSln())
41 {
42  MooseEnum temp_formulation = getParam<MooseEnum>("formulation");
43  if (temp_formulation == "penalty")
45  else if (temp_formulation == "kinematic")
47  else
48  mooseError("Formulation must be either Penalty or Kinematic");
49 }
50 
51 void
52 NodalConstraint::computeResidual(NumericVector<Number> & residual)
53 {
54  if ((_weights.size() == 0) && (_master_node_vector.size() == 1))
55  _weights.push_back(1.0);
56 
57  std::vector<dof_id_type> masterdof = _var.dofIndices();
58  std::vector<dof_id_type> slavedof = _var.dofIndicesNeighbor();
59  DenseVector<Number> re(masterdof.size());
60  DenseVector<Number> neighbor_re(slavedof.size());
61 
62  re.zero();
63  neighbor_re.zero();
64 
65  for (_i = 0; _i < slavedof.size(); ++_i)
66  {
67  for (_j = 0; _j < masterdof.size(); ++_j)
68  {
69  switch (_formulation)
70  {
71  case Moose::Penalty:
73  neighbor_re(_i) += computeQpResidual(Moose::Slave) * _var.scalingFactor();
74  break;
75  case Moose::Kinematic:
76  // Transfer the current residual of the slave node to the master nodes
77  Real res = residual(slavedof[_i]);
78  re(_j) += res * _weights[_j];
79  neighbor_re(_i) += -res / _master_node_vector.size() + computeQpResidual(Moose::Slave);
80  break;
81  }
82  }
83  }
84  _assembly.cacheResidualNodes(re, masterdof);
85  _assembly.cacheResidualNodes(neighbor_re, slavedof);
86 }
87 
88 void
89 NodalConstraint::computeJacobian(SparseMatrix<Number> & jacobian)
90 {
91  if ((_weights.size() == 0) && (_master_node_vector.size() == 1))
92  _weights.push_back(1.0);
93 
94  // Calculate Jacobian enteries and cache those entries along with the row and column indices
95  std::vector<dof_id_type> slavedof = _var.dofIndicesNeighbor();
96  std::vector<dof_id_type> masterdof = _var.dofIndices();
97 
98  DenseMatrix<Number> Kee(masterdof.size(), masterdof.size());
99  DenseMatrix<Number> Ken(masterdof.size(), slavedof.size());
100  DenseMatrix<Number> Kne(slavedof.size(), masterdof.size());
101  DenseMatrix<Number> Knn(slavedof.size(), slavedof.size());
102 
103  Kee.zero();
104  Ken.zero();
105  Kne.zero();
106  Knn.zero();
107 
108  for (_i = 0; _i < slavedof.size(); ++_i)
109  {
110  for (_j = 0; _j < masterdof.size(); ++_j)
111  {
112  switch (_formulation)
113  {
114  case Moose::Penalty:
119  break;
120  case Moose::Kinematic:
121  Kee(_j, _j) = 0.;
122  Ken(_j, _i) += jacobian(slavedof[_i], masterdof[_j]) * _weights[_j];
123  Kne(_i, _j) += -jacobian(slavedof[_i], masterdof[_j]) / masterdof.size() +
125  Knn(_i, _i) += -jacobian(slavedof[_i], slavedof[_i]) / masterdof.size() +
127  break;
128  }
129  }
130  }
131  _assembly.cacheJacobianBlock(Kee, masterdof, masterdof, _var.scalingFactor());
132  _assembly.cacheJacobianBlock(Ken, masterdof, slavedof, _var.scalingFactor());
133  _assembly.cacheJacobianBlock(Kne, slavedof, masterdof, _var.scalingFactor());
134  _assembly.cacheJacobianBlock(Knn, slavedof, slavedof, _var.scalingFactor());
135 }
136 
137 void
139 {
140 }
virtual void updateConnectivity()
Built the connectivity for this constraint.
void cacheJacobianBlock(DenseMatrix< Number > &jac_block, std::vector< dof_id_type > &idof_indices, std::vector< dof_id_type > &jdof_indices, Real scaling_factor)
Definition: Assembly.C:1522
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:42
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int _i
Counter for master and slave nodes.
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.
virtual void computeResidual(NumericVector< Number > &residual)
Computes the nodal residual.
std::vector< dof_id_type > & dofIndicesNeighbor()
Get DOF indices for currently selected element.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:37
InputParameters validParams< NodalConstraint >()
Assembly & _assembly
Definition: Constraint.h:79
std::vector< dof_id_type > _master_node_vector
node IDs of the master node
std::vector< dof_id_type > & dofIndices()
void cacheResidualNodes(DenseVector< Number > &res, std::vector< dof_id_type > &dof_index)
Lets an external class cache residual at a set of nodes.
Definition: Assembly.C:1396
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...
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
InputParameters validParams< Constraint >()
Definition: Constraint.C:20
NodalConstraint(const InputParameters &parameters)
MooseVariable & _var
Definition: Constraint.h:80
void scalingFactor(Real factor)
Set the scaling factor for this variable.