https://mooseframework.inl.gov
ArrayKernel.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
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 "ArrayKernel.h"
11 
12 #include "Assembly.h"
13 #include "MooseVariableFE.h"
14 #include "MooseVariableScalar.h"
15 #include "SubProblem.h"
16 #include "NonlinearSystem.h"
17 #include "FEProblemBase.h"
18 
19 #include "libmesh/threads.h"
20 #include "libmesh/quadrature.h"
21 
24 {
26  params.registerBase("ArrayKernel");
27  return params;
28 }
29 
31  : KernelBase(parameters),
33  false,
34  "variable",
37  _var(*mooseVariable()),
38  _test(_var.phi()),
39  _grad_test(_var.gradPhi()),
40  _array_grad_test(_var.arrayGradPhi()),
41  _phi(_assembly.phi(_var)),
42  _grad_phi(_assembly.gradPhi(_var)),
43  _u(_is_implicit ? _var.sln() : _var.slnOld()),
44  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld()),
45  _count(_var.count()),
46  _work_vector(_count)
47 {
49 
50  _save_in.resize(_save_in_strings.size());
51  _diag_save_in.resize(_diag_save_in_strings.size());
52 
53  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
54  {
56 
58  paramError("save_in", "cannot use solution variable as save-in variable");
59 
60  if (var->feType() != _var.feType())
61  paramError(
62  "save_in",
63  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
65 
66  _save_in[i] = var;
69  }
70 
71  _has_save_in = _save_in.size() > 0;
72 
73  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
74  {
76 
78  paramError("diag_save_in", "cannot use solution variable as diag save-in variable");
79 
80  if (var->feType() != _var.feType())
81  paramError(
82  "diag_save_in",
83  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
85 
86  _diag_save_in[i] = var;
89  }
90 
91  _has_diag_save_in = _diag_save_in.size() > 0;
92 }
93 
94 void
96 {
98 
100  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
101  {
102  initQpResidual();
103  for (_i = 0; _i < _test.size(); _i++)
104  {
105  _work_vector.setZero();
107  mooseAssert(_work_vector.size() == _count,
108  "Size of local residual is not equal to the number of array variable compoments");
109  _work_vector *= _JxW[_qp] * _coord[_qp];
111  }
112  }
113 
115 
116  if (_has_save_in)
117  {
119  for (const auto & var : _save_in)
120  {
121  auto * avar = dynamic_cast<ArrayMooseVariable *>(var);
122  if (avar)
123  avar->addSolution(_local_re);
124  else
125  mooseError("Save-in variable for an array kernel must be an array variable");
126  }
127  }
128 }
129 
130 void
132 {
134 
136  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
137  {
138  initQpJacobian();
139  for (_i = 0; _i < _test.size(); _i++)
140  for (_j = 0; _j < _phi.size(); _j++)
141  {
145  }
146  }
147 
149 
150  if (_has_diag_save_in)
151  {
154  for (const auto & var : _diag_save_in)
155  {
156  auto * avar = dynamic_cast<ArrayMooseVariable *>(var);
157  if (avar)
158  avar->addSolution(diag);
159  else
160  mooseError("Save-in variable for an array kernel must be an array variable");
161  }
162  }
163 }
164 
167 {
168  return RealEigenVector::Zero(_var.count());
169 }
170 
171 void
172 ArrayKernel::computeOffDiagJacobian(const unsigned int jvar_num)
173 {
174  const auto & jvar = getVariable(jvar_num);
175 
176  bool same_var = (jvar_num == _var.number());
177 
178  prepareMatrixTag(_assembly, _var.number(), jvar_num);
179 
180  // This (undisplaced) jvar could potentially yield the wrong phi size if this object is acting on
181  // the displaced mesh
182  auto phi_size = jvar.dofIndices().size();
183 
184  precalculateOffDiagJacobian(jvar_num);
185  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
186  {
187  initQpOffDiagJacobian(jvar);
188  for (_i = 0; _i < _test.size(); _i++)
189  for (_j = 0; _j < phi_size; _j++)
190  {
193  _local_ke, _i, _test.size(), _j, phi_size, _var.number(), jvar_num, _work_matrix);
194  }
195  }
196 
198 
199  if (_has_diag_save_in && same_var)
200  {
203  for (const auto & var : _diag_save_in)
204  {
205  auto * avar = dynamic_cast<ArrayMooseVariable *>(var);
206  if (avar)
207  avar->addSolution(diag);
208  else
209  mooseError("Save-in variable for an array kernel must be an array variable");
210  }
211  }
212 }
213 
216 {
217  if (jvar.number() == _var.number())
218  return computeQpJacobian().asDiagonal();
219  else
220  return RealEigenMatrix::Zero(_var.count(), jvar.count());
221 }
222 
223 void
225 {
228 
229  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
230  for (_i = 0; _i < _test.size(); _i++)
231  {
234  _local_ke, _i, _test.size(), 0, 1, _var.number(), jvar, _work_matrix);
235  }
236 
238 }
239 
242 {
243  return RealEigenMatrix::Zero(_var.count(), (unsigned int)jvar.order() + 1);
244 }
VarFieldType
Definition: MooseTypes.h:721
const libMesh::FEType & feType() const
Get the type of finite element object.
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
std::vector< MooseVariableFEBase * > _save_in
Definition: KernelBase.h:65
std::string incompatVarMsg(MooseVariableFieldBase &var1, MooseVariableFieldBase &var2)
Builds and returns a string of the form:
Definition: MooseError.C:26
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:52
std::vector< AuxVariableName > _diag_save_in_strings
Definition: KernelBase.h:71
virtual void precalculateOffDiagJacobian(unsigned int)
unsigned int number() const
Get variable number coming from libMesh.
void saveDiagLocalArrayJacobian(DenseMatrix< Number > &ke, unsigned int i, unsigned int ntest, unsigned int j, unsigned int nphi, unsigned int ivar, const RealEigenVector &v) const
Helper function for assembling diagonal Jacobian contriubutions on local quadrature points for an arr...
Definition: Assembly.h:1806
void addSolution(const DenseVector< libMesh::Number > &v)
Add passed in local DOF values onto the current solution.
unsigned int count() const
Get the number of components Note: For standard and vector variables, the number is one...
virtual void precalculateResidual()
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual void initQpOffDiagJacobian(const MooseVariableFEBase &jvar)
Put necessary evaluations depending on qp but independent on test and shape functions here for off-di...
Definition: ArrayKernel.h:78
This class provides an interface for common operations on field variables of both FE and FV types wit...
void saveLocalArrayResidual(DenseVector< Number > &re, unsigned int i, unsigned int ntest, const RealEigenVector &v) const
Helper function for assembling residual contriubutions on local quadrature points for an array kernel...
Definition: Assembly.h:1787
THREAD_ID _tid
The thread ID for this kernel.
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: KernelBase.h:69
DenseVector< Real > getJacobianDiagonal(DenseMatrix< Number > &ke)
Definition: Assembly.h:1857
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
const MooseVariableFieldBase & getVariable(unsigned int jvar_num) const
Retrieve the variable object from our system associated with jvar_num.
ArrayKernel(const InputParameters &parameters)
Definition: ArrayKernel.C:30
virtual void initQpResidual()
Put necessary evaluations depending on qp but independent on test functions here. ...
Definition: ArrayKernel.h:67
static InputParameters validParams()
Definition: ArrayKernel.C:23
SystemBase & _sys
Reference to the EquationSystem object.
void registerBase(const std::string &value)
This method must be called from every base "Moose System" to create linkage with the Action System...
MooseVariableFE< RealEigenVector > * mooseVariable() const
Return the MooseVariableFE object that this interface acts on.
virtual RealEigenVector computeQpJacobian()
Compute this Kernel&#39;s contribution to the diagonal Jacobian at the current quadrature point...
Definition: ArrayKernel.C:166
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:259
This is the common base class for the three main kernel types implemented in MOOSE, Kernel, VectorKernel and ArrayKernel.
Definition: KernelBase.h:23
virtual RealEigenMatrix computeQpOffDiagJacobianScalar(const MooseVariableScalar &jvar)
This is the virtual that derived classes should override for computing a full Jacobian component...
Definition: ArrayKernel.C:241
virtual void precalculateJacobian()
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: KernelBase.h:70
SubProblem & _subproblem
Reference to this kernel&#39;s SubProblem.
virtual void computeJacobian() override
Compute this ArrayKernel&#39;s contribution to the diagonal Jacobian entries.
Definition: ArrayKernel.C:131
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:714
const QBase *const & _qrule
active quadrature rule
Definition: KernelBase.h:49
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: KernelBase.h:64
FEProblemBase & _fe_problem
Reference to this kernel&#39;s FEProblemBase.
const ArrayVariableTestValue & _test
the current test function
Definition: ArrayKernel.h:88
unsigned int _i
current index for the test function
Definition: KernelBase.h:58
virtual ArrayMooseVariable & getArrayVariable(const THREAD_ID tid, const std::string &var_name)=0
Returns the variable reference for requested ArrayMooseVariable which may be in any system...
virtual void computeQpResidual(RealEigenVector &residual)=0
Compute this Kernel&#39;s contribution to the residual at the current quadrature point, to be filled in residual.
NonlinearSystemBase & getNonlinearSystemBase(const unsigned int sys_num)
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:1159
Eigen::Matrix< Real, Eigen::Dynamic, Eigen::Dynamic > RealEigenMatrix
Definition: MooseTypes.h:149
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
const unsigned int _count
Number of components of the array variable.
Definition: ArrayKernel.h:107
RealEigenVector _work_vector
Work vector for residual and diag jacobian.
Definition: ArrayKernel.h:111
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:55
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
virtual void initQpJacobian()
Put necessary evaluations depending on qp but independent on test and shape functions here...
Definition: ArrayKernel.h:72
virtual RealEigenMatrix computeQpOffDiagJacobian(const MooseVariableFEBase &jvar)
This is the virtual that derived classes should override for computing a full Jacobian component...
Definition: ArrayKernel.C:215
libMesh::Order order() const
Get the order of this variable Note: Order enum can be implicitly converted to unsigned int...
std::vector< AuxVariableName > _save_in_strings
Definition: KernelBase.h:66
unsigned int _j
current index for the shape function
Definition: KernelBase.h:61
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:144
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
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:173
virtual void computeResidual() override
Compute this ArrayKernel&#39;s contribution to the residual.
Definition: ArrayKernel.C:95
const ArrayVariablePhiValue & _phi
the current shape functions
Definition: ArrayKernel.h:95
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
Interface for objects that need to get values of MooseVariables.
virtual void computeOffDiagJacobianScalar(unsigned int jvar) override
Computes jacobian block with respect to a scalar variable.
Definition: ArrayKernel.C:224
Class for scalar variables (they are different).
void mooseError(Args &&... args) const
Emits an error prefixed with object name and 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:179
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealEigenVector
Definition: MooseTypes.h:146
RealEigenMatrix _work_matrix
Work vector for off diag jacobian.
Definition: ArrayKernel.h:114
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual according to active tags.
static InputParameters validParams()
Definition: KernelBase.C:21
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the active tags.
void saveFullLocalArrayJacobian(DenseMatrix< Number > &ke, unsigned int i, unsigned int ntest, unsigned int j, unsigned int nphi, unsigned int ivar, unsigned int jvar, const RealEigenMatrix &v) const
Helper function for assembling full Jacobian contriubutions on local quadrature points for an array k...
Definition: Assembly.h:1831
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes full Jacobian of jvar and the array variable this kernel operates on.
Definition: ArrayKernel.C:172
ArrayMooseVariable & _var
This is an array kernel so we cast to a ArrayMooseVariable.
Definition: ArrayKernel.h:85
SystemBase & sys()
Get the system this variable is part of.
void ErrorVector unsigned int
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:43
spin_mutex spin_mtx