www.mooseframework.org
Kernel.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 "Kernel.h"
11 
12 // MOOSE includes
13 #include "Assembly.h"
14 #include "MooseVariableFE.h"
15 #include "MooseVariableScalar.h"
16 #include "SubProblem.h"
17 #include "NonlinearSystem.h"
18 
19 #include "libmesh/threads.h"
20 #include "libmesh/quadrature.h"
21 
24 {
26  params.registerBase("Kernel");
27  return params;
28 }
29 
30 Kernel::Kernel(const InputParameters & parameters)
31  : KernelBase(parameters),
33  false,
34  "variable",
37  _var(*mooseVariable()),
38  _test(_var.phi()),
39  _grad_test(_var.gradPhi()),
40  _phi(_assembly.phi(_var)),
41  _grad_phi(_assembly.gradPhi(_var)),
42  _u(_is_implicit ? _var.sln() : _var.slnOld()),
43  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld())
44 {
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 
54  paramError("save_in", "cannot use solution variable as save-in variable");
55 
56  if (var->feType() != _var.feType())
57  paramError(
58  "save_in",
59  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
61 
62  _save_in[i] = var;
65  }
66 
67  _has_save_in = _save_in.size() > 0;
68 
69  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
70  {
72 
74  paramError("diag_save_in", "cannot use solution variable as diag save-in variable");
75 
76  if (var->feType() != _var.feType())
77  paramError(
78  "diag_save_in",
79  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
81 
82  _diag_save_in[i] = var;
85  }
86 
87  _has_diag_save_in = _diag_save_in.size() > 0;
88 }
89 
90 void
92 {
93  if (!hasVectorTags())
94  return;
95 
97 
99  for (_i = 0; _i < _test.size(); _i++)
100  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
102 
104 
105  if (_has_save_in)
106  {
107  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
108  for (const auto & var : _save_in)
109  var->sys().solution().add_vector(_local_re, var->dofIndices());
110  }
111 }
112 
113 void
115 {
117 
119  for (_i = 0; _i < _test.size(); _i++)
120  for (_j = 0; _j < _phi.size(); _j++)
121  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
123 
125 
127  {
129  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
130  for (const auto & var : _diag_save_in)
131  var->sys().solution().add_vector(diag, var->dofIndices());
132  }
133 }
134 
135 void
136 Kernel::computeOffDiagJacobian(const unsigned int jvar_num)
137 {
138  const auto & jvar = getVariable(jvar_num);
139 
140  if (jvar_num == _var.number())
141  computeJacobian();
142  else
143  {
144  prepareMatrixTag(_assembly, _var.number(), jvar_num);
145 
146  // This (undisplaced) jvar could potentially yield the wrong phi size if this object is acting
147  // on the displaced mesh
148  auto phi_size = jvar.dofIndices().size();
149  mooseAssert(
150  phi_size * jvar.count() == _local_ke.n(),
151  "The size of the phi container does not match the number of local Jacobian columns");
152 
153  if (_local_ke.m() != _test.size())
154  return;
155 
156  precalculateOffDiagJacobian(jvar_num);
157  if (jvar.count() == 1)
158  {
159  for (_i = 0; _i < _test.size(); _i++)
160  for (_j = 0; _j < phi_size; _j++)
161  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
162  _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar.number());
163  }
164  else
165  {
166  unsigned int n = phi_size;
167  for (_i = 0; _i < _test.size(); _i++)
168  for (_j = 0; _j < n; _j++)
169  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
170  {
171  RealEigenVector v = _JxW[_qp] * _coord[_qp] *
172  computeQpOffDiagJacobianArray(static_cast<ArrayMooseVariable &>(
173  const_cast<MooseVariableFieldBase &>(jvar)));
174  for (unsigned int k = 0; k < v.size(); ++k)
175  _local_ke(_i, _j + k * n) += v(k);
176  }
177  }
178 
180  }
181 }
182 
183 void
184 Kernel::computeOffDiagJacobianScalar(const unsigned int jvar)
185 {
188 
189  for (_i = 0; _i < _test.size(); _i++)
190  for (_j = 0; _j < jv.order(); _j++)
191  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
193 
195 }
196 
197 void
199 {
200  computeResidual();
201 
202  for (const auto & [ivariable, jvariable] : _fe_problem.couplingEntries(_tid, _sys.number()))
203  {
204  const unsigned int ivar = ivariable->number();
205  const unsigned int jvar = jvariable->number();
206 
207  if (ivar != _var.number())
208  continue;
209 
210  if (_is_implicit)
211  {
212  prepareShapes(jvar);
214  }
215  }
216 
218 }
VarFieldType
Definition: MooseTypes.h:634
virtual void computeResidual() override
Compute this Kernel&#39;s contribution to the residual.
Definition: Kernel.C:91
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase * > > & couplingEntries(const THREAD_ID tid, const unsigned int nl_sys_num)
static InputParameters validParams()
Definition: Kernel.C:23
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
virtual Real computeQpResidual()=0
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
bool hasVectorTags() const
virtual RealEigenVector computeQpOffDiagJacobianArray(const ArrayMooseVariable &jvar)
For coupling array variables.
Definition: Kernel.h:66
MooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
Definition: Kernel.h:72
Class for stuff related to variables.
Definition: Adaptivity.h:31
std::vector< MooseVariableFEBase * > _save_in
Definition: KernelBase.h:64
std::string incompatVarMsg(MooseVariableFieldBase &var1, MooseVariableFieldBase &var2)
Builds and returns a string of the form:
Definition: MooseError.C:23
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:51
std::vector< AuxVariableName > _diag_save_in_strings
Definition: KernelBase.h:70
virtual void precalculateOffDiagJacobian(unsigned int)
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpOffDiagJacobian(unsigned int)
For coupling standard variables.
Definition: Kernel.h:56
virtual void precalculateResidual()
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
bool computingScalingJacobian() const
Whether we are computing an initial Jacobian for automatic variable scaling.
Definition: SystemBase.C:1481
virtual void prepareShapes(unsigned int var_num)
Prepare shape functions.
virtual void computeResidualAndJacobian() override
Compute the residual and Jacobian together.
Definition: Kernel.C:198
unsigned int m() const
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:68
DenseVector< Real > getJacobianDiagonal(DenseMatrix< Number > &ke)
Definition: Assembly.h:1840
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.
Kernel(const InputParameters &parameters)
Definition: Kernel.C:30
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< Real > * mooseVariable() const
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:256
const VariableTestValue & _test
the current test function
Definition: Kernel.h:75
This is the common base class for the three main kernel types implemented in MOOSE, Kernel, VectorKernel and ArrayKernel.
Definition: KernelBase.h:22
virtual void precalculateJacobian()
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: KernelBase.h:69
SubProblem & _subproblem
Reference to this kernel&#39;s SubProblem.
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:627
const QBase *const & _qrule
active quadrature rule
Definition: KernelBase.h:48
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: KernelBase.h:63
FEProblemBase & _fe_problem
Reference to this kernel&#39;s FEProblemBase.
virtual Real computeQpJacobian()
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
Definition: Kernel.h:51
unsigned int _i
current index for the test function
Definition: KernelBase.h:57
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.
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:1125
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:54
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
virtual void computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: Kernel.C:114
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:65
unsigned int _j
current index for the shape function
Definition: KernelBase.h:60
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
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
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.
Interface for objects that need to get values of MooseVariables.
Class for scalar variables (they are different).
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
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealEigenVector
Definition: MooseTypes.h:138
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
virtual Real computeQpOffDiagJacobianScalar(unsigned int)
For coupling scalar variables.
Definition: Kernel.h:61
unsigned int n() const
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.
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: Kernel.C:136
const VariablePhiValue & _phi
the current shape functions
Definition: Kernel.h:81
SystemBase & sys()
Get the system this variable is part of.
bool _is_implicit
If the object is using implicit or explicit form.
virtual void computeOffDiagJacobianScalar(unsigned int jvar) override
Computes jacobian block with respect to a scalar variable.
Definition: Kernel.C:184
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:42