www.mooseframework.org
ElemElemConstraint.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 "ElemElemConstraint.h"
11 
12 // MOOSE includes
13 #include "Assembly.h"
14 #include "ElementPairInfo.h"
15 #include "FEProblem.h"
16 #include "MooseMesh.h"
17 #include "MooseVariableFE.h"
18 
19 #include "libmesh/quadrature.h"
20 
21 template <>
24 {
26  params.addParam<unsigned int>("interface_id", 0, "The id of the interface.");
27  params.addRequiredParam<NonlinearVariableName>(
28  "variable", "The name of the variable that this constraint is applied to.");
29  return params;
30 }
31 
33  : Constraint(parameters),
37  _fe_problem(*getCheckedPointerParam<FEProblemBase *>("_fe_problem_base")),
38  _dim(_mesh.dimension()),
39  _interface_id(getParam<unsigned int>("interface_id")),
40  _var(_sys.getFieldVariable<Real>(_tid, parameters.get<NonlinearVariableName>("variable"))),
41 
42  _current_elem(_assembly.elem()),
43 
44  _neighbor_elem(_assembly.neighbor()),
45 
46  _u(_var.sln()),
47  _grad_u(_var.gradSln()),
48 
49  _phi(_assembly.phi(_var)),
50  _grad_phi(_assembly.gradPhi(_var)),
51 
52  _test(_var.phi()),
53  _grad_test(_var.gradPhi()),
54 
55  _phi_neighbor(_assembly.phiNeighbor(_var)),
56  _grad_phi_neighbor(_assembly.gradPhiNeighbor(_var)),
57 
58  _test_neighbor(_var.phiNeighbor()),
59  _grad_test_neighbor(_var.gradPhiNeighbor()),
60 
61  _u_neighbor(_var.slnNeighbor()),
62  _grad_u_neighbor(_var.gradSlnNeighbor())
63 {
65 }
66 
67 void
68 ElemElemConstraint::reinit(const ElementPairInfo & element_pair_info)
69 {
70  reinitConstraintQuadrature(element_pair_info);
71 }
72 
73 void
75 {
76  _constraint_q_point.resize(element_pair_info._elem1_constraint_q_point.size());
77  _constraint_weight.resize(element_pair_info._elem1_constraint_JxW.size());
78  std::copy(element_pair_info._elem1_constraint_q_point.begin(),
79  element_pair_info._elem1_constraint_q_point.end(),
80  _constraint_q_point.begin());
81  std::copy(element_pair_info._elem1_constraint_JxW.begin(),
82  element_pair_info._elem1_constraint_JxW.end(),
83  _constraint_weight.begin());
84 }
85 
86 void
88 {
89  bool is_elem;
90  if (type == Moose::Element)
91  is_elem = true;
92  else
93  is_elem = false;
94 
95  const VariableTestValue & test_space = is_elem ? _test : _test_neighbor;
96  DenseVector<Number> & re = is_elem ? _assembly.residualBlock(_var.number())
98  for (_qp = 0; _qp < _constraint_q_point.size(); _qp++)
99  for (_i = 0; _i < test_space.size(); _i++)
101 }
102 
103 void
105 {
106  // Compute the residual for this element
108 
109  // Compute the residual for the neighbor
111 }
112 
113 void
115 {
116  const VariableTestValue & test_space =
118  const VariableTestValue & loc_phi =
120  DenseMatrix<Number> & Kxx =
131 
132  for (_qp = 0; _qp < _constraint_q_point.size(); _qp++)
133  for (_i = 0; _i < test_space.size(); _i++)
134  for (_j = 0; _j < loc_phi.size(); _j++)
136 }
137 
138 void
140 {
141  // Compute element-element Jacobian
143 
144  // Compute element-neighbor Jacobian
146 
147  // Compute neighbor-element Jacobian
149 
150  // Compute neighbor-neighbor Jacobian
152 }
VarFieldType
Definition: MooseTypes.h:488
MooseVariable & _var
std::vector< Real > _constraint_weight
Weights of quadrature points used in integration of constraint.
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
DenseVector< Number > & residualBlockNeighbor(unsigned int var_num, TagID tag_id=0)
Definition: Assembly.h:725
unsigned int number() const
Get variable number coming from libMesh.
virtual void computeJacobian()
Computes the jacobian for the current side.
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar, TagID tag=0)
Definition: Assembly.C:2019
std::vector< Point > _constraint_q_point
Quadrature points used in integration of constraint.
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...
DGResidualType
Definition: MooseTypes.h:509
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
ElemElemConstraint(const InputParameters &parameters)
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...
const std::string & type() const
Get the type of this object.
Definition: MooseObject.h:53
Enhances MooseVariableInterface interface provide values from neighbor elements.
const VariableTestValue & _test_neighbor
Neighbor test function.
virtual void computeElemNeighJacobian(Moose::DGJacobianType type)
Computes the element/neighbor-element/neighbor Jacobian.
std::vector< Point > _elem1_constraint_q_point
virtual Real computeQpJacobian(Moose::DGJacobianType type)=0
Compute the Jacobian for one of the constraint quadrature points.
OutputTools< Real >::VariableTestValue VariableTestValue
Definition: MooseTypes.h:205
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
virtual void reinitConstraintQuadrature(const ElementPairInfo &element_pair_info)
Set information needed for constraint integration.
const VariableTestValue & _test
Test function.
This is the ElementPairInfo class.
const VariablePhiValue & _phi_neighbor
Neighbor shape function.
Assembly & _assembly
Definition: Constraint.h:72
virtual Real computeQpResidual(Moose::DGResidualType type)=0
Compute the residual for one of the constraint quadrature points.
virtual void computeResidual()
Computes the residual for the current side.
DGJacobianType
Definition: MooseTypes.h:515
virtual void computeElemNeighResidual(Moose::DGResidualType type)
Computes the residual for this element or the neighbor.
MatType type
virtual void reinit(const ElementPairInfo &element_pair_info)
reinit element-element constraint
std::vector< Real > _elem1_constraint_JxW
InputParameters validParams< ElemElemConstraint >()
const VariablePhiValue & _phi
Shape function.
DenseMatrix< Number > & jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar, TagID tag=0)
Definition: Assembly.C:2033
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
DenseVector< Number > & residualBlock(unsigned int var_num, TagID tag_id=0)
Definition: Assembly.h:720
InputParameters validParams< Constraint >()
Definition: Constraint.C:16
unsigned int _qp
Definition: Constraint.h:76
unsigned int _i
Indices for looping over DOFs.