https://mooseframework.inl.gov
ArrayDGKernel.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 "ArrayDGKernel.h"
11 #include "Assembly.h"
12 #include "MooseVariable.h"
13 #include "Problem.h"
14 #include "SubProblem.h"
15 #include "SystemBase.h"
16 #include "MaterialData.h"
17 #include "ParallelUniqueId.h"
18 
19 #include "libmesh/dof_map.h"
20 #include "libmesh/dense_vector.h"
21 #include "libmesh/numeric_vector.h"
22 #include "libmesh/dense_subvector.h"
23 #include "libmesh/libmesh_common.h"
24 #include "libmesh/quadrature.h"
25 
28 {
30  return params;
31 }
32 
34  : DGKernelBase(parameters),
37  _var(*mooseVariable()),
38  _u(_is_implicit ? _var.sln() : _var.slnOld()),
39  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld()),
40 
41  _phi(_var.phiFace()),
42  _grad_phi(_var.gradPhiFace()),
43  _array_grad_phi(_var.arrayGradPhiFace()),
44 
45  _test(_var.phiFace()),
46  _grad_test(_var.gradPhiFace()),
47  _array_grad_test(_var.arrayGradPhiFace()),
48 
49  _phi_neighbor(_var.phiFaceNeighbor()),
50  _grad_phi_neighbor(_var.gradPhiFaceNeighbor()),
51  _array_grad_phi_neighbor(_var.arrayGradPhiFaceNeighbor()),
52 
53  _test_neighbor(_var.phiFaceNeighbor()),
54  _grad_test_neighbor(_var.gradPhiFaceNeighbor()),
55  _array_grad_test_neighbor(_var.arrayGradPhiFaceNeighbor()),
56 
57  _u_neighbor(_is_implicit ? _var.slnNeighbor() : _var.slnOldNeighbor()),
58  _grad_u_neighbor(_is_implicit ? _var.gradSlnNeighbor() : _var.gradSlnOldNeighbor()),
59 
60  _array_normals(_assembly.mappedNormals()),
61  _count(_var.count()),
62 
63  _work_vector(_count)
64 {
66 
67  _save_in.resize(_save_in_strings.size());
68  _diag_save_in.resize(_diag_save_in_strings.size());
69 
70  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
71  {
76 
78  mooseError("Trying to use solution variable " + _save_in_strings[i] +
79  " as a save_in variable in " + name());
80 
81  if (var->feType() != _var.feType())
82  paramError(
83  "save_in",
84  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
86 
87  _save_in[i] = var;
90  }
91 
92  _has_save_in = _save_in.size() > 0;
93 
94  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
95  {
100 
102  mooseError("Trying to use solution variable " + _diag_save_in_strings[i] +
103  " as a diag_save_in variable in " + name());
104 
105  if (var->feType() != _var.feType())
106  paramError(
107  "diag_save_in",
108  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
110 
111  _diag_save_in[i] = var;
114  }
115 
116  _has_diag_save_in = _diag_save_in.size() > 0;
117 }
118 
119 void
121 {
122  bool is_elem;
123  if (type == Moose::Element)
124  is_elem = true;
125  else
126  is_elem = false;
127 
128  const ArrayVariableTestValue & test_space = is_elem ? _test : _test_neighbor;
129 
130  if (is_elem)
132  else
134 
135  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
136  {
138  for (_i = 0; _i < test_space.size(); _i++)
139  {
140  _work_vector.setZero();
142  mooseAssert(_work_vector.size() == _count,
143  "Size of local residual is not equal to the number of array variable compoments");
144  _work_vector *= _JxW[_qp] * _coord[_qp];
146  }
147  }
148 
150 
151  if (_has_save_in)
152  {
153  Threads::spin_mutex::scoped_lock lock(_resid_vars_mutex);
154  for (const auto & var : _save_in)
155  {
156  auto * avar = dynamic_cast<ArrayMooseVariable *>(var);
157  if (!avar)
158  mooseError("Save-in variable for an array kernel must be an array variable");
159 
160  if (is_elem)
161  avar->addSolution(_local_re);
162  else
163  avar->addSolutionNeighbor(_local_re);
164  }
165  }
166 }
167 
168 void
170 {
171  const ArrayVariableTestValue & test_space =
173  const ArrayVariableTestValue & loc_phi =
175 
178  else
180 
181  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
182  {
184  for (_i = 0; _i < test_space.size(); _i++)
185  for (_j = 0; _j < loc_phi.size(); _j++)
186  {
189  _local_ke, _i, test_space.size(), _j, loc_phi.size(), _var.number(), v);
190  }
191  }
192 
194 
196  {
198  Threads::spin_mutex::scoped_lock lock(_jacoby_vars_mutex);
199  for (const auto & var : _diag_save_in)
200  {
201  auto * avar = dynamic_cast<ArrayMooseVariable *>(var);
202  if (!avar)
203  mooseError("Save-in variable for an array kernel must be an array variable");
204 
206  avar->addSolution(diag);
207  else
208  avar->addSolutionNeighbor(diag);
209  }
210  }
211 }
212 
214 {
215  return RealEigenVector::Zero(_count);
216 }
217 
218 void
219 ArrayDGKernel::computeOffDiagJacobian(const unsigned int jvar_num)
220 {
221  if (!excludeBoundary())
222  {
223  const auto & jvar = getVariable(jvar_num);
224 
225  precalculateOffDiagJacobian(jvar_num);
226 
227  // Compute element-element Jacobian
229 
230  // Compute element-neighbor Jacobian
232 
233  // Compute neighbor-element Jacobian
235 
236  // Compute neighbor-neighbor Jacobian
238  }
239 }
240 
241 void
243  const MooseVariableFEBase & jvar)
244 {
245  const ArrayVariableTestValue & test_space =
247 
250  else
252 
253  if (_local_ke.n() == 0 || _local_ke.m() == 0)
254  return;
255 
257  {
258  const auto & jv0 = static_cast<const MooseVariable &>(jvar);
259  const VariableTestValue & loc_phi =
261  : jv0.phiFaceNeighbor();
262 
263  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
264  {
266  for (_i = 0; _i < test_space.size(); _i++)
267  for (_j = 0; _j < loc_phi.size(); _j++)
268  {
271  _i,
272  test_space.size(),
273  _j,
274  loc_phi.size(),
275  _var.number(),
276  jvar.number(),
277  v);
278  }
279  }
280  }
282  {
283  const auto & jv1 = static_cast<const ArrayMooseVariable &>(jvar);
284  const ArrayVariableTestValue & loc_phi =
285  (type == Moose::ElementElement || type == Moose::NeighborElement) ? jv1.phiFace()
286  : jv1.phiFaceNeighbor();
287 
288  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
289  {
291  for (_i = 0; _i < test_space.size(); _i++)
292  for (_j = 0; _j < loc_phi.size(); _j++)
293  {
296  _i,
297  test_space.size(),
298  _j,
299  loc_phi.size(),
300  _var.number(),
301  jvar.number(),
302  v);
303  }
304  }
305  }
306  else
307  mooseError("Vector variable cannot be coupled into array DG kernel currently");
308 
310 
312  _var.number() == jvar.number())
313  {
315  Threads::spin_mutex::scoped_lock lock(_jacoby_vars_mutex);
316  for (const auto & var : _diag_save_in)
317  {
318  auto * avar = dynamic_cast<ArrayMooseVariable *>(var);
319  if (!avar)
320  mooseError("Save-in variable for an array kernel must be an array variable");
321 
323  avar->addSolution(diag);
324  else
325  avar->addSolutionNeighbor(diag);
326  }
327  }
328 }
329 
332  const MooseVariableFEBase & jvar)
333 {
334  if (jvar.number() == _var.number())
335  return computeQpJacobian(type).asDiagonal();
336  else
337  return RealEigenMatrix::Zero(_var.count(), jvar.count());
338 }
virtual void computeOffDiagJacobian(unsigned int jvar) override
Override this function to consider couplings of components of the array variable. ...
VarFieldType
Definition: MooseTypes.h:721
const QBase *const & _qrule
Quadrature rule.
Definition: DGKernelBase.h:107
static InputParameters validParams()
Factory constructor initializes all internal references needed for residual computation.
Definition: DGKernelBase.C:26
const libMesh::FEType & feType() const
Get the type of finite element object.
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
virtual void computeOffDiagElemNeighJacobian(Moose::DGJacobianType type, const MooseVariableFEBase &jvar) override
Computes the element-element off-diagonal Jacobian.
Class for stuff related to variables.
Definition: Adaptivity.h:31
std::string incompatVarMsg(MooseVariableFieldBase &var1, MooseVariableFieldBase &var2)
Builds and returns a string of the form:
Definition: MooseError.C:26
virtual void precalculateOffDiagJacobian(unsigned int)
unsigned int number() const
Get variable number coming from libMesh.
void saveDiagLocalArrayJacobian(DenseMatrix< Number > &ke, unsigned int i, unsigned int ntest, unsigned int j, unsigned int nphi, unsigned int ivar, const RealEigenVector &v) const
Helper function for assembling diagonal Jacobian contriubutions on local quadrature points for an arr...
Definition: Assembly.h:1806
std::vector< AuxVariableName > _save_in_strings
Definition: DGKernelBase.h:118
RealEigenVector _work_vector
Work vector for residual computation.
unsigned int count() const
Get the number of components Note: For standard and vector variables, the number is one...
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: DGKernelBase.h:121
virtual void initQpOffDiagJacobian(Moose::DGJacobianType, const MooseVariableFEBase &)
Put necessary evaluations depending on qp but independent on test and shape functions here for off-di...
Definition: ArrayDGKernel.h:93
unsigned int m() const
This class provides an interface for common operations on field variables of both FE and FV types wit...
unsigned int _i
Definition: DGKernelBase.h:136
virtual void initQpResidual(Moose::DGResidualType)
Put necessary evaluations depending on qp but independent on test functions here. ...
Definition: ArrayDGKernel.h:82
DGResidualType
Definition: MooseTypes.h:743
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: DGKernelBase.h:122
void saveLocalArrayResidual(DenseVector< Number > &re, unsigned int i, unsigned int ntest, const RealEigenVector &v) const
Helper function for assembling residual contriubutions on local quadrature points for an array kernel...
Definition: Assembly.h:1787
THREAD_ID _tid
The thread ID for this kernel.
DenseVector< Real > getJacobianDiagonal(DenseMatrix< Number > &ke)
Definition: Assembly.h:1857
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
virtual const std::string & name() const
Get the name of the class.
Definition: MooseBase.h:57
const MooseVariableFieldBase & getVariable(unsigned int jvar_num) const
Retrieve the variable object from our system associated with jvar_num.
Serves as a base class for DGKernel and ADDGKernel.
Definition: DGKernelBase.h:32
virtual void computeElemNeighJacobian(Moose::DGJacobianType type) override
Computes the element/neighbor-element/neighbor Jacobian.
static InputParameters validParams()
Factory constructor initializes all internal references needed for residual computation.
Definition: ArrayDGKernel.C:27
virtual void computeElemNeighResidual(Moose::DGResidualType type) override
Computes the residual for this element or the neighbor.
const ArrayVariableTestValue & _test
test functions
virtual RealEigenMatrix computeQpOffDiagJacobian(Moose::DGJacobianType type, const MooseVariableFEBase &jvar)
This is the virtual that derived classes should override for computing the off-diag Jacobian...
ArrayDGKernel(const InputParameters &parameters)
Definition: ArrayDGKernel.C:33
SystemBase & _sys
Reference to the EquationSystem object.
unsigned int _qp
Definition: DGKernelBase.h:134
void prepareMatrixTagNeighbor(Assembly &assembly, unsigned int ivar, unsigned int jvar, Moose::DGJacobianType type)
Prepare data for computing element jacobian according to the active tags for DG and interface kernels...
std::vector< MooseVariableFEBase * > _save_in
Definition: DGKernelBase.h:117
MooseVariableFE< RealEigenVector > * mooseVariable() const
Return the MooseVariableFE object that this interface acts on.
Enhances MooseVariableInterface interface provide values from neighbor elements.
const unsigned int _count
Number of components of the array variable.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:259
unsigned int _j
Definition: DGKernelBase.h:136
OutputTools< Real >::VariableTestValue VariableTestValue
Definition: MooseTypes.h:324
SubProblem & _subproblem
Reference to this kernel&#39;s SubProblem.
virtual void initQpJacobian(Moose::DGJacobianType)
Put necessary evaluations depending on qp but independent on test and shape functions here...
Definition: ArrayDGKernel.h:87
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:714
std::vector< AuxVariableName > _diag_save_in_strings
Definition: DGKernelBase.h:123
bool excludeBoundary() const
Check current element if it contains broken boundary.
Definition: DGKernelBase.C:188
const std::string & type() const
Get the type of this class.
Definition: MooseBase.h:51
const MooseArray< Real > & _coord
Coordinate transform mainly for curvilinear coordinates.
Definition: DGKernelBase.h:111
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 ...
Eigen::Matrix< Real, Eigen::Dynamic, Eigen::Dynamic > RealEigenMatrix
Definition: MooseTypes.h:149
virtual const MooseVariableFieldBase & getVariable(const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type=Moose::VarKindType::VAR_ANY, Moose::VarFieldType expected_var_field_type=Moose::VarFieldType::VAR_FIELD_ANY) const =0
Returns the variable reference for requested variable which must be of the expected_var_type (Nonline...
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
static Threads::spin_mutex _resid_vars_mutex
Definition: DGKernelBase.h:140
virtual bool hasVariable(const std::string &var_name) const
Query a system for a variable.
Definition: SystemBase.C:834
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
DGJacobianType
Definition: MooseTypes.h:749
forward declarations
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
void prepareVectorTagNeighbor(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual the according to active tags for DG and interface kernels...
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:173
virtual RealEigenVector computeQpJacobian(Moose::DGJacobianType)
This is the virtual that derived classes should override for computing the Jacobian on neighboring el...
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
const MooseArray< Real > & _JxW
Jacobian det times quadrature weighting on quadrature points.
Definition: DGKernelBase.h:109
virtual Moose::VarFieldType fieldType() const =0
Field type of this variable.
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
virtual void computeQpResidual(Moose::DGResidualType type, RealEigenVector &residual)=0
This is the virtual that derived classes should override for computing the residual on neighboring el...
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:179
const ArrayVariablePhiValue & _phi_neighbor
Side shape function.
const FieldVariablePhiValue & phiFace() const override final
Return the variable&#39;s shape functions on an element face.
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealEigenVector
Definition: MooseTypes.h:146
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
unsigned int n() const
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual according to active tags.
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the active tags.
void saveFullLocalArrayJacobian(DenseMatrix< Number > &ke, unsigned int i, unsigned int ntest, unsigned int j, unsigned int nphi, unsigned int ivar, unsigned int jvar, const RealEigenMatrix &v) const
Helper function for assembling full Jacobian contriubutions on local quadrature points for an array k...
Definition: Assembly.h:1831
SystemBase & sys()
Get the system this variable is part of.
const ArrayVariableTestValue & _test_neighbor
Side test function.
static Threads::spin_mutex _jacoby_vars_mutex
Definition: DGKernelBase.h:141
ArrayMooseVariable & _var
Variable this kernel operates on.
Definition: ArrayDGKernel.h:96
const ArrayVariablePhiValue & _phi
Shape functions.
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: DGKernelBase.h:116