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 
23 {
25  return params;
26 }
27 
29  : IntegratedBCBase(parameters),
31  false,
32  "variable",
35  _var(*mooseVariable()),
36  _normals(_assembly.normals()),
37  _phi(_assembly.phiFace(_var)),
38  _grad_phi(_assembly.gradPhiFace(_var)),
39  _test(_var.phiFace()),
40  _grad_test(_var.gradPhiFace()),
41  _u(_is_implicit ? _var.sln() : _var.slnOld()),
42  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld())
43 {
45 
46  _save_in.resize(_save_in_strings.size());
47  _diag_save_in.resize(_diag_save_in_strings.size());
48 
49  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
50  {
52 
53  if (var->feType() != _var.feType())
54  paramError(
55  "save_in",
56  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
58  _save_in[i] = var;
61  }
62 
63  _has_save_in = _save_in.size() > 0;
64 
65  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
66  {
68 
69  if (var->feType() != _var.feType())
70  paramError(
71  "diag_save_in",
72  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
74 
75  _diag_save_in[i] = var;
78  }
79 
80  _has_diag_save_in = _diag_save_in.size() > 0;
81 }
82 
83 void
85 {
87 
89 
90  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
91  {
93  for (_i = 0; _i < _test.size(); _i++)
95  }
96 
98 
99  if (_has_save_in)
100  {
101  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
102  for (unsigned int i = 0; i < _save_in.size(); i++)
103  _save_in[i]->sys().solution().add_vector(_local_re, _save_in[i]->dofIndices());
104  }
105 }
106 
107 void
109 {
111 
113 
114  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
115  {
117  for (_i = 0; _i < _test.size(); _i++)
118  for (_j = 0; _j < _phi.size(); _j++)
120  }
121 
123 
124  if (_has_diag_save_in)
125  {
126  unsigned int rows = _local_ke.m();
127  DenseVector<Number> diag(rows);
128  for (unsigned int i = 0; i < rows; i++)
129  diag(i) = _local_ke(i, i);
130 
131  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
132  for (unsigned int i = 0; i < _diag_save_in.size(); i++)
133  _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
134  }
135 }
136 
137 void
138 IntegratedBC::computeOffDiagJacobian(const unsigned int jvar_num)
139 {
140  const auto & jvar = getVariable(jvar_num);
141 
142  if (jvar_num == _var.number())
143  {
144  computeJacobian();
145  return;
146  }
147 
148  prepareMatrixTag(_assembly, _var.number(), jvar_num);
149 
150  precalculateOffDiagJacobian(jvar_num);
151 
152  // This (undisplaced) jvar could potentially yield the wrong phi size if this object is acting
153  // on the displaced mesh, so we obtain the variable on the proper system
154  auto phi_size = jvar.dofIndices().size();
155  mooseAssert(phi_size * jvar.count() == _local_ke.n(),
156  "The size of the phi container does not match the number of local Jacobian columns");
157 
158  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
159  {
161  for (_i = 0; _i < _test.size(); _i++)
162  for (_j = 0; _j < phi_size; _j++)
163  _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar_num);
164  }
165 
167 }
168 
169 void
171 {
173 
175  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
176  for (_i = 0; _i < _test.size(); _i++)
177  for (_j = 0; _j < jv.order(); _j++)
179 
181 }
182 
183 void
185 {
186  computeResidual();
187 
188  for (const auto & [ivariable, jvariable] : _fe_problem.couplingEntries(_tid, _sys.number()))
189  {
190  const unsigned int ivar = ivariable->number();
191  const unsigned int jvar = jvariable->number();
192 
193  if (ivar != _var.number())
194  continue;
195 
196  if (_is_implicit)
197  {
198  prepareShapes(jvar);
200  }
201  }
202 
204 }
const VariableTestValue & _test
test function values (in QPs)
Definition: IntegratedBC.h:97
VarFieldType
Definition: MooseTypes.h:634
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes d-ivar-residual / d-jvar...
Definition: IntegratedBC.C:138
virtual Real computeQpJacobian()
Method for computing the diagonal Jacobian at quadrature points.
Definition: IntegratedBC.h:49
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase * > > & couplingEntries(const THREAD_ID tid, const unsigned int nl_sys_num)
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
std::vector< MooseVariableFEBase * > _diag_save_in
Class for stuff related to variables.
Definition: Adaptivity.h:31
std::string incompatVarMsg(MooseVariableFieldBase &var1, MooseVariableFieldBase &var2)
Builds and returns a string of the form:
Definition: MooseError.C:23
std::vector< MooseVariableFEBase * > _save_in
virtual void precalculateOffDiagJacobian(unsigned int)
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpOffDiagJacobian(unsigned int)
Method for computing an off-diagonal jacobian component at quadrature points.
Definition: IntegratedBC.h:54
IntegratedBC(const InputParameters &parameters)
Definition: IntegratedBC.C:28
virtual void precalculateResidual()
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
unsigned int m() const
THREAD_ID _tid
The thread ID for this kernel.
const VariablePhiValue & _phi
shape function values (in QPs)
Definition: IntegratedBC.h:90
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
const FEType & feType() const
Get the type of finite element object.
const MooseVariableFieldBase & getVariable(unsigned int jvar_num) const
Retrieve the variable object from our system associated with jvar_num.
unsigned int _qp
quadrature point index
SystemBase & _sys
Reference to the EquationSystem object.
virtual void computeJacobian() override
Compute this object&#39;s contribution to the diagonal Jacobian entries.
Definition: IntegratedBC.C:108
virtual Real computeQpResidual()=0
Method for computing the residual at quadrature points.
MooseVariableFE< Real > * mooseVariable() const
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:256
virtual void precalculateJacobian()
std::vector< AuxVariableName > _save_in_strings
SubProblem & _subproblem
Reference to this kernel&#39;s SubProblem.
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:627
virtual void precalculateQpResidual()
Insertion point for evaluations that depend on qp but are independent of the test functions...
Definition: IntegratedBC.h:68
FEProblemBase & _fe_problem
Reference to this kernel&#39;s FEProblemBase.
void computeOffDiagJacobianScalar(unsigned int jvar) override
Computes jacobian block with respect to a scalar variable.
Definition: IntegratedBC.C:170
virtual MooseVariable & getStandardVariable(const THREAD_ID tid, const std::string &var_name)=0
Returns the variable reference for requested MooseVariable which may be in any system.
virtual Real computeQpOffDiagJacobianScalar(unsigned int jvar)
Method for computing an off-diagonal jacobian component from a scalar var.
Definition: IntegratedBC.h:59
void prepareShapes(unsigned int var_num) override final
Prepare shape functions.
void paramError(const std::string &param, Args... args) const
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
unsigned int number() const
Gets the number of this system.
Definition: SystemBase.C:1125
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.
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
virtual void computeResidualAndJacobian() override
Compute this object&#39;s contribution to the residual and Jacobian simultaneously.
Definition: IntegratedBC.C:184
MooseVariable & _var
Definition: IntegratedBC.h:82
virtual MooseVariableScalar & getScalarVariable(THREAD_ID tid, const std::string &var_name) const
Gets a reference to a scalar variable with specified number.
Definition: SystemBase.C:134
const QBase *const & _qrule
active quadrature rule
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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:163
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
static InputParameters validParams()
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:169
virtual void computeResidual() override
Compute this object&#39;s contribution to the residual.
Definition: IntegratedBC.C:84
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
unsigned int n() const
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual according to active tags.
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the active 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.
bool _is_implicit
If the object is using implicit or explicit form.
const MooseArray< Real > & _JxW
transformed Jacobian weights
virtual void precalculateQpJacobian()
Insertion point for evaluations that depend on qp but are independent of the test and shape functions...
Definition: IntegratedBC.h:74
virtual void precalculateQpOffDiagJacobian(const MooseVariableFEBase &)
Insertion point for evaluations that depend on qp but are independent of the test and shape functions...
Definition: IntegratedBC.h:80