www.mooseframework.org
DGKernel.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 "DGKernel.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(_assembly.phiFace(_var)),
42  _grad_phi(_assembly.gradPhiFace(_var)),
43 
44  _test(_var.phiFace()),
45  _grad_test(_var.gradPhiFace()),
46 
47  _phi_neighbor(_assembly.phiFaceNeighbor(_var)),
48  _grad_phi_neighbor(_assembly.gradPhiFaceNeighbor(_var)),
49 
50  _test_neighbor(_var.phiFaceNeighbor()),
51  _grad_test_neighbor(_var.gradPhiFaceNeighbor()),
52 
53  _u_neighbor(_is_implicit ? _var.slnNeighbor() : _var.slnOldNeighbor()),
54  _grad_u_neighbor(_is_implicit ? _var.gradSlnNeighbor() : _var.gradSlnOldNeighbor())
55 {
57 
58  _save_in.resize(_save_in_strings.size());
59  _diag_save_in.resize(_diag_save_in_strings.size());
60 
61  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
62  {
67 
69  mooseError("Trying to use solution variable " + _save_in_strings[i] +
70  " as a save_in variable in " + name());
71 
72  if (var->feType() != _var.feType())
73  paramError(
74  "save_in",
75  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
77 
78  _save_in[i] = var;
81  }
82 
83  _has_save_in = _save_in.size() > 0;
84 
85  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
86  {
91 
93  mooseError("Trying to use solution variable " + _diag_save_in_strings[i] +
94  " as a diag_save_in variable in " + name());
95 
96  if (var->feType() != _var.feType())
97  paramError(
98  "diag_save_in",
99  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
101 
102  _diag_save_in[i] = var;
105  }
106 
107  _has_diag_save_in = _diag_save_in.size() > 0;
108 }
109 
110 void
112 {
113  bool is_elem;
114  if (type == Moose::Element)
115  is_elem = true;
116  else
117  is_elem = false;
118 
119  const VariableTestValue & test_space = is_elem ? _test : _test_neighbor;
120 
121  if (is_elem)
123  else
125 
126  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
127  {
129  for (_i = 0; _i < test_space.size(); _i++)
131  }
132 
134 
135  if (_has_save_in)
136  {
137  Threads::spin_mutex::scoped_lock lock(_resid_vars_mutex);
138  for (const auto & var : _save_in)
139  {
140  const std::vector<dof_id_type> & dof_indices =
141  is_elem ? var->dofIndices() : var->dofIndicesNeighbor();
142  var->sys().solution().add_vector(_local_re, dof_indices);
143  }
144  }
145 }
146 
147 void
149 {
150  const VariableTestValue & test_space =
152  const VariableTestValue & loc_phi =
154 
157  else
159 
160  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
161  {
163  for (_i = 0; _i < test_space.size(); _i++)
164  for (_j = 0; _j < loc_phi.size(); _j++)
166  }
167 
169 
171  {
172  unsigned int rows = _local_ke.m();
173  DenseVector<Number> diag(rows);
174  for (unsigned int i = 0; i < rows; i++)
175  diag(i) = _local_ke(i, i);
176 
177  Threads::spin_mutex::scoped_lock lock(_jacoby_vars_mutex);
178  for (const auto & var : _diag_save_in)
179  {
181  var->sys().solution().add_vector(diag, var->dofIndices());
182  else
183  var->sys().solution().add_vector(diag, var->dofIndicesNeighbor());
184  }
185  }
186 }
187 
188 void
190  const MooseVariableFEBase & jvar)
191 {
192  const VariableTestValue & test_space =
194  const VariableTestValue & loc_phi =
196 
199  else
201 
202  if (_local_ke.n() == 0 || _local_ke.m() == 0)
203  return;
204 
205  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
206  {
208  for (_i = 0; _i < test_space.size(); _i++)
209  for (_j = 0; _j < loc_phi.size(); _j++)
210  _local_ke(_i, _j) +=
212  }
213 
215 }
216 
217 Real
219 {
220  return 0.;
221 }
VarFieldType
Definition: MooseTypes.h:634
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
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
std::string incompatVarMsg(MooseVariableFieldBase &var1, MooseVariableFieldBase &var2)
Builds and returns a string of the form:
Definition: MooseError.C:23
unsigned int number() const
Get variable number coming from libMesh.
std::vector< AuxVariableName > _save_in_strings
Definition: DGKernelBase.h:118
virtual void computeElemNeighResidual(Moose::DGResidualType type) override
Computes the residual for this element or the neighbor.
Definition: DGKernel.C:111
const VariablePhiValue & _phi_neighbor
Side shape function.
Definition: DGKernel.h:106
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual Real computeQpJacobian(Moose::DGJacobianType type)=0
This is the virtual that derived classes should override for computing the Jacobian on neighboring el...
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: DGKernelBase.h:121
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
DGResidualType
Definition: MooseTypes.h:656
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: DGKernelBase.h:122
THREAD_ID _tid
The thread ID for this kernel.
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:56
const FEType & feType() const
Get the type of finite element object.
Serves as a base class for DGKernel and ADDGKernel.
Definition: DGKernelBase.h:32
MooseVariable & _var
Variable this kernel operates on.
Definition: DGKernel.h:92
SystemBase & _sys
Reference to the EquationSystem object.
static InputParameters validParams()
Factory constructor initializes all internal references needed for residual computation.
Definition: DGKernel.C:27
const VariableTestValue & _test_neighbor
Side test function.
Definition: DGKernel.h:110
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< Real > * mooseVariable() const
Enhances MooseVariableInterface interface provide values from neighbor elements.
virtual Real computeQpOffDiagJacobian(Moose::DGJacobianType type, unsigned int jvar)
This is the virtual that derived classes should override for computing the off-diag Jacobian...
Definition: DGKernel.C:218
unsigned int _j
Definition: DGKernelBase.h:136
OutputTools< Real >::VariableTestValue VariableTestValue
Definition: MooseTypes.h:312
SubProblem & _subproblem
Reference to this kernel&#39;s SubProblem.
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:627
std::vector< AuxVariableName > _diag_save_in_strings
Definition: DGKernelBase.h:123
virtual void precalculateQpOffDiagJacobian(Moose::DGJacobianType, const MooseVariableFEBase &)
Insertion point for evaluations that depend on qp but are independent of the test and shape functions...
Definition: DGKernel.h:86
DGKernel(const InputParameters &parameters)
Definition: DGKernel.C:33
const std::string & type() const
Get the type of this class.
Definition: MooseBase.h:50
const MooseArray< Real > & _coord
Coordinate transform mainly for curvilinear coordinates.
Definition: DGKernelBase.h:111
virtual void precalculateQpJacobian(Moose::DGJacobianType)
Insertion point for evaluations that depend on qp but are independent of the test and shape functions...
Definition: DGKernel.h:80
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 ...
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:800
virtual void computeElemNeighJacobian(Moose::DGJacobianType type) override
Computes the element/neighbor-element/neighbor Jacobian.
Definition: DGKernel.C:148
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
const VariableTestValue & _test
test functions
Definition: DGKernel.h:102
DGJacobianType
Definition: MooseTypes.h:662
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...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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:163
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
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
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:169
const VariablePhiValue & _phi
Shape functions.
Definition: DGKernel.h:98
virtual void precalculateQpResidual(Moose::DGResidualType)
Insertion point for evaluations that depend on qp but are independent of the test functions...
Definition: DGKernel.h:74
virtual Real computeQpResidual(Moose::DGResidualType type)=0
This is the virtual that derived classes should override for computing the residual on neighboring el...
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.
SystemBase & sys()
Get the system this variable is part of.
virtual void computeOffDiagElemNeighJacobian(Moose::DGJacobianType type, const MooseVariableFEBase &jvar) override
Computes the element-element off-diagonal Jacobian.
Definition: DGKernel.C:189
static Threads::spin_mutex _jacoby_vars_mutex
Definition: DGKernelBase.h:141
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: DGKernelBase.h:116