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 
22 template <>
25 {
27  params.registerBase("Kernel");
28  return params;
29 }
30 
31 Kernel::Kernel(const InputParameters & parameters)
32  : KernelBase(parameters),
33  MooseVariableInterface<Real>(this,
34  false,
35  "variable",
38  _var(*mooseVariable()),
39  _test(_var.phi()),
40  _grad_test(_var.gradPhi()),
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 {
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 
55  paramError("save_in", "cannot use solution variable as save-in variable");
56 
57  if (var->feType() != _var.feType())
58  paramError(
59  "save_in",
60  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
62 
63  _save_in[i] = var;
66  }
67 
68  _has_save_in = _save_in.size() > 0;
69 
70  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
71  {
73 
75  paramError("diag_save_in", "cannot use solution variable as diag save-in variable");
76 
77  if (var->feType() != _var.feType())
78  paramError(
79  "diag_save_in",
80  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
82 
83  _diag_save_in[i] = var;
86  }
87 
88  _has_diag_save_in = _diag_save_in.size() > 0;
89 }
90 
91 void
93 {
95 
97  for (_i = 0; _i < _test.size(); _i++)
98  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
100 
102 
103  if (_has_save_in)
104  {
105  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
106  for (const auto & var : _save_in)
107  var->sys().solution().add_vector(_local_re, var->dofIndices());
108  }
109 }
110 
111 void
113 {
115 
117  for (_i = 0; _i < _test.size(); _i++)
118  for (_j = 0; _j < _phi.size(); _j++)
119  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
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 (const auto & var : _diag_save_in)
133  var->sys().solution().add_vector(diag, var->dofIndices());
134  }
135 }
136 
137 void
139 {
140  auto jvar_num = jvar.number();
141  if (jvar_num == _var.number())
142  computeJacobian();
143  else
144  {
145  prepareMatrixTag(_assembly, _var.number(), jvar_num);
146 
147  if (_local_ke.m() != _test.size() || _local_ke.n() != jvar.phiSize())
148  return;
149 
150  precalculateOffDiagJacobian(jvar_num);
151  for (_i = 0; _i < _test.size(); _i++)
152  for (_j = 0; _j < jvar.phiSize(); _j++)
153  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
154  _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar_num);
155 
157  }
158 }
159 
160 void
161 Kernel::computeOffDiagJacobian(unsigned int jvar)
162 {
163  mooseDeprecated("The computeOffDiagJacobian method signature has changed. Developers, please "
164  "pass in a MooseVariableFEBase reference instead of the variable number");
165  if (jvar == _var.number())
166  computeJacobian();
167  else
168  {
170 
172  for (_i = 0; _i < _test.size(); _i++)
173  for (_j = 0; _j < _phi.size(); _j++)
174  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
176 
178  }
179 }
180 
181 void
183 {
186 
187  for (_i = 0; _i < _test.size(); _i++)
188  for (_j = 0; _j < jv.order(); _j++)
189  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
191 
193 }
VarFieldType
Definition: MooseTypes.h:488
virtual Real computeQpJacobian()
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
Definition: KernelBase.h:111
virtual void computeResidual() override
Compute this Kernel&#39;s contribution to the residual.
Definition: Kernel.C:92
FEProblemBase & _fe_problem
Reference to this kernel&#39;s FEProblemBase.
Definition: KernelBase.h:130
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.
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.
MooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
Definition: Kernel.h:54
std::vector< MooseVariableFEBase * > _save_in
Definition: KernelBase.h:172
std::string incompatVarMsg(MooseVariableFEBase &var1, MooseVariableFEBase &var2)
Builds and returns a string of the form:
Definition: MooseError.C:22
NonlinearSystemBase & getNonlinearSystemBase()
virtual void precalculateOffDiagJacobian(unsigned int)
Definition: KernelBase.h:123
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:159
std::vector< AuxVariableName > _diag_save_in_strings
Definition: KernelBase.h:178
unsigned int number() const
Get variable number coming from libMesh.
SystemBase & _sys
Reference to the EquationSystem object.
Definition: KernelBase.h:133
InputParameters validParams< Kernel >()
Definition: Kernel.C:24
THREAD_ID _tid
The thread ID for this kernel.
Definition: KernelBase.h:136
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: KernelBase.h:176
DenseMatrix< Number > _local_ke
Holds residual entries as they are accumulated by this Kernel.
virtual void precalculateJacobian()
Definition: KernelBase.h:122
const FEType & feType() const
Get the type of finite element object.
void addMooseVariableDependency(MooseVariableFEBase *var)
Call this function to add the passed in MooseVariableFEBase as a variable that this object depends on...
Kernel(const InputParameters &parameters)
Definition: Kernel.C:31
virtual void computeOffDiagJacobian(MooseVariableFEBase &jvar) override
Computes d-residual / d-jvar... storing the result in Ke.
Definition: Kernel.C:138
virtual size_t phiSize() const =0
Return phi size.
void registerBase(const std::string &value)
This method must be called from every base "Moose System" to create linkage with the Action System...
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
Definition: KernelBase.h:139
MooseVariableFE< Real > * mooseVariable() const
Get the variable that this object is using.
const VariableTestValue & _test
the current test function
Definition: Kernel.h:57
This is the common base class for the two main kernel types implemented in MOOSE, EigenKernel and Ker...
Definition: KernelBase.h:44
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: KernelBase.h:177
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
const QBase *const & _qrule
active quadrature rule
Definition: KernelBase.h:156
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: KernelBase.h:171
unsigned int _i
current index for the test function
Definition: KernelBase.h:165
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
virtual bool hasVariable(const std::string &var_name) const
Query a system for a variable.
Definition: SystemBase.C:668
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:162
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
virtual void computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: Kernel.C:112
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:173
unsigned int _j
current index for the shape function
Definition: KernelBase.h:168
virtual void precalculateResidual()
Following methods are used for Kernels that need to perform a per-element calculation.
Definition: KernelBase.h:121
SubProblem & _subproblem
Reference to this kernel&#39;s SubProblem.
Definition: KernelBase.h:127
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
void mooseDeprecated(Args &&... args) const
Definition: MooseObject.h:161
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).
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
Definition: Moose.h:112
virtual Real computeQpOffDiagJacobian(unsigned int)
This is the virtual that derived classes should override for computing an off-diagonal Jacobian compo...
Definition: KernelBase.h:116
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.
const VariablePhiValue & _phi
the current shape functions
Definition: Kernel.h:63
SystemBase & sys()
Get the system this variable is part of.
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
InputParameters validParams< KernelBase >()
Definition: KernelBase.C:22
virtual void computeOffDiagJacobianScalar(unsigned int jvar) override
Computes jacobian block with respect to a scalar variable.
Definition: Kernel.C:182
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:150