www.mooseframework.org
DGKernel.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #include "DGKernel.h"
16 #include "Assembly.h"
17 #include "MooseVariable.h"
18 #include "Problem.h"
19 #include "SubProblem.h"
20 #include "SystemBase.h"
21 #include "MaterialData.h"
22 #include "ParallelUniqueId.h"
23 
24 #include "libmesh/dof_map.h"
25 #include "libmesh/dense_vector.h"
26 #include "libmesh/numeric_vector.h"
27 #include "libmesh/dense_subvector.h"
28 #include "libmesh/libmesh_common.h"
29 #include "libmesh/quadrature.h"
30 
31 template <>
34 {
41  params.addRequiredParam<NonlinearVariableName>(
42  "variable", "The name of the variable that this boundary condition applies to");
43  params.addParam<bool>("use_displaced_mesh",
44  false,
45  "Whether or not this object should use the "
46  "displaced mesh for computation. Note that in "
47  "the case this is true but no displacements "
48  "are provided in the Mesh block the "
49  "undisplaced mesh will still be used.");
50  params.addParamNamesToGroup("use_displaced_mesh", "Advanced");
51 
52  params.declareControllable("enable");
53  params.registerBase("DGKernel");
54  params.addParam<std::vector<AuxVariableName>>(
55  "save_in",
56  "The name of auxiliary variables to save this Kernel's residual contributions to. "
57  " Everything about that variable must match everything about this variable (the "
58  "type, what blocks it's on, etc.)");
59  params.addParam<std::vector<AuxVariableName>>(
60  "diag_save_in",
61  "The name of auxiliary variables to save this Kernel's diagonal Jacobian "
62  "contributions to. Everything about that variable must match everything "
63  "about this variable (the type, what blocks it's on, etc.)");
64 
65  return params;
66 }
67 
68 // Static mutex definitions
69 Threads::spin_mutex DGKernel::_resid_vars_mutex;
70 Threads::spin_mutex DGKernel::_jacoby_vars_mutex;
71 
73  : MooseObject(parameters),
74  BlockRestrictable(this),
75  BoundaryRestrictable(this, false), // false for _not_ nodal
76  SetupInterface(this),
77  TransientInterface(this),
78  FunctionInterface(this),
79  UserObjectInterface(this),
81  TwoMaterialPropertyInterface(this, blockIDs(), boundaryIDs()),
82  Restartable(parameters, "DGKernels"),
83  ZeroInterface(parameters),
84  MeshChangedInterface(parameters),
85  _subproblem(*parameters.get<SubProblem *>("_subproblem")),
86  _sys(*parameters.get<SystemBase *>("_sys")),
87  _tid(parameters.get<THREAD_ID>("_tid")),
88  _assembly(_subproblem.assembly(_tid)),
89  _var(_sys.getVariable(_tid, parameters.get<NonlinearVariableName>("variable"))),
90  _mesh(_subproblem.mesh()),
91 
92  _current_elem(_assembly.elem()),
93  _current_elem_volume(_assembly.elemVolume()),
94 
95  _neighbor_elem(_assembly.neighbor()),
96 
97  _current_side(_assembly.side()),
98  _current_side_elem(_assembly.sideElem()),
99  _current_side_volume(_assembly.sideElemVolume()),
100 
101  _coord_sys(_assembly.coordSystem()),
102  _q_point(_assembly.qPointsFace()),
103  _qrule(_assembly.qRuleFace()),
104  _JxW(_assembly.JxWFace()),
105  _coord(_assembly.coordTransformation()),
106 
107  _u(_is_implicit ? _var.sln() : _var.slnOld()),
108  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld()),
109 
110  _phi(_assembly.phiFace()),
111  _grad_phi(_assembly.gradPhiFace()),
112 
113  _test(_var.phiFace()),
114  _grad_test(_var.gradPhiFace()),
115 
116  _normals(_var.normals()),
117 
118  _phi_neighbor(_assembly.phiFaceNeighbor()),
119  _grad_phi_neighbor(_assembly.gradPhiFaceNeighbor()),
120 
121  _test_neighbor(_var.phiFaceNeighbor()),
122  _grad_test_neighbor(_var.gradPhiFaceNeighbor()),
123 
124  _u_neighbor(_is_implicit ? _var.slnNeighbor() : _var.slnOldNeighbor()),
125  _grad_u_neighbor(_is_implicit ? _var.gradSlnNeighbor() : _var.gradSlnOldNeighbor()),
126 
127  _save_in_strings(parameters.get<std::vector<AuxVariableName>>("save_in")),
128  _diag_save_in_strings(parameters.get<std::vector<AuxVariableName>>("diag_save_in"))
129 {
130  _save_in.resize(_save_in_strings.size());
131  _diag_save_in.resize(_diag_save_in_strings.size());
132 
133  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
134  {
136 
138  mooseError("Trying to use solution variable " + _save_in_strings[i] +
139  " as a save_in variable in " + name());
140 
141  if (var->feType() != _var.feType())
142  mooseError("Error in " + name() + ". When saving residual values in an Auxiliary variable "
143  "the AuxVariable must be the same type as the nonlinear "
144  "variable the object is acting on.");
145 
146  _save_in[i] = var;
149  }
150 
151  _has_save_in = _save_in.size() > 0;
152 
153  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
154  {
156 
158  mooseError("Trying to use solution variable " + _diag_save_in_strings[i] +
159  " as a diag_save_in variable in " + name());
160 
161  if (var->feType() != _var.feType())
162  mooseError("Error in " + name() + ". When saving diagonal Jacobian values in an Auxiliary "
163  "variable the AuxVariable must be the same type as the "
164  "nonlinear variable the object is acting on.");
165 
166  _diag_save_in[i] = var;
169  }
170 
171  _has_diag_save_in = _diag_save_in.size() > 0;
172 }
173 
175 
176 void
178 {
179  bool is_elem;
180  if (type == Moose::Element)
181  is_elem = true;
182  else
183  is_elem = false;
184 
185  const VariableTestValue & test_space = is_elem ? _test : _test_neighbor;
186  DenseVector<Number> & re = is_elem ? _assembly.residualBlock(_var.number())
188  _local_re.resize(re.size());
189  _local_re.zero();
190 
191  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
192  for (_i = 0; _i < test_space.size(); _i++)
193  _local_re(_i) += _JxW[_qp] * _coord[_qp] * computeQpResidual(type);
194 
195  re += _local_re;
196 
197  if (_has_save_in)
198  {
199  Threads::spin_mutex::scoped_lock lock(_resid_vars_mutex);
200  for (const auto & var : _save_in)
201  {
202  std::vector<dof_id_type> & dof_indices =
203  is_elem ? var->dofIndices() : var->dofIndicesNeighbor();
204  var->sys().solution().add_vector(_local_re, dof_indices);
205  }
206  }
207 }
208 
209 void
211 {
212  // Compute the residual for this element
214 
215  // Compute the residual for the neighbor
217 }
218 
219 void
221 {
222  const VariableTestValue & test_space =
224  const VariableTestValue & loc_phi =
226  DenseMatrix<Number> & Kxx =
227  type == Moose::ElementElement
229  : type == Moose::ElementNeighbor
232  : type == Moose::NeighborElement
237  _local_kxx.resize(Kxx.m(), Kxx.n());
238  _local_kxx.zero();
239 
240  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
241  for (_i = 0; _i < test_space.size(); _i++)
242  for (_j = 0; _j < loc_phi.size(); _j++)
243  _local_kxx(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpJacobian(type);
244 
245  Kxx += _local_kxx;
246 
248  {
249  unsigned int rows = _local_kxx.m();
250  DenseVector<Number> diag(rows);
251  for (unsigned int i = 0; i < rows; i++)
252  diag(i) = _local_kxx(i, i);
253 
254  Threads::spin_mutex::scoped_lock lock(_jacoby_vars_mutex);
255  for (const auto & var : _diag_save_in)
256  {
257  if (type == Moose::ElementElement)
258  var->sys().solution().add_vector(diag, var->dofIndices());
259  else
260  var->sys().solution().add_vector(diag, var->dofIndicesNeighbor());
261  }
262  }
263 }
264 void
266 {
267  // Compute element-element Jacobian
269 
270  // Compute element-neighbor Jacobian
272 
273  // Compute neighbor-element Jacobian
275 
276  // Compute neighbor-neighbor Jacobian
278 }
279 
280 void
282 {
283  const VariableTestValue & test_space =
285  const VariableTestValue & loc_phi =
287  DenseMatrix<Number> & Kxx =
288  type == Moose::ElementElement
290  : type == Moose::ElementNeighbor
292  : type == Moose::NeighborElement
296 
297  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
298  for (_i = 0; _i < test_space.size(); _i++)
299  for (_j = 0; _j < loc_phi.size(); _j++)
300  Kxx(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(type, jvar);
301 }
302 
303 void
305 {
306  if (jvar == _var.number())
307  computeJacobian();
308  else
309  {
310  // Compute element-element Jacobian
312 
313  // Compute element-neighbor Jacobian
315 
316  // Compute neighbor-element Jacobian
318 
319  // Compute neighbor-neighbor Jacobian
321  }
322 }
323 
324 Real
326 {
327  return 0.;
328 }
329 
332 {
333  return _var;
334 }
335 
336 SubProblem &
338 {
339  return _subproblem;
340 }
341 
342 const Real &
344 {
345  return _assembly.neighborVolume();
346 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
std::vector< AuxVariableName > _diag_save_in_strings
Definition: DGKernel.h:197
virtual void computeOffDiagElemNeighJacobian(Moose::DGJacobianType type, unsigned int jvar)
Computes the element-element off-diagonal Jacobian.
Definition: DGKernel.C:281
A class for creating restricted objects.
Definition: Restartable.h:31
THREAD_ID _tid
Definition: DGKernel.h:116
DenseMatrix< Number > _local_kxx
Holds residual entries as they are accumulated by this DGKernel.
Definition: DGKernel.h:187
const FEType & feType() const
Get the type of finite element object.
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
Class for stuff related to variables.
Definition: MooseVariable.h:43
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: DGKernel.h:190
InputParameters validParams< BlockRestrictable >()
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name)=0
Returns the variable reference for requested variable which may be in any system. ...
SubProblem & subProblem()
Return a reference to the subproblem.
Definition: DGKernel.C:337
std::vector< AuxVariableName > _save_in_strings
Definition: DGKernel.h:192
Assembly & _assembly
Definition: DGKernel.h:118
const MooseArray< Real > & _coord
Definition: DGKernel.h:145
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: DGKernel.h:195
const VariablePhiValue & _phi_neighbor
Side shape function.
Definition: DGKernel.h:169
virtual void computeResidual()
Computes the residual for the current side.
Definition: DGKernel.C:210
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...
void addMooseVariableDependency(MooseVariable *var)
Call this function to add the passed in MooseVariable as a variable that this object depends on...
InputParameters validParams< DGKernel >()
Definition: DGKernel.C:33
/class BoundaryRestrictable /brief Provides functionality for limiting the object to certain boundary...
const Real & neighborVolume()
Returns the reference to the current neighbor volume.
Definition: Assembly.C:353
DGResidualType
Definition: MooseTypes.h:184
Base class for a system (of equations)
Definition: SystemBase.h:91
DenseVector< Number > & residualBlock(unsigned int var_num, Moose::KernelType type=Moose::KT_NONTIME)
Definition: Assembly.h:504
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
MooseVariable & _var
Definition: DGKernel.h:119
const VariableTestValue & _test_neighbor
Side test function.
Definition: DGKernel.h:174
unsigned int _j
Definition: DGKernel.h:147
Interface for objects that needs transient capabilities.
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:325
Interface for notifications that the mesh has changed.
Every object that can be built by the factory should be derived from this class.
Definition: MooseObject.h:36
static Threads::spin_mutex _resid_vars_mutex
Definition: DGKernel.h:225
virtual void computeJacobian()
Computes the jacobian for the current side.
Definition: DGKernel.C:265
InputParameters validParams< MeshChangedInterface >()
DenseVector< Number > _local_re
Holds residual entries as they are accumulated by this DGKernel.
Definition: DGKernel.h:184
DGKernel(const InputParameters &parameters)
Factory constructor initializes all internal references needed for residual computation.
Definition: DGKernel.C:72
Interface for objects that need to use UserObjects.
virtual bool hasVariable(const std::string &var_name)
Query a system for a variable.
Definition: SystemBase.C:556
InputParameters validParams< MooseObject >()
Definition: MooseObject.C:22
SubProblem & _subproblem
Definition: DGKernel.h:113
const VariableTestValue & _test
Side shape function.
Definition: DGKernel.h:162
DGJacobianType
Definition: MooseTypes.h:190
const MooseArray< Real > & _JxW
Definition: DGKernel.h:144
MatType type
unsigned int number() const
Get variable number coming from libMesh.
unsigned int _i
Definition: DGKernel.h:147
std::vector< MooseVariable * > _diag_save_in
Definition: DGKernel.h:196
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:53
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
Interface to bring zero values inside objects.
Definition: ZeroInterface.h:35
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:150
InputParameters validParams< BoundaryRestrictable >()
virtual void computeOffDiagJacobian(unsigned int jvar)
Computes d-residual / d-jvar...
Definition: DGKernel.C:304
An interface that restricts an object to subdomains via the &#39;blocks&#39; input parameter.
InputParameters validParams< TwoMaterialPropertyInterface >()
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:156
SystemBase & _sys
Definition: DGKernel.h:114
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
virtual void computeElemNeighJacobian(Moose::DGJacobianType type)
Computes the element/neighbor-element/neighbor Jacobian.
Definition: DGKernel.C:220
const VariablePhiValue & _phi
Definition: DGKernel.h:157
virtual void computeElemNeighResidual(Moose::DGResidualType type)
Computes the residual for this element or the neighbor.
Definition: DGKernel.C:177
std::vector< MooseVariable * > _save_in
Definition: DGKernel.h:191
virtual Real computeQpResidual(Moose::DGResidualType type)=0
This is the virtual that derived classes should override for computing the residual on neighboring el...
const Real & getNeighborElemVolume()
The volume (or length) of the current neighbor.
Definition: DGKernel.C:343
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:887
SystemBase & sys()
Get the system this variable is part of.
MooseVariable & variable()
The variable number that this kernel operates on.
Definition: DGKernel.C:331
DenseVector< Number > & residualBlockNeighbor(unsigned int var_num, Moose::KernelType type=Moose::KT_NONTIME)
Definition: Assembly.h:509
DenseMatrix< Number > & jacobianBlockNeighbor(Moose::DGJacobianType type, unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:901
Interface for objects that need to use functions.
virtual ~DGKernel()
Definition: DGKernel.C:174
QBase *& _qrule
Definition: DGKernel.h:143
unsigned int _qp
Definition: DGKernel.h:141
InputParameters validParams< TransientInterface >()
unsigned int THREAD_ID
Definition: MooseTypes.h:79
static Threads::spin_mutex _jacoby_vars_mutex
Definition: DGKernel.h:226