www.mooseframework.org
TiedValueConstraint.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 "TiedValueConstraint.h"
11 
12 // MOOSE includes
13 #include "MooseVariableFE.h"
14 #include "SystemBase.h"
15 
16 #include "libmesh/sparse_matrix.h"
17 
19 
20 template <>
23 {
25  params.addParam<Real>("scaling", 1, "scaling factor to be applied to constraint equations");
26  params.set<bool>("use_displaced_mesh") = true;
27  return params;
28 }
29 
31  : NodeFaceConstraint(parameters),
32  _scaling(getParam<Real>("scaling")),
33  _residual_copy(_sys.residualGhosted())
34 {
35 }
36 
37 Real
39 {
40  return _u_master[_qp];
41 }
42 
43 Real
45 {
46  Real scaling_factor = _var.scalingFactor();
47  Real slave_resid = 0;
48  Real retVal = 0;
49  switch (type)
50  {
51  case Moose::Slave:
52  retVal = (_u_slave[_qp] - _u_master[_qp]) * _test_slave[_i][_qp] * _scaling;
53  break;
54  case Moose::Master:
55  slave_resid = _residual_copy(_current_node->dof_number(0, _var.number(), 0)) / scaling_factor;
56  retVal = slave_resid * _test_master[_i][_qp];
57  break;
58  default:
59  break;
60  }
61  return retVal;
62 }
63 
64 Real
66 {
67  Real scaling_factor = _var.scalingFactor();
68  Real slave_jac = 0;
69  Real retVal = 0;
70  switch (type)
71  {
72  case Moose::SlaveSlave:
73  retVal = _phi_slave[_j][_qp] * _test_slave[_i][_qp] * _scaling;
74  break;
75  case Moose::SlaveMaster:
76  retVal = -_phi_master[_j][_qp] * _test_slave[_i][_qp] * _scaling;
77  break;
78  case Moose::MasterSlave:
79  slave_jac =
80  (*_jacobian)(_current_node->dof_number(0, _var.number(), 0), _connected_dof_indices[_j]);
81  retVal = slave_jac * _test_master[_i][_qp] / scaling_factor;
82  break;
84  retVal = 0;
85  break;
86  default:
87  mooseError("Unsupported type");
88  break;
89  }
90  return retVal;
91 }
registerMooseObject("MooseApp", TiedValueConstraint)
ConstraintType
Definition: MooseTypes.h:523
MooseVariable & _var
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpResidual(Moose::ConstraintType type) override
This is the virtual that derived classes should override for computing the residual on neighboring el...
A TiedValueConstraint forces the value of a variable to be the same on both sides of an interface...
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:75
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
const VariableTestValue & _test_master
Side test function.
const Node *const & _current_node
current node being processed
const std::string & type() const
Get the type of this object.
Definition: MooseObject.h:53
unsigned int _j
Definition: Constraint.h:75
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.
InputParameters validParams< TiedValueConstraint >()
virtual Real computeQpSlaveValue() override
Compute the value the slave node should have at the beginning of a timestep.
MatType type
const VariablePhiValue & _phi_master
Side shape function.
ConstraintJacobianType
Definition: MooseTypes.h:543
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:76
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...