www.mooseframework.org
TiedValueConstraint.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 "TiedValueConstraint.h"
16 
17 // MOOSE includes
18 #include "MooseVariable.h"
19 #include "SystemBase.h"
20 
21 #include "libmesh/sparse_matrix.h"
22 
23 template <>
26 {
28  params.addParam<Real>("scaling", 1, "scaling factor to be applied to constraint equations");
29  params.set<bool>("use_displaced_mesh") = true;
30  return params;
31 }
32 
34  : NodeFaceConstraint(parameters),
35  _scaling(getParam<Real>("scaling")),
36  _residual_copy(_sys.residualGhosted())
37 {
38 }
39 
40 Real
42 {
43  return _u_master[_qp];
44 }
45 
46 Real
48 {
49  Real scaling_factor = _var.scalingFactor();
50  Real slave_resid = 0;
51  Real retVal = 0;
52  switch (type)
53  {
54  case Moose::Slave:
55  retVal = (_u_slave[_qp] - _u_master[_qp]) * _test_slave[_i][_qp] * _scaling;
56  break;
57  case Moose::Master:
58  slave_resid = _residual_copy(_current_node->dof_number(0, _var.number(), 0)) / scaling_factor;
59  retVal = slave_resid * _test_master[_i][_qp];
60  break;
61  default:
62  break;
63  }
64  return retVal;
65 }
66 
67 Real
69 {
70  Real scaling_factor = _var.scalingFactor();
71  Real slave_jac = 0;
72  Real retVal = 0;
73  switch (type)
74  {
75  case Moose::SlaveSlave:
76  retVal = _phi_slave[_j][_qp] * _test_slave[_i][_qp] * _scaling;
77  break;
78  case Moose::SlaveMaster:
79  retVal = -_phi_master[_j][_qp] * _test_slave[_i][_qp] * _scaling;
80  break;
81  case Moose::MasterSlave:
82  slave_jac =
83  (*_jacobian)(_current_node->dof_number(0, _var.number(), 0), _connected_dof_indices[_j]);
84  retVal = slave_jac * _test_master[_i][_qp] / scaling_factor;
85  break;
87  retVal = 0;
88  break;
89  }
90  return retVal;
91 }
ConstraintType
Definition: MooseTypes.h:198
virtual Real computeQpResidual(Moose::ConstraintType type) override
This is the virtual that derived classes should override for computing the residual on neighboring el...
T & set(const std::string &name, bool quiet_mode=false)
Returns a writable reference to the named parameters.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
NumericVector< Number > & _residual_copy
const VariableValue & _u_master
Holds the current solution at the current quadrature point.
unsigned int _i
Definition: Constraint.h:83
const VariableTestValue & _test_master
Side test function.
unsigned int _j
Definition: Constraint.h:83
VariableTestValue _test_slave
Shape function on the slave side. This will always only have one entry and that entry will always be ...
InputParameters validParams< NodeFaceConstraint >()
A NodeFaceConstraint is used when you need to create constraints between two surfaces in a mesh...
VariablePhiValue _phi_slave
Shape function on the slave side. This will always.
const VariableValue & _u_slave
Value of the unknown variable this BC is action on.
const Node *& _current_node
current node being processed
InputParameters validParams< TiedValueConstraint >()
virtual Real computeQpSlaveValue() override
Compute the value the slave node should have at the beginning of a timestep.
MatType type
unsigned int number() const
Get variable number coming from libMesh.
const VariablePhiValue & _phi_master
Side shape function.
ConstraintJacobianType
Definition: MooseTypes.h:204
TiedValueConstraint(const InputParameters &parameters)
std::vector< dof_id_type > _connected_dof_indices
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...
unsigned int _qp
Definition: Constraint.h:84
MooseVariable & _var
Definition: Constraint.h:80
void scalingFactor(Real factor)
Set the scaling factor for this variable.
virtual Real computeQpJacobian(Moose::ConstraintJacobianType type) override
This is the virtual that derived classes should override for computing the Jacobian on neighboring el...