www.mooseframework.org
DGFunctionDiffusionDirichletBC.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 
11 
12 // MOOSE includes
13 #include "Function.h"
14 #include "MooseVariableFE.h"
15 
16 #include "libmesh/numeric_vector.h"
17 #include "libmesh/utility.h"
18 
19 // C++ includes
20 #include <cmath>
21 
23 
26 {
28  params.addClassDescription(
29  "Diffusion Dirichlet boundary condition for discontinuous Galerkin method.");
30  params.addParam<Real>("value", 0.0, "The value the variable should have on the boundary");
31  params.addRequiredParam<FunctionName>("function", "The forcing function.");
32  params.addRequiredParam<Real>("epsilon", "Epsilon");
33  params.addRequiredParam<Real>("sigma", "Sigma");
34  params.addParam<MaterialPropertyName>(
35  "diff", 1, "The diffusion (or thermal conductivity or viscosity) coefficient.");
36 
37  return params;
38 }
39 
41  : IntegratedBC(parameters),
42  _func(getFunction("function")),
43  _epsilon(getParam<Real>("epsilon")),
44  _sigma(getParam<Real>("sigma")),
45  _diff(getMaterialProperty<Real>("diff"))
46 {
47 }
48 
49 Real
51 {
52  const int elem_b_order = std::max(libMesh::Order(1), _var.order());
53  const Real h_elem =
55 
56  const Real fn = _func.value(_t, _q_point[_qp]);
57  Real r = 0.0;
58  r -= (_diff[_qp] * _grad_u[_qp] * _normals[_qp] * _test[_i][_qp]);
59  r += _epsilon * (_u[_qp] - fn) * _diff[_qp] * _grad_test[_i][_qp] * _normals[_qp];
60  r += _sigma / h_elem * (_u[_qp] - fn) * _test[_i][_qp];
61 
62  return r;
63 }
64 
65 Real
67 {
68  const int elem_b_order = std::max(libMesh::Order(1), _var.order());
69  const Real h_elem =
71 
72  Real r = 0.0;
73  r -= (_diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test[_i][_qp]);
74  r += _epsilon * _phi[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] * _normals[_qp];
75  r += _sigma / h_elem * _phi[_j][_qp] * _test[_i][_qp];
76 
77  return r;
78 }
const VariableTestValue & _test
test function values (in QPs)
Definition: IntegratedBC.h:97
virtual Real computeQpResidual() override
Method for computing the residual at quadrature points.
const MooseArray< Point > & _normals
normals at quadrature points
Definition: IntegratedBC.h:85
DGFunctionDiffusionDirichletBC(const InputParameters &parameters)
virtual Real computeQpJacobian() override
Method for computing the diagonal Jacobian at quadrature points.
const VariableGradient & _grad_u
the gradient of the unknown variable this BC is acting on
Definition: IntegratedBC.h:106
static InputParameters validParams()
Definition: IntegratedBC.C:22
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int _i
i-th, j-th index for enumerating test and shape functions
registerMooseObject("MooseApp", DGFunctionDiffusionDirichletBC)
const VariablePhiValue & _phi
shape function values (in QPs)
Definition: IntegratedBC.h:90
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...
unsigned int _qp
quadrature point index
auto max(const L &left, const R &right)
const VariablePhiGradient & _grad_phi
gradients of shape functions (in QPs)
Definition: IntegratedBC.h:92
static InputParameters validParams()
Factory constructor, takes parameters so that all derived classes can be built using the same constru...
const MooseArray< Point > & _q_point
active quadrature points
const Real & _current_elem_volume
Volume of the current element.
const MaterialProperty< Real > & _diff
Base class for deriving any boundary condition of a integrated type.
Definition: IntegratedBC.h:18
const Real & _current_side_volume
Volume of the current side.
Order order() const
Get the order of this variable Note: Order enum can be implicitly converted to unsigned int...
MooseVariable & _var
Definition: IntegratedBC.h:82
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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 VariableTestGradient & _grad_test
gradients of test functions (in QPs)
Definition: IntegratedBC.h:99
virtual Real value(Real t, const Point &p) const
Override this to evaluate the scalar function at point (t,x,y,z), by default this returns zero...
Definition: Function.C:41
const VariableValue & _u
the values of the unknown variable this BC is acting on
Definition: IntegratedBC.h:104
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