Line data Source code
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 "DGKernelBase.h"
11 : #include "MooseVariable.h"
12 : #include "Problem.h"
13 : #include "SubProblem.h"
14 : #include "SystemBase.h"
15 : #include "MaterialData.h"
16 : #include "ParallelUniqueId.h"
17 :
18 : #include "libmesh/dof_map.h"
19 : #include "libmesh/dense_vector.h"
20 : #include "libmesh/numeric_vector.h"
21 : #include "libmesh/dense_subvector.h"
22 : #include "libmesh/libmesh_common.h"
23 : #include "libmesh/quadrature.h"
24 :
25 : InputParameters
26 202459 : DGKernelBase::validParams()
27 : {
28 202459 : InputParameters params = NeighborResidualObject::validParams();
29 202459 : params += TwoMaterialPropertyInterface::validParams();
30 202459 : params += BlockRestrictable::validParams();
31 202459 : params += BoundaryRestrictable::validParams();
32 607377 : params.addParam<bool>("use_displaced_mesh",
33 404918 : false,
34 : "Whether or not this object should use the "
35 : "displaced mesh for computation. Note that in "
36 : "the case this is true but no displacements "
37 : "are provided in the Mesh block the "
38 : "undisplaced mesh will still be used.");
39 202459 : params.addParamNamesToGroup("use_displaced_mesh", "Advanced");
40 :
41 202459 : params.addParam<std::vector<AuxVariableName>>(
42 : "save_in",
43 : {},
44 : "The name of auxiliary variables to save this Kernel's residual contributions to. "
45 : " Everything about that variable must match everything about this variable (the "
46 : "type, what blocks it's on, etc.)");
47 202459 : params.addParam<std::vector<AuxVariableName>>(
48 : "diag_save_in",
49 : {},
50 : "The name of auxiliary variables to save this Kernel's diagonal Jacobian "
51 : "contributions to. Everything about that variable must match everything "
52 : "about this variable (the type, what blocks it's on, etc.)");
53 202459 : params.addParamNamesToGroup("diag_save_in save_in", "Advanced");
54 :
55 : // DG Kernels always need one layer of ghosting.
56 202459 : params.addRelationshipManager("ElementSideNeighborLayers",
57 : Moose::RelationshipManagerType::GEOMETRIC |
58 : Moose::RelationshipManagerType::ALGEBRAIC |
59 : Moose::RelationshipManagerType::COUPLING);
60 :
61 202459 : params.addParam<std::vector<BoundaryName>>(
62 : "exclude_boundary", "The internal side sets to be excluded from this kernel.");
63 202459 : params.registerBase("DGKernel");
64 202459 : params.registerSystemAttributeName("DGKernel");
65 :
66 202459 : return params;
67 0 : }
68 :
69 : // Static mutex definitions
70 : Threads::spin_mutex DGKernelBase::_resid_vars_mutex;
71 : Threads::spin_mutex DGKernelBase::_jacoby_vars_mutex;
72 :
73 1433 : DGKernelBase::DGKernelBase(const InputParameters & parameters)
74 : : NeighborResidualObject(parameters),
75 : BlockRestrictable(this),
76 : BoundaryRestrictable(this, false), // false for _not_ nodal
77 : NeighborCoupleableMooseVariableDependencyIntermediateInterface(this, false, false),
78 : TwoMaterialPropertyInterface(this, blockIDs(), boundaryIDs()),
79 : ElementIDInterface(this),
80 2866 : _current_elem(_assembly.elem()),
81 1433 : _current_elem_volume(_assembly.elemVolume()),
82 :
83 1433 : _neighbor_elem(_assembly.neighbor()),
84 1433 : _neighbor_elem_volume(_assembly.neighborVolume()),
85 :
86 1433 : _current_side(_assembly.side()),
87 1433 : _current_side_elem(_assembly.sideElem()),
88 1433 : _current_side_volume(_assembly.sideElemVolume()),
89 :
90 1433 : _coord_sys(_assembly.coordSystem()),
91 1433 : _q_point(_assembly.qPointsFace()),
92 1433 : _qrule(_assembly.qRuleFace()),
93 1433 : _JxW(_assembly.JxWFace()),
94 1433 : _coord(_assembly.coordTransformation()),
95 1433 : _normals(_assembly.normals()),
96 :
97 1433 : _save_in_strings(parameters.get<std::vector<AuxVariableName>>("save_in")),
98 4299 : _diag_save_in_strings(parameters.get<std::vector<AuxVariableName>>("diag_save_in"))
99 : {
100 : // Gather information on broken boundaries
101 2866 : std::vector<BoundaryName> bnd = isParamValid("exclude_boundary")
102 1433 : ? getParam<std::vector<BoundaryName>>("exclude_boundary")
103 1433 : : std::vector<BoundaryName>(0);
104 1433 : auto bnd_ids = _mesh.getBoundaryIDs(bnd);
105 :
106 : // check if the broken boundary ids are valid
107 1433 : auto & valid_ids = _mesh.meshSidesetIds();
108 1433 : std::vector<BoundaryName> diff;
109 1433 : for (unsigned int i = 0; i < bnd_ids.size(); ++i)
110 0 : if (valid_ids.find(bnd_ids[i]) == valid_ids.end())
111 0 : diff.push_back(bnd[i]);
112 1433 : if (!diff.empty())
113 : {
114 0 : auto msg = "DGKernel '" + name() +
115 : "' contains the following boundary names that do not exist on the mesh: " +
116 0 : Moose::stringify(diff, ", ");
117 0 : paramError("exclude_boundary", msg);
118 0 : }
119 :
120 1433 : _excluded_boundaries.insert(bnd_ids.begin(), bnd_ids.end());
121 1433 : }
122 :
123 : void
124 2545108 : DGKernelBase::computeResidual()
125 : {
126 2545108 : if (!excludeBoundary())
127 : {
128 2545108 : precalculateResidual();
129 :
130 : // Compute the residual for this element
131 2545108 : computeElemNeighResidual(Moose::Element);
132 :
133 : // Compute the residual for the neighbor
134 2545108 : computeElemNeighResidual(Moose::Neighbor);
135 : }
136 2545108 : }
137 :
138 : void
139 113132 : DGKernelBase::computeJacobian()
140 : {
141 113132 : if (!excludeBoundary())
142 : {
143 113132 : precalculateJacobian();
144 :
145 : // Compute element-element Jacobian
146 113132 : computeElemNeighJacobian(Moose::ElementElement);
147 :
148 : // Compute element-neighbor Jacobian
149 113132 : computeElemNeighJacobian(Moose::ElementNeighbor);
150 :
151 : // Compute neighbor-element Jacobian
152 113132 : computeElemNeighJacobian(Moose::NeighborElement);
153 :
154 : // Compute neighbor-neighbor Jacobian
155 113132 : computeElemNeighJacobian(Moose::NeighborNeighbor);
156 : }
157 113132 : }
158 :
159 : void
160 18262 : DGKernelBase::computeOffDiagJacobian(const unsigned int jvar_num)
161 : {
162 18262 : if (!excludeBoundary())
163 : {
164 18262 : const auto & jvar = getVariable(jvar_num);
165 :
166 18262 : if (jvar_num == variable().number())
167 15260 : computeJacobian();
168 : else
169 : {
170 3002 : precalculateOffDiagJacobian(jvar_num);
171 :
172 : // Compute element-element Jacobian
173 3002 : computeOffDiagElemNeighJacobian(Moose::ElementElement, jvar);
174 :
175 : // Compute element-neighbor Jacobian
176 3002 : computeOffDiagElemNeighJacobian(Moose::ElementNeighbor, jvar);
177 :
178 : // Compute neighbor-element Jacobian
179 3002 : computeOffDiagElemNeighJacobian(Moose::NeighborElement, jvar);
180 :
181 : // Compute neighbor-neighbor Jacobian
182 3002 : computeOffDiagElemNeighJacobian(Moose::NeighborNeighbor, jvar);
183 : }
184 : }
185 18262 : }
186 :
187 : bool
188 2723040 : DGKernelBase::excludeBoundary() const
189 : {
190 2723040 : if (_excluded_boundaries.empty())
191 2723040 : return false;
192 :
193 0 : auto boundary_ids = _mesh.getBoundaryIDs(_current_elem, _current_side);
194 0 : for (auto bid : boundary_ids)
195 0 : if (_excluded_boundaries.find(bid) != _excluded_boundaries.end())
196 0 : return true;
197 :
198 : // make sure we will also break on the neighboring side
199 0 : unsigned int neighbor_side = _neighbor_elem->which_neighbor_am_i(_current_elem);
200 0 : boundary_ids = _mesh.getBoundaryIDs(_neighbor_elem, neighbor_side);
201 0 : for (auto bid : boundary_ids)
202 0 : if (_excluded_boundaries.find(bid) != _excluded_boundaries.end())
203 0 : return true;
204 :
205 0 : return false;
206 0 : }
207 :
208 : void
209 154780 : DGKernelBase::prepareShapes(const unsigned int var_num)
210 : {
211 154780 : _subproblem.prepareFaceShapes(var_num, _tid);
212 154780 : }
|