www.mooseframework.org
DGConvection.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 "DGConvection.h"
16 
17 template <>
20 {
22  params.addRequiredParam<RealVectorValue>("velocity", "Velocity vector");
23  params.addClassDescription("DG upwinding for the convection");
24  return params;
25 }
26 
28  : DGKernel(parameters), _velocity(getParam<RealVectorValue>("velocity"))
29 {
30 }
31 
32 Real
34 {
35  Real r = 0;
36 
37  Real vdotn = _velocity * _normals[_qp];
38 
39  switch (type)
40  {
41  case Moose::Element:
42  if (vdotn >= 0)
43  r += vdotn * _u[_qp] * _test[_i][_qp];
44  else
45  r += vdotn * _u_neighbor[_qp] * _test[_i][_qp];
46  break;
47 
48  case Moose::Neighbor:
49  if (vdotn >= 0)
50  r -= vdotn * _u[_qp] * _test_neighbor[_i][_qp];
51  else
52  r -= vdotn * _u_neighbor[_qp] * _test_neighbor[_i][_qp];
53  break;
54  }
55 
56  return r;
57 }
58 
59 Real
61 {
62  Real r = 0;
63 
64  Real vdotn = _velocity * _normals[_qp];
65 
66  switch (type)
67  {
69  if (vdotn >= 0)
70  r += vdotn * _phi[_j][_qp] * _test[_i][_qp];
71  break;
72 
74  if (vdotn < 0)
75  r += vdotn * _phi_neighbor[_j][_qp] * _test[_i][_qp];
76  break;
77 
79  if (vdotn >= 0)
80  r -= vdotn * _phi[_j][_qp] * _test_neighbor[_i][_qp];
81  break;
82 
84  if (vdotn < 0)
85  r -= vdotn * _phi_neighbor[_j][_qp] * _test_neighbor[_i][_qp];
86  break;
87  }
88 
89  return r;
90 }
VectorValue< Real > RealVectorValue
Definition: Assembly.h:40
const VariableValue & _u_neighbor
Holds the current solution at the current quadrature point.
Definition: DGKernel.h:179
InputParameters validParams< DGKernel >()
Definition: DGKernel.C:33
const VariablePhiValue & _phi_neighbor
Side shape function.
Definition: DGKernel.h:169
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
DGResidualType
Definition: MooseTypes.h:184
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...
RealVectorValue _velocity
Definition: DGConvection.h:34
DGConvection(const InputParameters &parameters)
Definition: DGConvection.C:27
const VariableTestValue & _test_neighbor
Side test function.
Definition: DGKernel.h:174
unsigned int _j
Definition: DGKernel.h:147
The DGKernel class is responsible for calculating the residuals for various physics on internal sides...
Definition: DGKernel.h:47
virtual Real computeQpResidual(Moose::DGResidualType type)
This is the virtual that derived classes should override for computing the residual on neighboring el...
Definition: DGConvection.C:33
const VariableTestValue & _test
Side shape function.
Definition: DGKernel.h:162
DGJacobianType
Definition: MooseTypes.h:190
MatType type
unsigned int _i
Definition: DGKernel.h:147
InputParameters validParams< DGConvection >()
Definition: DGConvection.C:19
const MooseArray< Point > & _normals
Normal vectors at the quadrature points.
Definition: DGKernel.h:166
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
const VariablePhiValue & _phi
Definition: DGKernel.h:157
const VariableValue & _u
Holds the current solution at the current quadrature point on the face.
Definition: DGKernel.h:152
unsigned int _qp
Definition: DGKernel.h:141
virtual Real computeQpJacobian(Moose::DGJacobianType type)
This is the virtual that derived classes should override for computing the Jacobian on neighboring el...
Definition: DGConvection.C:60