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 #include "SystemBase.h"
19 
20 #include "libmesh/quadrature.h"
21 
24 {
26  params.addParam<unsigned int>("interface_id", 0, "The id of the interface.");
27  return params;
28 }
29 
31  : Constraint(parameters),
35  _fe_problem(*getCheckedPointerParam<FEProblemBase *>("_fe_problem_base")),
36  _dim(_mesh.dimension()),
37  _interface_id(getParam<unsigned int>("interface_id")),
38  _var(_sys.getFieldVariable<Real>(_tid, parameters.get<NonlinearVariableName>("variable"))),
39 
40  _current_elem(_assembly.elem()),
41 
42  _neighbor_elem(_assembly.neighbor()),
43 
44  _u(_var.sln()),
45  _grad_u(_var.gradSln()),
46 
47  _phi(_assembly.phi(_var)),
48  _grad_phi(_assembly.gradPhi(_var)),
49 
50  _test(_var.phi()),
51  _grad_test(_var.gradPhi()),
52 
53  _phi_neighbor(_assembly.phiNeighbor(_var)),
54  _grad_phi_neighbor(_assembly.gradPhiNeighbor(_var)),
55 
56  _test_neighbor(_var.phiNeighbor()),
57  _grad_test_neighbor(_var.gradPhiNeighbor()),
58 
59  _u_neighbor(_var.slnNeighbor()),
60  _grad_u_neighbor(_var.gradSlnNeighbor())
61 {
63 }
64 
65 void
66 ElemElemConstraint::reinit(const ElementPairInfo & element_pair_info)
67 {
68  reinitConstraintQuadrature(element_pair_info);
69 }
70 
71 void
73 {
74  _constraint_q_point.resize(element_pair_info._elem1_constraint_q_point.size());
75  _constraint_weight.resize(element_pair_info._elem1_constraint_JxW.size());
76  std::copy(element_pair_info._elem1_constraint_q_point.begin(),
77  element_pair_info._elem1_constraint_q_point.end(),
78  _constraint_q_point.begin());
79  std::copy(element_pair_info._elem1_constraint_JxW.begin(),
80  element_pair_info._elem1_constraint_JxW.end(),
81  _constraint_weight.begin());
82 }
83 
84 void
86 {
87  bool is_elem;
88  if (type == Moose::Element)
89  is_elem = true;
90  else
91  is_elem = false;
92 
93  const VariableTestValue & test_space = is_elem ? _test : _test_neighbor;
94  if (is_elem)
96  else
98  for (_qp = 0; _qp < _constraint_q_point.size(); _qp++)
99  for (_i = 0; _i < test_space.size(); _i++)
102 }
103 
104 void
106 {
107  // Compute the residual for this element
109 
110  // Compute the residual for the neighbor
112 }
113 
114 void
116 {
117  const VariableTestValue & test_space =
119  const VariableTestValue & loc_phi =
122 
123  for (_qp = 0; _qp < _constraint_q_point.size(); _qp++)
124  for (_i = 0; _i < test_space.size(); _i++)
125  for (_j = 0; _j < loc_phi.size(); _j++)
127 
129 }
130 
131 void
133 {
134  // Compute element-element Jacobian
136 
137  // Compute element-neighbor Jacobian
139 
140  // Compute neighbor-element Jacobian
142 
143  // Compute neighbor-neighbor Jacobian
145 }
VarFieldType
Definition: MooseTypes.h:634
virtual void computeJacobian() override
Computes the jacobian for the current side.
MooseVariable & _var
std::vector< Real > _constraint_weight
Weights of quadrature points used in integration of constraint.
virtual void computeResidual() override
Computes the residual for the current side.
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
unsigned int number() const
Get variable number coming from libMesh.
std::vector< Point > _constraint_q_point
Quadrature points used in integration of constraint.
T * get(const std::unique_ptr< T > &u)
The MooseUtils::get() specializations are used to support making forwards-compatible code changes fro...
Definition: MooseUtils.h:1147
Base class for all Constraint types.
Definition: Constraint.h:19
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
DGResidualType
Definition: MooseTypes.h:656
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
ElemElemConstraint(const InputParameters &parameters)
void prepareMatrixTagNeighbor(Assembly &assembly, unsigned int ivar, unsigned int jvar, Moose::DGJacobianType type)
Prepare data for computing element jacobian according to the active tags for DG and interface kernels...
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.
static InputParameters validParams()
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:312
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:627
virtual void reinitConstraintQuadrature(const ElementPairInfo &element_pair_info)
Set information needed for constraint integration.
const VariableTestValue & _test
Test function.
const std::string & type() const
Get the type of this class.
Definition: MooseBase.h:51
This is the ElementPairInfo class.
const VariablePhiValue & _phi_neighbor
Neighbor shape function.
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
virtual Real computeQpResidual(Moose::DGResidualType type)=0
Compute the residual for one of the constraint quadrature points.
DGJacobianType
Definition: MooseTypes.h:662
virtual void computeElemNeighResidual(Moose::DGResidualType type)
Computes the residual for this element or the neighbor.
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
void prepareVectorTagNeighbor(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual the according to active tags for DG and interface kernels...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void reinit(const ElementPairInfo &element_pair_info)
reinit element-element constraint
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
std::vector< Real > _elem1_constraint_JxW
const VariablePhiValue & _phi
Shape function.
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...
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual according to active tags.
static InputParameters validParams()
Definition: Constraint.C:15
unsigned int _qp
Definition: Constraint.h:36
void ErrorVector unsigned int
unsigned int _i
Indices for looping over DOFs.