www.mooseframework.org
DGKernel.h
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 #ifndef DGKERNEL_H
16 #define DGKERNEL_H
17 
18 // local includes
19 #include "MooseArray.h"
20 #include "MooseObject.h"
21 #include "BlockRestrictable.h"
22 #include "BoundaryRestrictable.h"
23 #include "SetupInterface.h"
24 #include "TransientInterface.h"
25 #include "UserObjectInterface.h"
27 #include "FunctionInterface.h"
29 #include "Restartable.h"
30 #include "ZeroInterface.h"
31 #include "MeshChangedInterface.h"
32 
33 // Forward Declarations
34 class MooseMesh;
35 class SubProblem;
36 class Assembly;
37 
38 class DGKernel;
39 
40 template <>
42 
47 class DGKernel : public MooseObject,
48  public BlockRestrictable,
49  public BoundaryRestrictable,
50  public SetupInterface,
51  public TransientInterface,
52  public FunctionInterface,
53  public UserObjectInterface,
56  public Restartable,
57  public ZeroInterface,
59 {
60 public:
69 
70  virtual ~DGKernel();
71 
76 
81 
86 
90  virtual void computeResidual();
91 
96 
100  virtual void computeJacobian();
101 
105  virtual void computeOffDiagElemNeighJacobian(Moose::DGJacobianType type, unsigned int jvar);
106 
110  virtual void computeOffDiagJacobian(unsigned int jvar);
111 
112 protected:
115 
117 
121  // unsigned int _dim;
122 
123  const Elem *& _current_elem;
124 
126  const Real & _current_elem_volume;
127 
129  const Elem *& _neighbor_elem;
130 
132  unsigned int & _current_side;
134  const Elem *& _current_side_elem;
135 
137  const Real & _current_side_volume;
138 
141  unsigned int _qp;
143  QBase *& _qrule;
146 
147  unsigned int _i, _j;
148 
150 
152  const VariableValue & _u;
153 
156  // shape functions
159  // test functions
160 
167 
172 
177 
182 
184  DenseVector<Number> _local_re;
185 
187  DenseMatrix<Number> _local_kxx;
188 
191  std::vector<MooseVariable *> _save_in;
192  std::vector<AuxVariableName> _save_in_strings;
193 
196  std::vector<MooseVariable *> _diag_save_in;
197  std::vector<AuxVariableName> _diag_save_in_strings;
198 
203  virtual Real computeQpResidual(Moose::DGResidualType type) = 0;
204 
209  virtual Real computeQpJacobian(Moose::DGJacobianType type) = 0;
210 
214  virtual Real computeQpOffDiagJacobian(Moose::DGJacobianType type, unsigned int jvar);
215 
217  const Real & getNeighborElemVolume();
218 
219 public:
220  // boundary id used for internal edges (all DG kernels lives on this boundary id -- a made-up
221  // number)
222  static const BoundaryID InternalBndId;
223 
224 protected:
225  static Threads::spin_mutex _resid_vars_mutex;
226  static Threads::spin_mutex _jacoby_vars_mutex;
227 };
228 
229 #endif // DGKERNEL_H
const VariablePhiGradient & _grad_phi_neighbor
Gradient of side shape function.
Definition: DGKernel.h:171
std::vector< AuxVariableName > _diag_save_in_strings
Definition: DGKernel.h:197
virtual void computeOffDiagElemNeighJacobian(Moose::DGJacobianType type, unsigned int jvar)
Computes the element-element off-diagonal Jacobian.
Definition: DGKernel.C:281
const VariableTestGradient & _grad_test
Gradient of side shape function.
Definition: DGKernel.h:164
A class for creating restricted objects.
Definition: Restartable.h:31
THREAD_ID _tid
Definition: DGKernel.h:116
DenseMatrix< Number > _local_kxx
Holds residual entries as they are accumulated by this DGKernel.
Definition: DGKernel.h:187
Keeps track of stuff related to assembling.
Definition: Assembly.h:63
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
Class for stuff related to variables.
Definition: MooseVariable.h:43
const VariableValue & _u_neighbor
Holds the current solution at the current quadrature point.
Definition: DGKernel.h:179
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: DGKernel.h:190
InputParameters validParams< DGKernel >()
Definition: DGKernel.C:33
SubProblem & subProblem()
Return a reference to the subproblem.
Definition: DGKernel.C:337
static const BoundaryID InternalBndId
Definition: DGKernel.h:222
BoundaryID _boundary_id
Definition: DGKernel.h:149
std::vector< AuxVariableName > _save_in_strings
Definition: DGKernel.h:192
Assembly & _assembly
Definition: DGKernel.h:118
const MooseArray< Real > & _coord
Definition: DGKernel.h:145
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: DGKernel.h:195
const VariablePhiValue & _phi_neighbor
Side shape function.
Definition: DGKernel.h:169
virtual void computeResidual()
Computes the residual for the current side.
Definition: DGKernel.C:210
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual Real computeQpJacobian(Moose::DGJacobianType type)=0
This is the virtual that derived classes should override for computing the Jacobian on neighboring el...
const Elem *& _neighbor_elem
The neighboring element.
Definition: DGKernel.h:129
/class BoundaryRestrictable /brief Provides functionality for limiting the object to certain boundary...
const MooseArray< Point > & _q_point
Definition: DGKernel.h:142
const VariableTestGradient & _grad_test_neighbor
Gradient of side shape function.
Definition: DGKernel.h:176
const Real & _current_elem_volume
The volume (or length) of the current element.
Definition: DGKernel.h:126
unsigned int & _current_side
Current side.
Definition: DGKernel.h:132
DGResidualType
Definition: MooseTypes.h:184
Base class for a system (of equations)
Definition: SystemBase.h:91
const Elem *& _current_side_elem
Current side element.
Definition: DGKernel.h:134
const VariablePhiGradient & _grad_phi
Definition: DGKernel.h:158
MooseMesh & _mesh
Definition: DGKernel.h:120
MooseVariable & _var
Definition: DGKernel.h:119
const VariableTestValue & _test_neighbor
Side test function.
Definition: DGKernel.h:174
const Moose::CoordinateSystemType & _coord_sys
Coordinate system.
Definition: DGKernel.h:140
const Real & _current_side_volume
The volume (or length) of the current side.
Definition: DGKernel.h:137
unsigned int _j
Definition: DGKernel.h:147
Interface for objects that needs transient capabilities.
const VariableGradient & _grad_u_neighbor
Holds the current solution gradient at the current quadrature point.
Definition: DGKernel.h:181
virtual Real computeQpOffDiagJacobian(Moose::DGJacobianType type, unsigned int jvar)
This is the virtual that derived classes should override for computing the off-diag Jacobian...
Definition: DGKernel.C:325
Interface for notifications that the mesh has changed.
The DGKernel class is responsible for calculating the residuals for various physics on internal sides...
Definition: DGKernel.h:47
Every object that can be built by the factory should be derived from this class.
Definition: MooseObject.h:36
static Threads::spin_mutex _resid_vars_mutex
Definition: DGKernel.h:225
virtual void computeJacobian()
Computes the jacobian for the current side.
Definition: DGKernel.C:265
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
DenseVector< Number > _local_re
Holds residual entries as they are accumulated by this DGKernel.
Definition: DGKernel.h:184
DGKernel(const InputParameters &parameters)
Factory constructor initializes all internal references needed for residual computation.
Definition: DGKernel.C:72
Interface for objects that need to use UserObjects.
CoordinateSystemType
Definition: MooseTypes.h:212
SubProblem & _subproblem
Definition: DGKernel.h:113
const VariableTestValue & _test
Side shape function.
Definition: DGKernel.h:162
DGJacobianType
Definition: MooseTypes.h:190
const MooseArray< Real > & _JxW
Definition: DGKernel.h:144
MatType type
unsigned int _i
Definition: DGKernel.h:147
std::vector< MooseVariable * > _diag_save_in
Definition: DGKernel.h:196
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:53
const VariableGradient & _grad_u
Holds the current solution gradient at the current quadrature point on the face.
Definition: DGKernel.h:155
Interface to bring zero values inside objects.
Definition: ZeroInterface.h:35
virtual void computeOffDiagJacobian(unsigned int jvar)
Computes d-residual / d-jvar...
Definition: DGKernel.C:304
An interface that restricts an object to subdomains via the &#39;blocks&#39; input parameter.
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseObject.h:53
const MooseArray< Point > & _normals
Normal vectors at the quadrature points.
Definition: DGKernel.h:166
SystemBase & _sys
Definition: DGKernel.h:114
virtual void computeElemNeighJacobian(Moose::DGJacobianType type)
Computes the element/neighbor-element/neighbor Jacobian.
Definition: DGKernel.C:220
const Elem *& _current_elem
Definition: DGKernel.h:123
const VariablePhiValue & _phi
Definition: DGKernel.h:157
virtual void computeElemNeighResidual(Moose::DGResidualType type)
Computes the residual for this element or the neighbor.
Definition: DGKernel.C:177
std::vector< MooseVariable * > _save_in
Definition: DGKernel.h:191
virtual Real computeQpResidual(Moose::DGResidualType type)=0
This is the virtual that derived classes should override for computing the residual on neighboring el...
const Real & getNeighborElemVolume()
The volume (or length) of the current neighbor.
Definition: DGKernel.C:343
const VariableValue & _u
Holds the current solution at the current quadrature point on the face.
Definition: DGKernel.h:152
MooseVariable & variable()
The variable number that this kernel operates on.
Definition: DGKernel.C:331
Interface for objects that need to use functions.
virtual ~DGKernel()
Definition: DGKernel.C:174
QBase *& _qrule
Definition: DGKernel.h:143
unsigned int _qp
Definition: DGKernel.h:141
boundary_id_type BoundaryID
Definition: MooseTypes.h:75
unsigned int THREAD_ID
Definition: MooseTypes.h:79
static Threads::spin_mutex _jacoby_vars_mutex
Definition: DGKernel.h:226