https://mooseframework.inl.gov
AuxKernel.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 "AuxKernel.h"
11 
12 // local includes
13 #include "FEProblem.h"
14 #include "SubProblem.h"
15 #include "AuxiliarySystem.h"
16 #include "MooseTypes.h"
17 #include "Assembly.h"
18 #include "MortarNodalAuxKernel.h"
19 
20 #include "libmesh/numeric_vector.h"
21 #include "libmesh/dof_map.h"
22 #include "libmesh/quadrature.h"
23 #include "libmesh/boundary_info.h"
24 
25 template <typename ComputeValueType>
28 {
30 
31  params.registerBase("AuxKernel");
32  if (typeid(AuxKernelTempl<ComputeValueType>).name() == typeid(VectorAuxKernel).name())
33  params.registerBase("VectorAuxKernel");
34  if (typeid(AuxKernelTempl<ComputeValueType>).name() == typeid(ArrayAuxKernel).name())
35  params.registerBase("ArrayAuxKernel");
36 
37  return params;
38 }
39 
40 template <typename ComputeValueType>
42  : AuxKernelBase(parameters),
43  MooseVariableInterface<ComputeValueType>(
44  this,
45  parameters.getCheckedPointerParam<SystemBase *>("_sys")
46  ->getVariable(parameters.get<THREAD_ID>("_tid"),
47  parameters.get<AuxVariableName>("variable"))
48  .isNodal(),
49  "variable",
51  std::is_same<Real, ComputeValueType>::value
53  : (std::is_same<RealVectorValue, ComputeValueType>::value
56 
57  _var(_aux_sys.getActualFieldVariable<ComputeValueType>(
58  _tid, parameters.get<AuxVariableName>("variable"))),
59  _nodal(_var.isNodal()),
60  _u(_nodal ? _var.nodalValueArray() : _var.sln()),
61 
62  _test(_bnd ? _var.phiFace() : _var.phi()),
63  _q_point(_bnd ? _assembly.qPointsFace() : _assembly.qPoints()),
64  _qrule(_bnd ? _assembly.qRuleFace() : _assembly.qRule()),
65  _JxW(_bnd ? _assembly.JxWFace() : _assembly.JxW()),
66  _coord(_assembly.coordTransformation()),
67 
68  _current_elem(_assembly.elem()),
69  _current_side(_assembly.side()),
70  _current_elem_volume(_assembly.elemVolume()),
71  _current_side_volume(_assembly.sideElemVolume()),
72 
73  _current_node(_assembly.node()),
74  _current_boundary_id(_assembly.currentBoundaryID()),
75  _solution(_aux_sys.solution()),
76 
77  _current_lower_d_elem(_assembly.lowerDElem())
78 {
79 }
80 
81 template <typename ComputeValueType>
82 const VariableValue &
83 AuxKernelTempl<ComputeValueType>::coupledDot(const std::string & var_name, unsigned int comp) const
84 {
85  auto var = getVar(var_name, comp);
86  if (var->kind() == Moose::VAR_AUXILIARY)
87  mooseError(
88  name(),
89  ": Unable to couple time derivative of an auxiliary variable into the auxiliary system.");
90 
91  return Coupleable::coupledDot(var_name, comp);
92 }
93 
94 template <typename ComputeValueType>
95 const VariableValue &
97  unsigned int comp) const
98 {
99  auto var = getVar(var_name, comp);
100  if (var->kind() == Moose::VAR_AUXILIARY)
101  mooseError(
102  name(),
103  ": Unable to couple time derivative of an auxiliary variable into the auxiliary system.");
104 
105  return Coupleable::coupledDotDu(var_name, comp);
106 }
107 
108 template <>
109 void
110 AuxKernelTempl<Real>::setDofValueHelper(const Real & value)
111 {
112  mooseAssert(_n_local_dofs == 1,
113  "Should only be calling setDofValue if there is one dof for the aux var");
114  _var.setDofValue(value, 0);
115 }
116 
117 template <>
118 void
120 {
121  mooseError("Not implemented");
122 }
123 
124 template <typename ComputeValueType>
125 void
127 {
128  if (_coincident_lower_d_calc)
129  _var.insertLower(_aux_sys.solution());
130  else
131  _var.insert(_aux_sys.solution());
132 }
133 
134 template <typename ComputeValueType>
135 void
137 {
138  precalculateValue();
139 
140  if (isNodal()) /* nodal variables */
141  {
142  mooseAssert(!_coincident_lower_d_calc,
143  "Nodal evaluations are point evaluations. We don't have to concern ourselves with "
144  "coincidence of lower-d blocks and higher-d faces because they share nodes");
145  if (_var.isNodalDefined())
146  {
147  _qp = 0;
148  ComputeValueType value = computeValue();
149  // update variable data, which is referenced by other kernels, so the value is up-to-date
150  _var.setNodalValue(value);
151  }
152  }
153  else /* elemental variables */
154  {
155  _n_local_dofs = _coincident_lower_d_calc ? _var.dofIndicesLower().size() : _var.numberOfDofs();
156 
157  if (_coincident_lower_d_calc)
158  {
159  static const std::string lower_error = "Make sure that the lower-d variable lives on a "
160  "lower-d block that is a superset of the boundary";
161  if (!_current_lower_d_elem)
162  mooseError("No lower-dimensional element. ", lower_error);
163  if (!_n_local_dofs)
164  mooseError("No degrees of freedom. ", lower_error);
165  }
166 
167  if (_n_local_dofs == 1) /* p0 */
168  {
169  ComputeValueType value = 0;
170  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
171  value += _JxW[_qp] * _coord[_qp] * computeValue();
172  value /= (_bnd ? _current_side_volume : _current_elem_volume);
173  if (_var.isFV())
174  setDofValueHelper(value);
175  else
176  {
177  // update the variable data referenced by other kernels.
178  // Note that this will update the values at the quadrature points too
179  // (because this is an Elemental variable)
180  if (_coincident_lower_d_calc)
181  {
182  _local_sol.resize(1);
183  if constexpr (std::is_same<Real, ComputeValueType>::value)
184  _local_sol(0) = value;
185  else
186  mooseAssert(false, "We should not enter the single dof branch with a vector variable");
187  _var.setLowerDofValues(_local_sol);
188  }
189  else
190  _var.setNodalValue(value);
191  }
192  }
193  else /* high-order */
194  {
195  _local_re.resize(_n_local_dofs);
196  _local_re.zero();
197  _local_ke.resize(_n_local_dofs, _n_local_dofs);
198  _local_ke.zero();
199 
200  const auto & test = _coincident_lower_d_calc ? _var.phiLower() : _test;
201 
202  // assemble the local mass matrix and the load
203  for (unsigned int i = 0; i < test.size(); i++)
204  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
205  {
206  ComputeValueType t = _JxW[_qp] * _coord[_qp] * test[i][_qp];
207  _local_re(i) += t * computeValue();
208  for (unsigned int j = 0; j < test.size(); j++)
209  _local_ke(i, j) += t * test[j][_qp];
210  }
211  // mass matrix is always SPD but in case of boundary restricted, it will be rank deficient
212  _local_sol.resize(_n_local_dofs);
213  if (_bnd)
214  _local_ke.svd_solve(_local_re, _local_sol);
215  else
216  _local_ke.cholesky_solve(_local_re, _local_sol);
217 
218  _coincident_lower_d_calc ? _var.setLowerDofValues(_local_sol) : _var.setDofValues(_local_sol);
219  }
220  }
221 }
222 
223 template <>
224 void
226 {
227  precalculateValue();
228 
229  if (isNodal()) /* nodal variables */
230  {
231  if (_var.isNodalDefined())
232  {
233  _qp = 0;
234  RealEigenVector value = computeValue();
235  // update variable data, which is referenced by other kernels, so the value is up-to-date
236  _var.setNodalValue(value);
237  }
238  }
239  else /* elemental variables */
240  {
241  _n_local_dofs = _var.numberOfDofs();
242  if (_n_local_dofs == 1) /* p0 */
243  {
244  RealEigenVector value = RealEigenVector::Zero(_var.count());
245  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
246  value += _JxW[_qp] * _coord[_qp] * computeValue();
247  value /= (_bnd ? _current_side_volume : _current_elem_volume);
248  // update the variable data referenced by other kernels.
249  // Note that this will update the values at the quadrature points too
250  // (because this is an Elemental variable)
251  _var.setNodalValue(value);
252  }
253  else /* high-order */
254  {
255  _local_re.resize(_n_local_dofs);
256  for (unsigned int i = 0; i < _local_re.size(); ++i)
257  _local_re(i) = RealEigenVector::Zero(_var.count());
258  _local_ke.resize(_n_local_dofs, _n_local_dofs);
259  _local_ke.zero();
260 
261  // assemble the local mass matrix and the load
262  for (unsigned int i = 0; i < _test.size(); i++)
263  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
264  {
265  Real t = _JxW[_qp] * _coord[_qp] * _test[i][_qp];
266  _local_re(i) += t * computeValue();
267  for (unsigned int j = 0; j < _test.size(); j++)
268  _local_ke(i, j) += t * _test[j][_qp];
269  }
270 
271  // mass matrix is always SPD
272  _local_sol.resize(_n_local_dofs);
273  for (unsigned int i = 0; i < _local_re.size(); ++i)
274  _local_sol(i) = RealEigenVector::Zero(_var.count());
275  DenseVector<Number> re(_n_local_dofs);
276  DenseVector<Number> sol(_n_local_dofs);
277  for (unsigned int i = 0; i < _var.count(); ++i)
278  {
279  for (unsigned int j = 0; j < _n_local_dofs; ++j)
280  re(j) = _local_re(j)(i);
281 
282  if (_bnd)
283  _local_ke.svd_solve(re, sol);
284  else
285  _local_ke.cholesky_solve(re, sol);
286 
287  for (unsigned int j = 0; j < _n_local_dofs; ++j)
288  _local_sol(j)(i) = sol(j);
289  }
290 
291  _var.setDofValues(_local_sol);
292  }
293  }
294 }
295 
296 template <typename ComputeValueType>
299 {
300  if (_sys.solutionStatesInitialized())
301  mooseError("The solution states have already been initialized when calling ",
302  type(),
303  "::uOld().\n\n",
304  "Make sure to call uOld() within the object constructor.");
305 
306  return _nodal ? _var.nodalValueOldArray() : _var.slnOld();
307 }
308 
309 template <typename ComputeValueType>
312 {
313  if (_sys.solutionStatesInitialized())
314  mooseError("The solution states have already been initialized when calling ",
315  type(),
316  "::uOlder().\n\n",
317  "Make sure to call uOlder() within the object constructor.");
318 
319  return _nodal ? _var.nodalValueOlderArray() : _var.slnOlder();
320 }
321 
322 template <typename ComputeValueType>
323 bool
325 {
326  return dynamic_cast<MortarNodalAuxKernelTempl<ComputeValueType> *>(this) != nullptr;
327 }
328 
329 // Explicitly instantiates the three versions of the AuxKernelTempl class
330 template class AuxKernelTempl<Real>;
331 template class AuxKernelTempl<RealVectorValue>;
332 template class AuxKernelTempl<RealEigenVector>;
std::string name(const ElemQuality q)
VarFieldType
Definition: MooseTypes.h:722
Base class for auxiliary kernels.
Definition: AuxKernelBase.h:42
const OutputTools< ComputeValueType >::VariableValue & uOlder() const
Retrieves the older value of the variable that this AuxKernel operates on.
Definition: AuxKernel.C:311
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
T * get(const std::unique_ptr< T > &u)
The MooseUtils::get() specializations are used to support making forwards-compatible code changes fro...
Definition: MooseUtils.h:1135
virtual void compute() override
Computes the value and stores it in the solution vector.
Definition: AuxKernel.C:136
virtual const VariableValue & coupledDot(const std::string &var_name, unsigned int comp=0) const
Time derivative of a coupled variable.
Definition: Coupleable.C:1170
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual void setDofValue(const OutputData &value, unsigned int index)=0
Degree of freedom value setters.
Base class for a system (of equations)
Definition: SystemBase.h:84
void registerBase(const std::string &value)
This method must be called from every base "Moose System" to create linkage with the Action System...
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:715
bool isMortar()
Definition: AuxKernel.C:324
AuxKernelTempl(const InputParameters &parameters)
Definition: AuxKernel.C:41
virtual const OutputTools< ComputeValueType >::VariableValue & value()
The value of the variable this object is operating on.
Base class for creating new nodally-based mortar auxiliary kernels.
virtual const VariableValue & coupledDot(const std::string &var_name, unsigned int comp=0) const override
Time derivative of a coupled variable.
Definition: AuxKernel.C:83
static InputParameters validParams()
Definition: AuxKernelBase.C:20
OutputTools< Real >::VariableValue VariableValue
Definition: MooseTypes.h:315
virtual const VariableValue & coupledDotDu(const std::string &var_name, unsigned int comp=0) const
Time derivative of a coupled variable with respect to the coefficients.
Definition: Coupleable.C:1468
const OutputTools< ComputeValueType >::VariableValue & uOld() const
Retrieves the old value of the variable that this AuxKernel operates on.
Definition: AuxKernel.C:298
void insert()
Insert the just computed values into the auxiliary solution vector.
Definition: AuxKernel.C:126
MooseVariableField< ComputeValueType > & _var
This is a regular kernel so we cast to a regular MooseVariable, hides base _var.
Definition: AuxKernel.h:103
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Interface for objects that need to get values of MooseVariables.
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type and optionally a file path to the top-level block p...
Definition: MooseBase.h:271
static InputParameters validParams()
Definition: AuxKernel.C:27
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealEigenVector
Definition: MooseTypes.h:146
Base class for creating new auxiliary kernels and auxiliary boundary conditions.
Definition: AuxKernel.h:17
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
virtual const VariableValue & coupledDotDu(const std::string &var_name, unsigned int comp=0) const override
Time derivative of a coupled variable with respect to the coefficients.
Definition: AuxKernel.C:96
void setDofValueHelper(const ComputeValueType &dof_value)
Currently only used when the auxiliary variable is a finite volume variable, this helps call through ...
unsigned int _n_local_dofs
number of local dofs for elemental variables
Definition: AuxKernel.h:148
unsigned int THREAD_ID
Definition: MooseTypes.h:209