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 #include "MooseVariableFE.h"
12 
13 #include "libmesh/utility.h"
14 
15 registerMooseObject("MooseApp", DGDiffusion);
16 
19 {
21  params.addClassDescription("Computes residual contribution for the diffusion operator using "
22  "discontinous Galerkin method.");
23  // See header file for sigma and epsilon
24  params.addRequiredParam<Real>("sigma", "sigma");
25  params.addRequiredParam<Real>("epsilon", "epsilon");
26  params.addParam<MaterialPropertyName>(
27  "diff", 1., "The diffusion (or thermal conductivity or viscosity) coefficient.");
28  return params;
29 }
30 
32  : DGKernel(parameters),
33  _epsilon(getParam<Real>("epsilon")),
34  _sigma(getParam<Real>("sigma")),
35  _diff(getMaterialProperty<Real>("diff")),
36  _diff_neighbor(getNeighborMaterialProperty<Real>("diff"))
37 {
38 }
39 
40 Real
42 {
43  Real r = 0.0;
44 
45  const int elem_b_order = std::max(libMesh::Order(1), _var.order());
46  const Real h_elem =
48 
49  switch (type)
50  {
51  case Moose::Element:
52  r -= 0.5 *
53  (_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 *
63  (_diff[_qp] * _grad_u[_qp] * _normals[_qp] +
66  r += _epsilon * 0.5 * (_u[_qp] - _u_neighbor[_qp]) * _diff_neighbor[_qp] *
68  r -= _sigma / h_elem * (_u[_qp] - _u_neighbor[_qp]) * _test_neighbor[_i][_qp];
69  break;
70  }
71 
72  return r;
73 }
74 
75 Real
77 {
78  Real r = 0.0;
79 
80  const int elem_b_order = std::max(libMesh::Order(1), _var.order());
81  const Real h_elem =
83 
84  switch (type)
85  {
87  r -= 0.5 * _diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test[_i][_qp];
88  r += _epsilon * 0.5 * _phi[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] * _normals[_qp];
89  r += _sigma / h_elem * _phi[_j][_qp] * _test[_i][_qp];
90  break;
91 
94  r += _epsilon * 0.5 * -_phi_neighbor[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] *
95  _normals[_qp];
96  r += _sigma / h_elem * -_phi_neighbor[_j][_qp] * _test[_i][_qp];
97  break;
98 
100  r += 0.5 * _diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test_neighbor[_i][_qp];
101  r += _epsilon * 0.5 * _phi[_j][_qp] * _diff_neighbor[_qp] * _grad_test_neighbor[_i][_qp] *
102  _normals[_qp];
103  r -= _sigma / h_elem * _phi[_j][_qp] * _test_neighbor[_i][_qp];
104  break;
105 
107  r += 0.5 * _diff_neighbor[_qp] * _grad_phi_neighbor[_j][_qp] * _normals[_qp] *
109  r += _epsilon * 0.5 * -_phi_neighbor[_j][_qp] * _diff_neighbor[_qp] *
111  r -= _sigma / h_elem * -_phi_neighbor[_j][_qp] * _test_neighbor[_i][_qp];
112  break;
113  }
114 
115  return r;
116 }
DG kernel for diffusion.
Definition: DGDiffusion.h:24
const VariablePhiGradient & _grad_phi_neighbor
Gradient of side shape function.
Definition: DGKernel.h:108
const VariableTestGradient & _grad_test
Gradient of side shape function.
Definition: DGKernel.h:104
Real _epsilon
Definition: DGDiffusion.h:35
const Real & _current_side_volume
The volume (or length) of the current side.
Definition: DGKernelBase.h:100
const VariableValue & _u_neighbor
Holds the current solution at the current quadrature point.
Definition: DGKernel.h:114
const VariablePhiValue & _phi_neighbor
Side shape function.
Definition: DGKernel.h:106
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:112
unsigned int _i
Definition: DGKernelBase.h:136
DGResidualType
Definition: MooseTypes.h:656
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
Gradient of shape function.
Definition: DGKernel.h:100
auto max(const L &left, const R &right)
MooseVariable & _var
Variable this kernel operates on.
Definition: DGKernel.h:92
static InputParameters validParams()
Factory constructor initializes all internal references needed for residual computation.
Definition: DGKernel.C:27
const VariableTestValue & _test_neighbor
Side test function.
Definition: DGKernel.h:110
unsigned int _qp
Definition: DGKernelBase.h:134
const VariableGradient & _grad_u_neighbor
Holds the current solution gradient at the current quadrature point.
Definition: DGKernel.h:116
unsigned int _j
Definition: DGKernelBase.h:136
The DGKernel class is responsible for calculating the residuals for various physics on internal sides...
Definition: DGKernel.h:18
const std::string & type() const
Get the type of this class.
Definition: MooseBase.h:50
const MooseArray< Point > & _normals
Normal vectors at the quadrature points.
Definition: DGKernelBase.h:113
const MaterialProperty< Real > & _diff_neighbor
Definition: DGDiffusion.h:38
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:41
const VariableTestValue & _test
test functions
Definition: DGKernel.h:102
Order order() const
Get the order of this variable Note: Order enum can be implicitly converted to unsigned int...
DGJacobianType
Definition: MooseTypes.h:662
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const VariableGradient & _grad_u
Holds the current solution gradient at the current quadrature point on the face.
Definition: DGKernel.h:96
DGDiffusion(const InputParameters &parameters)
Definition: DGDiffusion.C:31
registerMooseObject("MooseApp", DGDiffusion)
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...
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 VariablePhiValue & _phi
Shape functions.
Definition: DGKernel.h:98
const VariableValue & _u
Holds the current solution at the current quadrature point on the face.
Definition: DGKernel.h:94
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:76
const MaterialProperty< Real > & _diff
Definition: DGDiffusion.h:37
static InputParameters validParams()
Definition: DGDiffusion.C:18
const Real & _current_elem_volume
The volume (or length) of the current element.
Definition: DGKernelBase.h:87
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template pow< 2 >(tan(_arg))+1.0) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(sqrt