www.mooseframework.org
IntegratedBC.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 
15 #include "IntegratedBC.h"
16 
17 // MOOSE includes
18 #include "Assembly.h"
19 #include "SubProblem.h"
20 #include "SystemBase.h"
21 #include "MooseVariable.h"
22 #include "MooseVariableScalar.h"
23 
24 #include "libmesh/quadrature.h"
25 
26 template <>
29 {
31  params += validParams<RandomInterface>();
33 
34  params.addParam<std::vector<AuxVariableName>>(
35  "save_in",
36  "The name of auxiliary variables to save this BC's residual contributions to. "
37  "Everything about that variable must match everything about this variable (the "
38  "type, what blocks it's on, etc.)");
39  params.addParam<std::vector<AuxVariableName>>(
40  "diag_save_in",
41  "The name of auxiliary variables to save this BC's diagonal jacobian "
42  "contributions to. Everything about that variable must match everything "
43  "about this variable (the type, what blocks it's on, etc.)");
44 
45  params.addParamNamesToGroup("diag_save_in save_in", "Advanced");
46 
47  // Integrated BCs always rely on Boundary MaterialData
48  params.set<Moose::MaterialDataType>("_material_data_type") = Moose::BOUNDARY_MATERIAL_DATA;
49 
50  return params;
51 }
52 
54  : BoundaryCondition(parameters, false), // False is because this is NOT nodal
55  RandomInterface(parameters, _fe_problem, _tid, false),
57  MaterialPropertyInterface(this, boundaryIDs()),
58  _current_elem(_assembly.elem()),
59  _current_elem_volume(_assembly.elemVolume()),
60  _current_side(_assembly.side()),
61  _current_side_elem(_assembly.sideElem()),
62  _current_side_volume(_assembly.sideElemVolume()),
63 
64  _normals(_var.normals()),
65 
66  _qrule(_assembly.qRuleFace()),
67  _q_point(_assembly.qPointsFace()),
68  _JxW(_assembly.JxWFace()),
69  _coord(_assembly.coordTransformation()),
70 
71  _phi(_assembly.phiFace()),
72  _grad_phi(_assembly.gradPhiFace()),
73 
74  _test(_var.phiFace()),
75  _grad_test(_var.gradPhiFace()),
76 
77  _u(_is_implicit ? _var.sln() : _var.slnOld()),
78  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld()),
79 
80  _save_in_strings(parameters.get<std::vector<AuxVariableName>>("save_in")),
81  _diag_save_in_strings(parameters.get<std::vector<AuxVariableName>>("diag_save_in"))
82 {
83  _save_in.resize(_save_in_strings.size());
84  _diag_save_in.resize(_diag_save_in_strings.size());
85 
86  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
87  {
89 
90  if (var->feType() != _var.feType())
91  mooseError("Error in " + name() + ". When saving residual values in an Auxiliary variable "
92  "the AuxVariable must be the same type as the nonlinear "
93  "variable the object is acting on.");
94 
95  _save_in[i] = var;
98  }
99 
100  _has_save_in = _save_in.size() > 0;
101 
102  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
103  {
105 
106  if (var->feType() != _var.feType())
107  mooseError("Error in " + name() + ". When saving diagonal Jacobian values in an Auxiliary "
108  "variable the AuxVariable must be the same type as the "
109  "nonlinear variable the object is acting on.");
110 
111  _diag_save_in[i] = var;
114  }
115 
116  _has_diag_save_in = _diag_save_in.size() > 0;
117 }
118 
120 
121 void
123 {
124  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
125  _local_re.resize(re.size());
126  _local_re.zero();
127 
128  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
129  for (_i = 0; _i < _test.size(); _i++)
131 
132  re += _local_re;
133 
134  if (_has_save_in)
135  {
136  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
137  for (unsigned int i = 0; i < _save_in.size(); i++)
138  _save_in[i]->sys().solution().add_vector(_local_re, _save_in[i]->dofIndices());
139  }
140 }
141 
142 void
144 {
145  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
146  _local_ke.resize(ke.m(), ke.n());
147  _local_ke.zero();
148 
149  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
150  for (_i = 0; _i < _test.size(); _i++)
151  for (_j = 0; _j < _phi.size(); _j++)
153 
154  ke += _local_ke;
155 
156  if (_has_diag_save_in)
157  {
158  unsigned int rows = ke.m();
159  DenseVector<Number> diag(rows);
160  for (unsigned int i = 0; i < rows; i++)
161  diag(i) = _local_ke(i, i);
162 
163  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
164  for (unsigned int i = 0; i < _diag_save_in.size(); i++)
165  _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
166  }
167 }
168 
169 void
171 {
172  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);
173 
174  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
175  for (_i = 0; _i < _test.size(); _i++)
176  for (_j = 0; _j < _phi.size(); _j++)
177  {
178  if (_var.number() == jvar)
179  ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpJacobian();
180  else
181  ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar);
182  }
183 }
184 
185 void
187 {
188  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);
189 
191  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
192  for (_i = 0; _i < _test.size(); _i++)
193  for (_j = 0; _j < jv.order(); _j++)
194  ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar);
195 }
196 
197 Real
199 {
200  return 0;
201 }
202 
203 Real
205 {
206  return 0;
207 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
InputParameters validParams< MaterialPropertyInterface >()
const VariableTestValue & _test
test function values (in QPs)
Definition: IntegratedBC.h:105
Interface for objects that need parallel consistent random numbers without patterns over the course o...
SubProblem & _subproblem
Reference to SubProblem.
std::vector< MooseVariable * > _save_in
Definition: IntegratedBC.h:124
const FEType & feType() const
Get the type of finite element object.
Class for stuff related to variables.
Definition: MooseVariable.h:43
InputParameters validParams< IntegratedBC >()
Definition: IntegratedBC.C:28
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name)=0
Returns the variable reference for requested variable which may be in any system. ...
const MooseArray< Real > & _coord
coordinate transformation
Definition: IntegratedBC.h:91
IntegratedBC(const InputParameters &parameters)
Definition: IntegratedBC.C:53
MaterialDataType
MaterialData types.
Definition: MooseTypes.h:129
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
DenseVector< Number > _local_re
Holds residual entries as their accumulated by this Kernel.
Definition: IntegratedBC.h:117
void addMooseVariableDependency(MooseVariable *var)
Call this function to add the passed in MooseVariable as a variable that this object depends on...
const MooseArray< Real > & _JxW
transformed Jacobian weights
Definition: IntegratedBC.h:89
virtual Real computeQpJacobian()
Definition: IntegratedBC.C:198
const VariablePhiValue & _phi
shape function values (in QPs)
Definition: IntegratedBC.h:98
DenseVector< Number > & residualBlock(unsigned int var_num, Moose::KernelType type=Moose::KT_NONTIME)
Definition: Assembly.h:504
virtual Real computeQpResidual()=0
virtual void computeJacobianBlock(unsigned int jvar)
Computes d-ivar-residual / d-jvar...
Definition: IntegratedBC.C:170
MooseVariable & _var
variable this BC works on
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: IntegratedBC.h:128
std::vector< AuxVariableName > _diag_save_in_strings
Definition: IntegratedBC.h:130
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: IntegratedBC.h:123
InputParameters validParams< RandomInterface >()
std::vector< AuxVariableName > _save_in_strings
Definition: IntegratedBC.h:125
virtual void computeJacobian()
Definition: IntegratedBC.C:143
Assembly & _assembly
Reference to assembly.
Order order() const
Get the order of this variable Note: Order enum can be implicitly converted to unsigned int...
virtual Real computeQpOffDiagJacobian(unsigned int jvar)
This is the virtual that derived classes should override for computing an off-diagonal jacobian compo...
Definition: IntegratedBC.C:204
virtual ~IntegratedBC()
Definition: IntegratedBC.C:119
Base class for creating new types of boundary conditions.
unsigned int number() const
Get variable number coming from libMesh.
virtual void computeResidual()
Definition: IntegratedBC.C:122
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
An interface for accessing Materials.
DenseMatrix< Number > _local_ke
Holds residual entries as they are accumulated by this Kernel.
Definition: IntegratedBC.h:120
virtual void addVariableToZeroOnResidual(std::string var_name)
Adds this variable to the list of variables to be zeroed during each residual evaluation.
Definition: SystemBase.C:150
unsigned int _qp
quadrature point index
Definition: IntegratedBC.h:83
unsigned int _j
Definition: IntegratedBC.h:93
SystemBase & _sys
Reference to SystemBase.
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
Class for scalar variables (they are different).
std::vector< MooseVariable * > _diag_save_in
Definition: IntegratedBC.h:129
virtual void addVariableToZeroOnJacobian(std::string var_name)
Adds this variable to the list of variables to be zeroed during each Jacobian evaluation.
Definition: SystemBase.C:156
void computeJacobianBlockScalar(unsigned int jvar)
Computes jacobian block with respect to a scalar variable.
Definition: IntegratedBC.C:186
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...
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:887
unsigned int _i
i-th, j-th index for enumerating test and shape functions
Definition: IntegratedBC.h:93
SystemBase & sys()
Get the system this variable is part of.
THREAD_ID _tid
Thread id.
virtual MooseVariableScalar & getScalarVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a scalar variable with specified number.
Definition: SystemBase.C:121
InputParameters validParams< BoundaryCondition >()
QBase *& _qrule
active quadrature rule
Definition: IntegratedBC.h:85