www.mooseframework.org
DGDiffusion.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 "DGDiffusion.h"
11 
12 // MOOSE includes
13 #include "MooseVariableFE.h"
14 
15 #include "libmesh/utility.h"
16 
17 registerMooseObject("MooseApp", DGDiffusion);
18 
19 template <>
22 {
24  // See header file for sigma and epsilon
25  params.addRequiredParam<Real>("sigma", "sigma");
26  params.addRequiredParam<Real>("epsilon", "epsilon");
27  params.addParam<MaterialPropertyName>(
28  "diff", 1., "The diffusion (or thermal conductivity or viscosity) coefficient.");
29  return params;
30 }
31 
33  : DGKernel(parameters),
34  _epsilon(getParam<Real>("epsilon")),
35  _sigma(getParam<Real>("sigma")),
36  _diff(getMaterialProperty<Real>("diff")),
37  _diff_neighbor(getNeighborMaterialProperty<Real>("diff"))
38 {
39 }
40 
41 Real
43 {
44  Real r = 0;
45 
46  const unsigned int elem_b_order = _var.order();
47  const double h_elem =
48  _current_elem->volume() / _current_side_elem->volume() * 1. / Utility::pow<2>(elem_b_order);
49 
50  switch (type)
51  {
52  case Moose::Element:
53  r -= 0.5 * (_diff[_qp] * _grad_u[_qp] * _normals[_qp] +
55  _test[_i][_qp];
56  r += _epsilon * 0.5 * (_u[_qp] - _u_neighbor[_qp]) * _diff[_qp] * _grad_test[_i][_qp] *
57  _normals[_qp];
58  r += _sigma / h_elem * (_u[_qp] - _u_neighbor[_qp]) * _test[_i][_qp];
59  break;
60 
61  case Moose::Neighbor:
62  r += 0.5 * (_diff[_qp] * _grad_u[_qp] * _normals[_qp] +
65  r += _epsilon * 0.5 * (_u[_qp] - _u_neighbor[_qp]) * _diff_neighbor[_qp] *
67  r -= _sigma / h_elem * (_u[_qp] - _u_neighbor[_qp]) * _test_neighbor[_i][_qp];
68  break;
69  }
70 
71  return r;
72 }
73 
74 Real
76 {
77  Real r = 0;
78 
79  const unsigned int elem_b_order = _var.order();
80  const double h_elem =
81  _current_elem->volume() / _current_side_elem->volume() * 1. / Utility::pow<2>(elem_b_order);
82 
83  switch (type)
84  {
86  r -= 0.5 * _diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test[_i][_qp];
87  r += _epsilon * 0.5 * _phi[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] * _normals[_qp];
88  r += _sigma / h_elem * _phi[_j][_qp] * _test[_i][_qp];
89  break;
90 
93  r += _epsilon * 0.5 * -_phi_neighbor[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] *
94  _normals[_qp];
95  r += _sigma / h_elem * -_phi_neighbor[_j][_qp] * _test[_i][_qp];
96  break;
97 
99  r += 0.5 * _diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test_neighbor[_i][_qp];
100  r += _epsilon * 0.5 * _phi[_j][_qp] * _diff_neighbor[_qp] * _grad_test_neighbor[_i][_qp] *
101  _normals[_qp];
102  r -= _sigma / h_elem * _phi[_j][_qp] * _test_neighbor[_i][_qp];
103  break;
104 
106  r += 0.5 * _diff_neighbor[_qp] * _grad_phi_neighbor[_j][_qp] * _normals[_qp] *
108  r += _epsilon * 0.5 * -_phi_neighbor[_j][_qp] * _diff_neighbor[_qp] *
110  r -= _sigma / h_elem * -_phi_neighbor[_j][_qp] * _test_neighbor[_i][_qp];
111  break;
112  }
113 
114  return r;
115 }
DG kernel for diffusion.
Definition: DGDiffusion.h:30
const VariableTestValue & _test_neighbor
Side test function.
Definition: DGKernelBase.h:171
const VariablePhiGradient & _grad_phi_neighbor
Gradient of side shape function.
Definition: DGKernelBase.h:168
const VariablePhiValue & _phi_neighbor
Side shape function.
Definition: DGKernelBase.h:166
Real _epsilon
Definition: DGDiffusion.h:39
InputParameters validParams< DGDiffusion >()
Definition: DGDiffusion.C:21
const VariableValue & _u_neighbor
Holds the current solution at the current quadrature point.
Definition: DGKernelBase.h:176
InputParameters validParams< DGKernel >()
Definition: DGKernel.C:28
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
const VariableValue & _u
Holds the current solution at the current quadrature point on the face.
Definition: DGKernelBase.h:149
unsigned int _i
Definition: DGKernelBase.h:144
DGResidualType
Definition: MooseTypes.h:509
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 std::string & type() const
Get the type of this object.
Definition: MooseObject.h:53
unsigned int _qp
Definition: DGKernelBase.h:138
const VariableGradient & _grad_u_neighbor
Holds the current solution gradient at the current quadrature point.
Definition: DGKernelBase.h:178
unsigned int _j
Definition: DGKernelBase.h:144
The DGKernel class is responsible for calculating the residuals for various physics on internal sides...
Definition: DGKernel.h:23
const VariablePhiGradient & _grad_phi
Definition: DGKernelBase.h:155
const MooseArray< Point > & _normals
Normal vectors at the quadrature points.
Definition: DGKernelBase.h:163
const MaterialProperty< Real > & _diff_neighbor
Definition: DGDiffusion.h:42
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:42
Order order() const
Get the order of this variable Note: Order enum can be implicitly converted to unsigned int...
DGJacobianType
Definition: MooseTypes.h:515
MatType type
const VariablePhiValue & _phi
Definition: DGKernelBase.h:154
const VariableTestGradient & _grad_test
Gradient of side shape function.
Definition: DGKernelBase.h:161
const Elem *const & _current_elem
Definition: DGKernelBase.h:120
const VariableTestValue & _test
Side shape function.
Definition: DGKernelBase.h:159
DGDiffusion(const InputParameters &parameters)
Definition: DGDiffusion.C:32
const VariableGradient & _grad_u
Holds the current solution gradient at the current quadrature point on the face.
Definition: DGKernelBase.h:152
registerMooseObject("MooseApp", DGDiffusion)
const Elem *& _current_side_elem
Current side element.
Definition: DGKernelBase.h:131
MooseVariable & _var
Definition: DGKernelBase.h:117
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 VariableTestGradient & _grad_test_neighbor
Gradient of side shape function.
Definition: DGKernelBase.h:173
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:75
const MaterialProperty< Real > & _diff
Definition: DGDiffusion.h:41