www.mooseframework.org
ElemElemConstraint.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 "ElemElemConstraint.h"
16 
17 // MOOSE includes
18 #include "Assembly.h"
19 #include "ElementPairInfo.h"
20 #include "FEProblem.h"
21 #include "MooseMesh.h"
22 #include "MooseVariable.h"
23 
24 #include "libmesh/quadrature.h"
25 
26 template <>
29 {
31  params.addParam<unsigned int>("interface_id", 0, "The id of the interface.");
32  return params;
33 }
34 
36  : Constraint(parameters),
38  _fe_problem(*parameters.get<FEProblemBase *>("_fe_problem_base")),
39  _dim(_mesh.dimension()),
40 
41  _current_elem(_assembly.elem()),
42 
43  _neighbor_elem(_assembly.neighbor()),
44 
45  _u(_var.sln()),
46  _grad_u(_var.gradSln()),
47 
48  _phi(_assembly.phi()),
49  _grad_phi(_assembly.gradPhi()),
50 
51  _test(_var.phi()),
52  _grad_test(_var.gradPhi()),
53 
54  _phi_neighbor(_assembly.phiNeighbor()),
55  _grad_phi_neighbor(_assembly.gradPhiNeighbor()),
56 
57  _test_neighbor(_var.phiNeighbor()),
58  _grad_test_neighbor(_var.gradPhiNeighbor()),
59 
60  _u_neighbor(_var.slnNeighbor()),
61  _grad_u_neighbor(_var.gradSlnNeighbor())
62 {
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  DenseVector<Number> & re = is_elem ? _assembly.residualBlock(_var.number())
96  for (_qp = 0; _qp < _constraint_q_point.size(); _qp++)
97  for (_i = 0; _i < test_space.size(); _i++)
99 }
100 
101 void
103 {
104  // Compute the residual for this element
106 
107  // Compute the residual for the neighbor
109 }
110 
111 void
113 {
114  const VariableTestValue & test_space =
116  const VariableTestValue & loc_phi =
118  DenseMatrix<Number> & Kxx =
119  type == Moose::ElementElement
121  : type == Moose::ElementNeighbor
124  : type == Moose::NeighborElement
129 
130  for (_qp = 0; _qp < _constraint_q_point.size(); _qp++)
131  for (_i = 0; _i < test_space.size(); _i++)
132  for (_j = 0; _j < loc_phi.size(); _j++)
133  Kxx(_i, _j) += _constraint_weight[_qp] * computeQpJacobian(type);
134 }
135 
136 void
138 {
139  // Compute element-element Jacobian
141 
142  // Compute element-neighbor Jacobian
144 
145  // Compute neighbor-element Jacobian
147 
148  // Compute neighbor-neighbor Jacobian
150 }
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...
virtual void computeJacobian()
Computes the jacobian for the current side.
std::vector< Point > _constraint_q_point
Quadrature points used in integration of constraint.
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...
DGResidualType
Definition: MooseTypes.h:184
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
DenseVector< Number > & residualBlock(unsigned int var_num, Moose::KernelType type=Moose::KT_NONTIME)
Definition: Assembly.h:504
ElemElemConstraint(const InputParameters &parameters)
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.
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:79
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:190
virtual void computeElemNeighResidual(Moose::DGResidualType type)
Computes the residual for this element or the neighbor.
MatType type
unsigned int number() const
Get variable number coming from libMesh.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
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.
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...
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:887
InputParameters validParams< Constraint >()
Definition: Constraint.C:20
DenseVector< Number > & residualBlockNeighbor(unsigned int var_num, Moose::KernelType type=Moose::KT_NONTIME)
Definition: Assembly.h:509
DenseMatrix< Number > & jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:901
unsigned int _qp
Definition: Constraint.h:84
MooseVariable & _var
Definition: Constraint.h:80
unsigned int _i
Indices for looping over DOFs.