www.mooseframework.org
DGDiffusion.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 "DGDiffusion.h"
16 
17 // MOOSE includes
18 #include "MooseVariable.h"
19 
20 #include "libmesh/utility.h"
21 
22 template <>
25 {
27  // See header file for sigma and epsilon
28  params.addRequiredParam<Real>("sigma", "sigma");
29  params.addRequiredParam<Real>("epsilon", "epsilon");
30  params.addParam<MaterialPropertyName>(
31  "diff", 1., "The diffusion (or thermal conductivity or viscosity) coefficient.");
32  return params;
33 }
34 
36  : DGKernel(parameters),
37  _epsilon(getParam<Real>("epsilon")),
38  _sigma(getParam<Real>("sigma")),
39  _diff(getMaterialProperty<Real>("diff")),
40  _diff_neighbor(getNeighborMaterialProperty<Real>("diff"))
41 {
42 }
43 
44 Real
46 {
47  Real r = 0;
48 
49  const unsigned int elem_b_order = _var.order();
50  const double h_elem =
51  _current_elem->volume() / _current_side_elem->volume() * 1. / Utility::pow<2>(elem_b_order);
52 
53  switch (type)
54  {
55  case Moose::Element:
56  r -= 0.5 * (_diff[_qp] * _grad_u[_qp] * _normals[_qp] +
58  _test[_i][_qp];
59  r += _epsilon * 0.5 * (_u[_qp] - _u_neighbor[_qp]) * _diff[_qp] * _grad_test[_i][_qp] *
60  _normals[_qp];
61  r += _sigma / h_elem * (_u[_qp] - _u_neighbor[_qp]) * _test[_i][_qp];
62  break;
63 
64  case Moose::Neighbor:
65  r += 0.5 * (_diff[_qp] * _grad_u[_qp] * _normals[_qp] +
68  r += _epsilon * 0.5 * (_u[_qp] - _u_neighbor[_qp]) * _diff_neighbor[_qp] *
70  r -= _sigma / h_elem * (_u[_qp] - _u_neighbor[_qp]) * _test_neighbor[_i][_qp];
71  break;
72  }
73 
74  return r;
75 }
76 
77 Real
79 {
80  Real r = 0;
81 
82  const unsigned int elem_b_order = _var.order();
83  const double h_elem =
84  _current_elem->volume() / _current_side_elem->volume() * 1. / Utility::pow<2>(elem_b_order);
85 
86  switch (type)
87  {
89  r -= 0.5 * _diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test[_i][_qp];
90  r += _epsilon * 0.5 * _phi[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] * _normals[_qp];
91  r += _sigma / h_elem * _phi[_j][_qp] * _test[_i][_qp];
92  break;
93 
95  r -= 0.5 * _diff_neighbor[_qp] * _grad_phi_neighbor[_j][_qp] * _normals[_qp] * _test[_i][_qp];
96  r += _epsilon * 0.5 * -_phi_neighbor[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] *
97  _normals[_qp];
98  r += _sigma / h_elem * -_phi_neighbor[_j][_qp] * _test[_i][_qp];
99  break;
100 
102  r += 0.5 * _diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test_neighbor[_i][_qp];
103  r += _epsilon * 0.5 * _phi[_j][_qp] * _diff_neighbor[_qp] * _grad_test_neighbor[_i][_qp] *
104  _normals[_qp];
105  r -= _sigma / h_elem * _phi[_j][_qp] * _test_neighbor[_i][_qp];
106  break;
107 
109  r += 0.5 * _diff_neighbor[_qp] * _grad_phi_neighbor[_j][_qp] * _normals[_qp] *
110  _test_neighbor[_i][_qp];
111  r += _epsilon * 0.5 * -_phi_neighbor[_j][_qp] * _diff_neighbor[_qp] *
112  _grad_test_neighbor[_i][_qp] * _normals[_qp];
113  r -= _sigma / h_elem * -_phi_neighbor[_j][_qp] * _test_neighbor[_i][_qp];
114  break;
115  }
116 
117  return r;
118 }
const VariablePhiGradient & _grad_phi_neighbor
Gradient of side shape function.
Definition: DGKernel.h:171
const VariableTestGradient & _grad_test
Gradient of side shape function.
Definition: DGKernel.h:164
Real _epsilon
Definition: DGDiffusion.h:45
InputParameters validParams< DGDiffusion >()
Definition: DGDiffusion.C:24
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...
const VariableTestGradient & _grad_test_neighbor
Gradient of side shape function.
Definition: DGKernel.h:176
DGResidualType
Definition: MooseTypes.h:184
const Elem *& _current_side_elem
Current side element.
Definition: DGKernel.h:134
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...
const VariablePhiGradient & _grad_phi
Definition: DGKernel.h:158
MooseVariable & _var
Definition: DGKernel.h:119
const VariableTestValue & _test_neighbor
Side test function.
Definition: DGKernel.h:174
unsigned int _j
Definition: DGKernel.h:147
const VariableGradient & _grad_u_neighbor
Holds the current solution gradient at the current quadrature point.
Definition: DGKernel.h:181
The DGKernel class is responsible for calculating the residuals for various physics on internal sides...
Definition: DGKernel.h:47
Order order() const
Get the order of this variable Note: Order enum can be implicitly converted to unsigned int...
const MaterialProperty< Real > & _diff_neighbor
Definition: DGDiffusion.h:48
virtual Real computeQpResidual(Moose::DGResidualType type) override
This is the virtual that derived classes should override for computing the residual on neighboring el...
Definition: DGDiffusion.C:45
const VariableTestValue & _test
Side shape function.
Definition: DGKernel.h:162
DGJacobianType
Definition: MooseTypes.h:190
MatType type
unsigned int _i
Definition: DGKernel.h:147
const VariableGradient & _grad_u
Holds the current solution gradient at the current quadrature point on the face.
Definition: DGKernel.h:155
DGDiffusion(const InputParameters &parameters)
Definition: DGDiffusion.C:35
const MooseArray< Point > & _normals
Normal vectors at the quadrature points.
Definition: DGKernel.h:166
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...
const Elem *& _current_elem
Definition: DGKernel.h:123
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
virtual Real computeQpJacobian(Moose::DGJacobianType type) override
This is the virtual that derived classes should override for computing the Jacobian on neighboring el...
Definition: DGDiffusion.C:78
const MaterialProperty< Real > & _diff
Definition: DGDiffusion.h:47
unsigned int _qp
Definition: DGKernel.h:141