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 
24 template <>
27 {
29  params.addParam<Real>("value", 0.0, "The value the variable should have on the boundary");
30  params.addRequiredParam<FunctionName>("function", "The forcing function.");
31  params.addRequiredParam<Real>("epsilon", "Epsilon");
32  params.addRequiredParam<Real>("sigma", "Sigma");
33  params.addParam<MaterialPropertyName>(
34  "diff", 1, "The diffusion (or thermal conductivity or viscosity) coefficient.");
35 
36  return params;
37 }
38 
40  : IntegratedBC(parameters),
41  _func(getFunction("function")),
42  _epsilon(getParam<Real>("epsilon")),
43  _sigma(getParam<Real>("sigma")),
44  _diff(getMaterialProperty<Real>("diff"))
45 {
46 }
47 
48 Real
50 {
51  const unsigned int elem_b_order = _var.order();
52  const double h_elem =
53  _current_elem->volume() / _current_side_elem->volume() * 1. / Utility::pow<2>(elem_b_order);
54 
55  Real fn = _func.value(_t, _q_point[_qp]);
56  Real r = 0;
57  r -= (_diff[_qp] * _grad_u[_qp] * _normals[_qp] * _test[_i][_qp]);
58  r += _epsilon * (_u[_qp] - fn) * _diff[_qp] * _grad_test[_i][_qp] * _normals[_qp];
59  r += _sigma / h_elem * (_u[_qp] - fn) * _test[_i][_qp];
60 
61  return r;
62 }
63 
64 Real
66 {
67  const unsigned int elem_b_order = _var.order();
68  const double h_elem =
69  _current_elem->volume() / _current_side_elem->volume() * 1. / Utility::pow<2>(elem_b_order);
70 
71  Real r = 0;
72  r -= (_diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test[_i][_qp]);
73  r += _epsilon * _phi[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] * _normals[_qp];
74  r += _sigma / h_elem * _phi[_j][_qp] * _test[_i][_qp];
75 
76  return r;
77 }
const VariableTestValue & _test
test function values (in QPs)
Definition: IntegratedBC.h:68
InputParameters validParams< DGFunctionDiffusionDirichletBC >()
virtual Real computeQpResidual() override
method for computing the residual at quadrature points
virtual Real value(Real t, const Point &p)
Override this to evaluate the scalar function at point (t,x,y,z), by default this returns zero...
Definition: Function.C:38
InputParameters validParams< IntegratedBC >()
Definition: IntegratedBC.C:23
const MooseArray< Point > & _normals
normals at quadrature points
Definition: IntegratedBC.h:56
DGFunctionDiffusionDirichletBC(const InputParameters &parameters)
Factory constructor, takes parameters so that all derived classes can be built using the same constru...
const Elem *const & _current_elem
current element
const Elem *const & _current_side_elem
current side element
const VariableGradient & _grad_u
the gradient of the unknown variable this BC is acting on
Definition: IntegratedBC.h:77
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:61
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
const VariablePhiGradient & _grad_phi
gradients of shape functions (in QPs)
Definition: IntegratedBC.h:63
const MooseArray< Point > & _q_point
active quadrature points
const MaterialProperty< Real > & _diff
Base class for deriving any boundary condition of a integrated type.
Definition: IntegratedBC.h:24
Order order() const
Get the order of this variable Note: Order enum can be implicitly converted to unsigned int...
MooseVariable & _var
Definition: IntegratedBC.h:53
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:70
const VariableValue & _u
the values of the unknown variable this BC is acting on
Definition: IntegratedBC.h:75