www.mooseframework.org
IntegratedBC.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 "IntegratedBC.h"
11 
12 // MOOSE includes
13 #include "Assembly.h"
14 #include "SubProblem.h"
15 #include "SystemBase.h"
16 #include "MooseVariableFE.h"
17 #include "MooseVariableScalar.h"
18 
19 #include "libmesh/quadrature.h"
20 
21 template <>
24 {
26  return params;
27 }
28 
30  : IntegratedBCBase(parameters),
31  MooseVariableInterface<Real>(this,
32  false,
33  "variable",
36  _var(*mooseVariable()),
37  _normals(_assembly.normals()),
38  _phi(_assembly.phiFace(_var)),
39  _grad_phi(_assembly.gradPhiFace(_var)),
40  _test(_var.phiFace()),
41  _grad_test(_var.gradPhiFace()),
42  _u(_is_implicit ? _var.sln() : _var.slnOld()),
43  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld())
44 {
46 
47  _save_in.resize(_save_in_strings.size());
48  _diag_save_in.resize(_diag_save_in_strings.size());
49 
50  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
51  {
53 
54  if (var->feType() != _var.feType())
55  paramError(
56  "save_in",
57  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
59  _save_in[i] = var;
62  }
63 
64  _has_save_in = _save_in.size() > 0;
65 
66  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
67  {
69 
70  if (var->feType() != _var.feType())
71  paramError(
72  "diag_save_in",
73  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
75 
76  _diag_save_in[i] = var;
79  }
80 
81  _has_diag_save_in = _diag_save_in.size() > 0;
82 }
83 
84 void
86 {
88 
89  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
90  for (_i = 0; _i < _test.size(); _i++)
92 
94 
95  if (_has_save_in)
96  {
97  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
98  for (unsigned int i = 0; i < _save_in.size(); i++)
99  _save_in[i]->sys().solution().add_vector(_local_re, _save_in[i]->dofIndices());
100  }
101 }
102 
103 void
105 {
107 
108  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
109  for (_i = 0; _i < _test.size(); _i++)
110  for (_j = 0; _j < _phi.size(); _j++)
112 
114 
115  if (_has_diag_save_in)
116  {
117  unsigned int rows = _local_ke.m();
118  DenseVector<Number> diag(rows);
119  for (unsigned int i = 0; i < rows; i++)
120  diag(i) = _local_ke(i, i);
121 
122  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
123  for (unsigned int i = 0; i < _diag_save_in.size(); i++)
124  _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
125  }
126 }
127 
128 void
130 {
131  size_t jvar_num = jvar.number();
132  prepareMatrixTag(_assembly, _var.number(), jvar_num);
133 
134  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
135  for (_i = 0; _i < _test.size(); _i++)
136  for (_j = 0; _j < jvar.phiFaceSize(); _j++)
137  {
138  if (_var.number() == jvar_num)
140  else
141  _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar_num);
142  }
143 
145 }
146 
147 void
148 IntegratedBC::computeJacobianBlock(unsigned int jvar)
149 {
150  mooseDeprecated("The computeJacobianBlock method signature has changed. Developers, please "
151  "pass in a MooseVariableFEBase reference instead of the variable number");
153 
154  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
155  for (_i = 0; _i < _test.size(); _i++)
156  for (_j = 0; _j < _phi.size(); _j++)
157  {
158  if (_var.number() == jvar)
160  else
162  }
163 
165 }
166 
167 void
169 {
171 
173  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
174  for (_i = 0; _i < _test.size(); _i++)
175  for (_j = 0; _j < jv.order(); _j++)
177 
179 }
const VariableTestValue & _test
test function values (in QPs)
Definition: IntegratedBC.h:68
VarFieldType
Definition: MooseTypes.h:488
SubProblem & _subproblem
Reference to SubProblem.
InputParameters validParams< IntegratedBCBase >()
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
std::vector< MooseVariableFEBase * > _diag_save_in
virtual MooseVariable & getStandardVariable(THREAD_ID tid, const std::string &var_name)=0
Returns the variable reference for requested MooseVariable which may be in any system.
std::string incompatVarMsg(MooseVariableFEBase &var1, MooseVariableFEBase &var2)
Builds and returns a string of the form:
Definition: MooseError.C:22
std::vector< MooseVariableFEBase * > _save_in
InputParameters validParams< IntegratedBC >()
Definition: IntegratedBC.C:23
virtual Real computeQpJacobian()
unsigned int number() const
Get variable number coming from libMesh.
IntegratedBC(const InputParameters &parameters)
Definition: IntegratedBC.C:29
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
const VariablePhiValue & _phi
shape function values (in QPs)
Definition: IntegratedBC.h:61
DenseMatrix< Number > _local_ke
Holds residual entries as they are accumulated by this Kernel.
const FEType & feType() const
Get the type of finite element object.
unsigned int _qp
quadrature point index
void addMooseVariableDependency(MooseVariableFEBase *var)
Call this function to add the passed in MooseVariableFEBase as a variable that this object depends on...
virtual void computeJacobian() override
Definition: IntegratedBC.C:104
virtual Real computeQpResidual()=0
method for computing the residual at quadrature points
MooseVariableFE< Real > * mooseVariable() const
Get the variable that this object is using.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:259
std::vector< AuxVariableName > _save_in_strings
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
virtual void computeJacobianBlock(MooseVariableFEBase &jvar) override
Computes d-ivar-residual / d-jvar...
Definition: IntegratedBC.C:129
Assembly & _assembly
Reference to assembly.
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
const MooseArray< Real > & _coord
coordinate transformation
bool _has_save_in
The aux variables to save the residual contributions to.
virtual size_t phiFaceSize() const =0
Return phiFace size.
void paramError(const std::string &param, Args... args)
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
Definition: MooseObject.h:108
MooseVariable & _var
Definition: IntegratedBC.h:53
const QBase *const & _qrule
active quadrature rule
std::vector< AuxVariableName > _diag_save_in_strings
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:178
SystemBase & _sys
Reference to SystemBase.
void mooseDeprecated(Args &&... args) const
Definition: MooseObject.h:161
virtual Real computeQpOffDiagJacobian(unsigned int)
This is the virtual that derived classes should override for computing an off-diagonal jacobian compo...
void computeJacobianBlockScalar(unsigned int jvar) override
Computes jacobian block with respect to a scalar variable.
Definition: IntegratedBC.C:168
DenseVector< Number > _local_re
Holds residual entries as they are accumulated by this Kernel.
Interface for objects that need to get values of MooseVariables.
Class for scalar variables (they are different).
Base class for deriving any boundary condition of a integrated type.
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:184
virtual void computeResidual() override
Definition: IntegratedBC.C:85
Definition: Moose.h:112
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual the according to active tags.
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the ative tags.
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
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:149
const MooseArray< Real > & _JxW
transformed Jacobian weights