20 #include "libmesh/numeric_vector.h" 21 #include "libmesh/dof_map.h" 22 #include "libmesh/quadrature.h" 23 #include "libmesh/boundary_info.h" 25 template <
typename ComputeValueType>
40 template <
typename ComputeValueType>
45 parameters.getCheckedPointerParam<
SystemBase *>(
"_sys")
47 parameters.
get<AuxVariableName>(
"variable"))
57 _var(_aux_sys.getActualFieldVariable<ComputeValueType>(
58 _tid, parameters.
get<AuxVariableName>(
"variable"))),
59 _nodal(_var.isNodal()),
60 _u(_nodal ? _var.nodalValueArray() : _var.sln()),
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()),
68 _current_elem(_assembly.elem()),
69 _current_side(_assembly.side()),
70 _current_elem_volume(_assembly.elemVolume()),
71 _current_side_volume(_assembly.sideElemVolume()),
73 _current_node(_assembly.node()),
74 _current_boundary_id(_assembly.currentBoundaryID()),
75 _solution(_aux_sys.solution()),
77 _current_lower_d_elem(_assembly.lowerDElem())
81 template <
typename ComputeValueType>
85 auto var = getVar(var_name, comp);
89 ": Unable to couple time derivative of an auxiliary variable into the auxiliary system.");
94 template <
typename ComputeValueType>
97 unsigned int comp)
const 99 auto var = getVar(var_name, comp);
103 ": Unable to couple time derivative of an auxiliary variable into the auxiliary system.");
113 "Should only be calling setDofValue if there is one dof for the aux var");
124 template <
typename ComputeValueType>
128 if (_coincident_lower_d_calc)
129 _var.insertLower(_aux_sys.solution());
131 _var.insert(_aux_sys.solution());
134 template <
typename ComputeValueType>
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())
148 ComputeValueType
value = computeValue();
150 _var.setNodalValue(
value);
155 _n_local_dofs = _coincident_lower_d_calc ? _var.dofIndicesLower().size() : _var.numberOfDofs();
157 if (_coincident_lower_d_calc)
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);
164 mooseError(
"No degrees of freedom. ", lower_error);
167 if (_n_local_dofs == 1)
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);
174 setDofValueHelper(
value);
180 if (_coincident_lower_d_calc)
182 _local_sol.resize(1);
183 if constexpr (std::is_same<Real, ComputeValueType>::value)
184 _local_sol(0) =
value;
186 mooseAssert(
false,
"We should not enter the single dof branch with a vector variable");
187 _var.setLowerDofValues(_local_sol);
190 _var.setNodalValue(
value);
195 _local_re.resize(_n_local_dofs);
197 _local_ke.resize(_n_local_dofs, _n_local_dofs);
200 const auto & test = _coincident_lower_d_calc ? _var.phiLower() : _test;
203 for (
unsigned int i = 0; i < test.size(); i++)
204 for (_qp = 0; _qp < _qrule->n_points(); _qp++)
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];
212 _local_sol.resize(_n_local_dofs);
214 _local_ke.svd_solve(_local_re, _local_sol);
216 _local_ke.cholesky_solve(_local_re, _local_sol);
218 _coincident_lower_d_calc ? _var.setLowerDofValues(_local_sol) : _var.setDofValues(_local_sol);
231 if (_var.isNodalDefined())
236 _var.setNodalValue(
value);
241 _n_local_dofs = _var.numberOfDofs();
242 if (_n_local_dofs == 1)
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);
251 _var.setNodalValue(
value);
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);
262 for (
unsigned int i = 0; i < _test.size(); i++)
263 for (_qp = 0; _qp < _qrule->n_points(); _qp++)
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];
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());
277 for (
unsigned int i = 0; i < _var.count(); ++i)
279 for (
unsigned int j = 0; j < _n_local_dofs; ++j)
280 re(j) = _local_re(j)(i);
283 _local_ke.svd_solve(re, sol);
285 _local_ke.cholesky_solve(re, sol);
287 for (
unsigned int j = 0; j < _n_local_dofs; ++j)
288 _local_sol(j)(i) = sol(j);
291 _var.setDofValues(_local_sol);
296 template <
typename ComputeValueType>
300 if (_sys.solutionStatesInitialized())
301 mooseError(
"The solution states have already been initialized when calling ",
304 "Make sure to call uOld() within the object constructor.");
306 return _nodal ? _var.nodalValueOldArray() : _var.slnOld();
309 template <
typename ComputeValueType>
313 if (_sys.solutionStatesInitialized())
314 mooseError(
"The solution states have already been initialized when calling ",
317 "Make sure to call uOlder() within the object constructor.");
319 return _nodal ? _var.nodalValueOlderArray() : _var.slnOlder();
322 template <
typename ComputeValueType>
std::string name(const ElemQuality q)
Base class for auxiliary kernels.
const OutputTools< ComputeValueType >::VariableValue & uOlder() const
Retrieves the older value of the variable that this AuxKernel operates on.
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
T * get(const std::unique_ptr< T > &u)
The MooseUtils::get() specializations are used to support making forwards-compatible code changes fro...
virtual void compute() override
Computes the value and stores it in the solution vector.
virtual const VariableValue & coupledDot(const std::string &var_name, unsigned int comp=0) const
Time derivative of a coupled variable.
virtual void setDofValue(const OutputData &value, unsigned int index)=0
Degree of freedom value setters.
Base class for a system (of equations)
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
VarKindType
Framework-wide stuff.
AuxKernelTempl(const InputParameters ¶meters)
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.
static InputParameters validParams()
OutputTools< Real >::VariableValue VariableValue
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.
const OutputTools< ComputeValueType >::VariableValue & uOld() const
Retrieves the old value of the variable that this AuxKernel operates on.
void insert()
Insert the just computed values into the auxiliary solution vector.
MooseVariableField< ComputeValueType > & _var
This is a regular kernel so we cast to a regular MooseVariable, hides base _var.
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...
static InputParameters validParams()
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealEigenVector
Base class for creating new auxiliary kernels and auxiliary boundary conditions.
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.
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