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 
26 template <>
29 {
31  return params;
32 }
33 
34 DGKernel::DGKernel(const InputParameters & parameters) : DGKernelBase(parameters) {}
35 
37 
38 void
40 {
41  bool is_elem;
42  if (type == Moose::Element)
43  is_elem = true;
44  else
45  is_elem = false;
46 
47  const VariableTestValue & test_space = is_elem ? _test : _test_neighbor;
48 
49  if (is_elem)
51  else
53 
54  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
55  for (_i = 0; _i < test_space.size(); _i++)
57 
59 
60  if (_has_save_in)
61  {
62  Threads::spin_mutex::scoped_lock lock(_resid_vars_mutex);
63  for (const auto & var : _save_in)
64  {
65  const std::vector<dof_id_type> & dof_indices =
66  is_elem ? var->dofIndices() : var->dofIndicesNeighbor();
67  var->sys().solution().add_vector(_local_re, dof_indices);
68  }
69  }
70 }
71 
72 void
74 {
75  const VariableTestValue & test_space =
77  const VariableTestValue & loc_phi =
79 
82  else
84 
85  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
86  for (_i = 0; _i < test_space.size(); _i++)
87  for (_j = 0; _j < loc_phi.size(); _j++)
89 
91 
93  {
94  unsigned int rows = _local_ke.m();
95  DenseVector<Number> diag(rows);
96  for (unsigned int i = 0; i < rows; i++)
97  diag(i) = _local_ke(i, i);
98 
99  Threads::spin_mutex::scoped_lock lock(_jacoby_vars_mutex);
100  for (const auto & var : _diag_save_in)
101  {
103  var->sys().solution().add_vector(diag, var->dofIndices());
104  else
105  var->sys().solution().add_vector(diag, var->dofIndicesNeighbor());
106  }
107  }
108 }
109 
110 void
112 {
113  const VariableTestValue & test_space =
115  const VariableTestValue & loc_phi =
117 
120  else
122 
123  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
124  for (_i = 0; _i < test_space.size(); _i++)
125  for (_j = 0; _j < loc_phi.size(); _j++)
127 
129 }
130 
131 Real
133 {
134  return 0.;
135 }
const QBase *const & _qrule
Definition: DGKernelBase.h:140
const VariableTestValue & _test_neighbor
Side test function.
Definition: DGKernelBase.h:171
const VariablePhiValue & _phi_neighbor
Side shape function.
Definition: DGKernelBase.h:166
InputParameters validParams< DGKernelBase >()
Definition: DGKernelBase.C:28
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
unsigned int number() const
Get variable number coming from libMesh.
virtual void computeElemNeighResidual(Moose::DGResidualType type) override
Computes the residual for this element or the neighbor.
Definition: DGKernel.C:39
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:186
InputParameters validParams< DGKernel >()
Definition: DGKernel.C:28
unsigned int _i
Definition: DGKernelBase.h:144
DGResidualType
Definition: MooseTypes.h:509
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: DGKernelBase.h:187
DenseMatrix< Number > _local_ke
Holds residual entries as they are accumulated by this Kernel.
Serves as a base class for DGKernel and ADDGKernel.
Definition: DGKernelBase.h:46
const std::string & type() const
Get the type of this object.
Definition: MooseObject.h:53
unsigned int _qp
Definition: DGKernelBase.h:138
void prepareMatrixTagNeighbor(Assembly &assembly, unsigned int ivar, unsigned int jvar, Moose::DGJacobianType type)
Prepare data for computing element jacobian according to the ative tags for DG and interface kernels...
std::vector< MooseVariableFEBase * > _save_in
Definition: DGKernelBase.h:182
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:132
unsigned int _j
Definition: DGKernelBase.h:144
OutputTools< Real >::VariableTestValue VariableTestValue
Definition: MooseTypes.h:205
DGKernel(const InputParameters &parameters)
Factory constructor initializes all internal references needed for residual computation.
Definition: DGKernel.C:34
const MooseArray< Real > & _coord
Definition: DGKernelBase.h:142
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:199
virtual void computeElemNeighJacobian(Moose::DGJacobianType type) override
Computes the element/neighbor-element/neighbor Jacobian.
Definition: DGKernel.C:73
Assembly & _assembly
Definition: DGKernelBase.h:116
DGJacobianType
Definition: MooseTypes.h:515
MatType type
void prepareVectorTagNeighbor(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual the according to active tags for DG and interface kernels...
const VariablePhiValue & _phi
Definition: DGKernelBase.h:154
const VariableTestValue & _test
Side shape function.
Definition: DGKernelBase.h:159
DenseVector< Number > _local_re
Holds residual entries as they are accumulated by this Kernel.
const MooseArray< Real > & _JxW
Definition: DGKernelBase.h:141
MooseVariable & _var
Definition: DGKernelBase.h:117
virtual Real computeQpResidual(Moose::DGResidualType type)=0
This is the virtual that derived classes should override for computing the residual on neighboring el...
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual the according to active tags.
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the ative tags.
virtual void computeOffDiagElemNeighJacobian(Moose::DGJacobianType type, unsigned int jvar) override
Computes the element-element off-diagonal Jacobian.
Definition: DGKernel.C:111
virtual ~DGKernel()
Definition: DGKernel.C:36
static Threads::spin_mutex _jacoby_vars_mutex
Definition: DGKernelBase.h:200
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: DGKernelBase.h:181