www.mooseframework.org
ADKernel.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 "ADKernel.h"
11 #include "Assembly.h"
12 #include "MooseVariable.h"
13 #include "Problem.h"
14 #include "SubProblem.h"
15 #include "NonlinearSystemBase.h"
16 
17 // libmesh includes
18 #include "libmesh/threads.h"
19 
20 defineADBaseValidParams(ADKernel, KernelBase, params.registerBase("Kernel"););
21 defineADBaseValidParams(ADVectorKernel, KernelBase, params.registerBase("VectorKernel"););
22 
23 template <typename T, ComputeStage compute_stage>
25  : KernelBase(parameters),
26  MooseVariableInterface<T>(this,
27  false,
28  "variable",
30  std::is_same<T, Real>::value ? Moose::VarFieldType::VAR_FIELD_STANDARD
32  _var(*this->mooseVariable()),
33  _test(_var.phi()),
34  _grad_test(_var.template adGradPhi<compute_stage>()),
35  _u(_var.template adSln<compute_stage>()),
36  _grad_u(_var.template adGradSln<compute_stage>()),
37  _ad_JxW(_assembly.adJxW<compute_stage>()),
38  _ad_coord(_assembly.adCoordTransformation<compute_stage>()),
39  _ad_q_point(_assembly.template adQPoints<compute_stage>()),
40  _phi(_assembly.phi(_var)),
41  _grad_phi(_assembly.template adGradPhi<T, compute_stage>(_var))
42 {
44  _save_in.resize(_save_in_strings.size());
45  _diag_save_in.resize(_diag_save_in_strings.size());
46 
47  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
48  {
50 
52  paramError("save_in", "cannot use solution variable as save-in variable");
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 
60  _save_in[i] = var;
63  }
64 
65  _has_save_in = _save_in.size() > 0;
66 
67  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
68  {
70 
72  paramError("diag_save_in", "cannot use solution variable as diag save-in variable");
73 
74  if (var->feType() != _var.feType())
75  paramError(
76  "diag_save_in",
77  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
79 
80  _diag_save_in[i] = var;
83  }
84 
85  _has_diag_save_in = _diag_save_in.size() > 0;
86 }
87 
88 template <typename T, ComputeStage compute_stage>
89 void
91 {
92  prepareVectorTag(_assembly, _var.number());
93 
94  precalculateResidual();
95  for (_i = 0; _i < _test.size(); _i++)
96  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
97  _local_re(_i) += _ad_JxW[_qp] * _ad_coord[_qp] * computeQpResidual();
98 
99  accumulateTaggedLocalResidual();
100 
101  if (_has_save_in)
102  {
103  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
104  for (unsigned int i = 0; i < _save_in.size(); i++)
105  _save_in[i]->sys().solution().add_vector(_local_re, _save_in[i]->dofIndices());
106  }
107 }
108 
109 template <>
110 void
112 {
113 }
114 
115 template <>
116 void
118 {
119 }
120 
121 template <typename T, ComputeStage compute_stage>
122 void
124 {
125  prepareMatrixTag(_assembly, _var.number(), _var.number());
126 
127  size_t ad_offset = _var.number() * _sys.getMaxVarNDofsPerElem();
128 
129  precalculateResidual();
130  for (_i = 0; _i < _test.size(); _i++)
131  {
132  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
133  {
134  DualReal residual = _ad_JxW[_qp] * _ad_coord[_qp] * computeQpResidual();
135  for (_j = 0; _j < _var.phiSize(); _j++)
136  _local_ke(_i, _j) += residual.derivatives()[ad_offset + _j];
137  }
138  }
139 
140  accumulateTaggedLocalMatrix();
141 
142  if (_has_diag_save_in)
143  {
144  unsigned int rows = _local_ke.m();
145  DenseVector<Number> diag(rows);
146  for (unsigned int i = 0; i < rows; i++)
147  diag(i) = _local_ke(i, i);
148 
149  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
150  for (unsigned int i = 0; i < _diag_save_in.size(); i++)
151  _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
152  }
153 }
154 
155 template <>
156 void
158 {
159 }
160 template <>
161 void
163 {
164 }
165 
166 template <typename T, ComputeStage compute_stage>
167 void
169 {
170  std::vector<DualReal> residuals(_test.size(), 0);
171 
172  precalculateResidual();
173  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
174  for (_i = 0; _i < _test.size(); _i++)
175  residuals[_i] += _ad_JxW[_qp] * _ad_coord[_qp] * computeQpResidual();
176 
177  std::vector<std::pair<MooseVariableFEBase *, MooseVariableFEBase *>> & ce =
178  _assembly.couplingEntries();
179  for (const auto & it : ce)
180  {
181  MooseVariableFEBase & ivariable = *(it.first);
182  MooseVariableFEBase & jvariable = *(it.second);
183 
184  unsigned int ivar = ivariable.number();
185  unsigned int jvar = jvariable.number();
186 
187  if (ivar != _var.number())
188  continue;
189 
190  size_t ad_offset = jvar * _sys.getMaxVarNDofsPerElem();
191 
192  prepareMatrixTag(_assembly, ivar, jvar);
193 
194  if (_local_ke.m() != _test.size() || _local_ke.n() != jvariable.phiSize())
195  continue;
196 
197  precalculateResidual();
198  for (_i = 0; _i < _test.size(); _i++)
199  for (_j = 0; _j < jvariable.phiSize(); _j++)
200  _local_ke(_i, _j) += residuals[_i].derivatives()[ad_offset + _j];
201 
202  accumulateTaggedLocalMatrix();
203  }
204 }
205 
206 template <typename T, ComputeStage compute_stage>
207 void
209 {
210 }
211 
212 template class ADKernelTempl<Real, RESIDUAL>;
213 template class ADKernelTempl<Real, JACOBIAN>;
VarFieldType
Definition: MooseTypes.h:488
FEProblemBase & _fe_problem
Reference to this kernel&#39;s FEProblemBase.
Definition: KernelBase.h:130
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::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
virtual void computeOffDiagJacobianScalar(unsigned int jvar) override
Computes jacobian block with respect to a scalar variable.
Definition: ADKernel.C:208
NonlinearSystemBase & getNonlinearSystemBase()
std::vector< AuxVariableName > _diag_save_in_strings
Definition: KernelBase.h:178
unsigned int number() const
Get variable number coming from libMesh.
defineADBaseValidParams(ADKernel, KernelBase, params.registerBase("Kernel");)
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...
DualNumber< Real, NumberArray< AD_MAX_DOFS_PER_ELEM, Real > > DualReal
Definition: DualReal.h:29
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: KernelBase.h:176
virtual void computeADOffDiagJacobian() override
Definition: ADKernel.C:168
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...
virtual size_t phiSize() const =0
Return phi size.
MooseVariableFE< T > * mooseVariable() const
Get the variable that this object is using.
ADKernelTempl(const InputParameters &parameters)
Definition: ADKernel.C:24
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
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: KernelBase.h:171
MooseVariableFE< T > & _var
This is a regular kernel so we cast to a regular MooseVariable.
Definition: ADKernel.h:91
virtual bool hasVariable(const std::string &var_name) const
Query a system for a variable.
Definition: SystemBase.C:668
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
std::vector< AuxVariableName > _save_in_strings
Definition: KernelBase.h:173
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
Interface for objects that need to get values of MooseVariables.
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 computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: ADKernel.C:123
Definition: Moose.h:112
SystemBase & sys()
Get the system this variable is part of.
virtual void computeResidual() override
Compute this Kernel&#39;s contribution to the residual.
Definition: ADKernel.C:90