www.mooseframework.org
DGFunctionDiffusionDirichletBC.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 
16 
17 // MOOSE includes
18 #include "Function.h"
19 #include "MooseVariable.h"
20 
21 #include "libmesh/numeric_vector.h"
22 #include "libmesh/utility.h"
23 
24 // C++ includes
25 #include <cmath>
26 
27 template <>
30 {
32  params.addParam<Real>("value", 0.0, "The value the variable should have on the boundary");
33  params.addRequiredParam<FunctionName>("function", "The forcing function.");
34  params.addRequiredParam<Real>("epsilon", "Epsilon");
35  params.addRequiredParam<Real>("sigma", "Sigma");
36  params.addParam<MaterialPropertyName>(
37  "diff", 1, "The diffusion (or thermal conductivity or viscosity) coefficient.");
38 
39  return params;
40 }
41 
43  : IntegratedBC(parameters),
44  _func(getFunction("function")),
45  _epsilon(getParam<Real>("epsilon")),
46  _sigma(getParam<Real>("sigma")),
47  _diff(getMaterialProperty<Real>("diff"))
48 {
49 }
50 
51 Real
53 {
54  const unsigned int elem_b_order = _var.order();
55  const double h_elem =
56  _current_elem->volume() / _current_side_elem->volume() * 1. / Utility::pow<2>(elem_b_order);
57 
58  Real fn = _func.value(_t, _q_point[_qp]);
59  Real r = 0;
60  r -= (_diff[_qp] * _grad_u[_qp] * _normals[_qp] * _test[_i][_qp]);
61  r += _epsilon * (_u[_qp] - fn) * _diff[_qp] * _grad_test[_i][_qp] * _normals[_qp];
62  r += _sigma / h_elem * (_u[_qp] - fn) * _test[_i][_qp];
63 
64  return r;
65 }
66 
67 Real
69 {
70  const unsigned int elem_b_order = _var.order();
71  const double h_elem =
72  _current_elem->volume() / _current_side_elem->volume() * 1. / Utility::pow<2>(elem_b_order);
73 
74  Real r = 0;
75  r -= (_diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test[_i][_qp]);
76  r += _epsilon * _phi[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] * _normals[_qp];
77  r += _sigma / h_elem * _phi[_j][_qp] * _test[_i][_qp];
78 
79  return r;
80 }
const VariableTestValue & _test
test function values (in QPs)
Definition: IntegratedBC.h:105
InputParameters validParams< DGFunctionDiffusionDirichletBC >()
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:43
InputParameters validParams< IntegratedBC >()
Definition: IntegratedBC.C:28
const MooseArray< Point > & _normals
normals at quadrature points
Definition: IntegratedBC.h:80
DGFunctionDiffusionDirichletBC(const InputParameters &parameters)
Factory constructor, takes parameters so that all derived classes can be built using the same constru...
const Elem *& _current_side_elem
current side element
Definition: IntegratedBC.h:75
const VariableGradient & _grad_u
the gradient of the unknown variable this BC is acting on
Definition: IntegratedBC.h:114
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
const VariablePhiValue & _phi
shape function values (in QPs)
Definition: IntegratedBC.h:98
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
gradients of shape functions (in QPs)
Definition: IntegratedBC.h:100
MooseVariable & _var
variable this BC works on
const MaterialProperty< Real > & _diff
Base class for deriving any boundary condition of a integrated type.
Definition: IntegratedBC.h:33
Order order() const
Get the order of this variable Note: Order enum can be implicitly converted to unsigned int...
unsigned int _qp
quadrature point index
Definition: IntegratedBC.h:83
unsigned int _j
Definition: IntegratedBC.h:93
const MooseArray< Point > & _q_point
active quadrature points
Definition: IntegratedBC.h:87
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:107
const Elem *& _current_elem
current element
Definition: IntegratedBC.h:69
unsigned int _i
i-th, j-th index for enumerating test and shape functions
Definition: IntegratedBC.h:93
const VariableValue & _u
the values of the unknown variable this BC is acting on
Definition: IntegratedBC.h:112