www.mooseframework.org
DGKernelBase.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 "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 
27 {
32  params.addParam<bool>("use_displaced_mesh",
33  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  params.addParamNamesToGroup("use_displaced_mesh", "Advanced");
40 
41  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  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  params.addParamNamesToGroup("diag_save_in save_in", "Advanced");
54 
55  // DG Kernels always need one layer of ghosting.
56  params.addRelationshipManager("ElementSideNeighborLayers",
60 
61  params.addParam<std::vector<BoundaryName>>(
62  "exclude_boundary", "The internal side sets to be excluded from this kernel.");
63  params.registerBase("DGKernel");
64 
65  return params;
66 }
67 
68 // Static mutex definitions
69 Threads::spin_mutex DGKernelBase::_resid_vars_mutex;
70 Threads::spin_mutex DGKernelBase::_jacoby_vars_mutex;
71 
73  : NeighborResidualObject(parameters),
74  BlockRestrictable(this),
75  BoundaryRestrictable(this, false), // false for _not_ nodal
77  TwoMaterialPropertyInterface(this, blockIDs(), boundaryIDs()),
78  ElementIDInterface(this),
79  _current_elem(_assembly.elem()),
80  _current_elem_volume(_assembly.elemVolume()),
81 
82  _neighbor_elem(_assembly.neighbor()),
83  _neighbor_elem_volume(_assembly.neighborVolume()),
84 
85  _current_side(_assembly.side()),
86  _current_side_elem(_assembly.sideElem()),
87  _current_side_volume(_assembly.sideElemVolume()),
88 
89  _coord_sys(_assembly.coordSystem()),
90  _q_point(_assembly.qPointsFace()),
91  _qrule(_assembly.qRuleFace()),
92  _JxW(_assembly.JxWFace()),
93  _coord(_assembly.coordTransformation()),
94  _normals(_assembly.normals()),
95 
96  _save_in_strings(parameters.get<std::vector<AuxVariableName>>("save_in")),
97  _diag_save_in_strings(parameters.get<std::vector<AuxVariableName>>("diag_save_in"))
98 {
99  // Gather information on broken boundaries
100  std::vector<BoundaryName> bnd = isParamValid("exclude_boundary")
101  ? getParam<std::vector<BoundaryName>>("exclude_boundary")
102  : std::vector<BoundaryName>(0);
103  auto bnd_ids = _mesh.getBoundaryIDs(bnd);
104 
105  // check if the broken boundary ids are valid
106  auto & valid_ids = _mesh.meshSidesetIds();
107  std::vector<BoundaryName> diff;
108  for (unsigned int i = 0; i < bnd_ids.size(); ++i)
109  if (valid_ids.find(bnd_ids[i]) == valid_ids.end())
110  diff.push_back(bnd[i]);
111  if (!diff.empty())
112  {
113  auto msg = "DGKernel '" + name() +
114  "' contains the following boundary names that do not exist on the mesh: " +
115  Moose::stringify(diff, ", ");
116  paramError("exclude_boundary", msg);
117  }
118 
119  _excluded_boundaries.insert(bnd_ids.begin(), bnd_ids.end());
120 }
121 
122 void
124 {
125  if (!excludeBoundary())
126  {
128 
129  // Compute the residual for this element
131 
132  // Compute the residual for the neighbor
134  }
135 }
136 
137 void
139 {
140  if (!excludeBoundary())
141  {
143 
144  // Compute element-element Jacobian
146 
147  // Compute element-neighbor Jacobian
149 
150  // Compute neighbor-element Jacobian
152 
153  // Compute neighbor-neighbor Jacobian
155  }
156 }
157 
158 void
159 DGKernelBase::computeOffDiagJacobian(const unsigned int jvar_num)
160 {
161  if (!excludeBoundary())
162  {
163  const auto & jvar = getVariable(jvar_num);
164 
165  if (jvar_num == variable().number())
166  computeJacobian();
167  else
168  {
169  precalculateOffDiagJacobian(jvar_num);
170 
171  // Compute element-element Jacobian
173 
174  // Compute element-neighbor Jacobian
176 
177  // Compute neighbor-element Jacobian
179 
180  // Compute neighbor-neighbor Jacobian
182  }
183  }
184 }
185 
186 bool
188 {
189  if (_excluded_boundaries.empty())
190  return false;
191 
192  auto boundary_ids = _mesh.getBoundaryIDs(_current_elem, _current_side);
193  for (auto bid : boundary_ids)
194  if (_excluded_boundaries.find(bid) != _excluded_boundaries.end())
195  return true;
196 
197  // make sure we will also break on the neighboring side
198  unsigned int neighbor_side = _neighbor_elem->which_neighbor_am_i(_current_elem);
199  boundary_ids = _mesh.getBoundaryIDs(_neighbor_elem, neighbor_side);
200  for (auto bid : boundary_ids)
201  if (_excluded_boundaries.find(bid) != _excluded_boundaries.end())
202  return true;
203 
204  return false;
205 }
206 
207 void
208 DGKernelBase::prepareShapes(const unsigned int var_num)
209 {
211 }
MooseMesh & _mesh
Reference to this Kernel&#39;s mesh object.
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar...
Definition: DGKernelBase.C:159
std::set< BoundaryID > _excluded_boundaries
Broken boundaries.
Definition: DGKernelBase.h:148
static InputParameters validParams()
Factory constructor initializes all internal references needed for residual computation.
Definition: DGKernelBase.C:26
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
static InputParameters validParams()
const std::set< BoundaryID > & meshSidesetIds() const
Returns a read-only reference to the set of sidesets currently present in the Mesh.
Definition: MooseMesh.C:2933
virtual void precalculateOffDiagJacobian(unsigned int)
T * get(const std::unique_ptr< T > &u)
The MooseUtils::get() specializations are used to support making forwards-compatible code changes fro...
Definition: MooseUtils.h:1147
virtual void precalculateResidual()
DGKernelBase(const InputParameters &parameters)
Definition: DGKernelBase.C:72
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
/class BoundaryRestrictable /brief Provides functionality for limiting the object to certain boundary...
virtual void computeJacobian() override
Computes the jacobian for the current side.
Definition: DGKernelBase.C:138
void addRelationshipManager(const std::string &name, Moose::RelationshipManagerType rm_type, Moose::RelationshipManagerInputParameterCallback input_parameter_callback=nullptr)
Tells MOOSE about a RelationshipManager that this object needs.
THREAD_ID _tid
The thread ID for this kernel.
virtual void computeElemNeighResidual(Moose::DGResidualType type)=0
Computes the residual for this element or the neighbor.
virtual const std::string & name() const
Get the name of the class.
Definition: MooseBase.h:56
const MooseVariableFieldBase & getVariable(unsigned int jvar_num) const
Retrieve the variable object from our system associated with jvar_num.
static InputParameters validParams()
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
void registerBase(const std::string &value)
This method must be called from every base "Moose System" to create linkage with the Action System...
virtual void precalculateJacobian()
void prepareShapes(unsigned int var_num) override final
Prepare shape functions.
Definition: DGKernelBase.C:208
SubProblem & _subproblem
Reference to this kernel&#39;s SubProblem.
bool excludeBoundary() const
Check current element if it contains broken boundary.
Definition: DGKernelBase.C:187
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 ...
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:62
static Threads::spin_mutex _resid_vars_mutex
Definition: DGKernelBase.h:140
static InputParameters validParams()
This is a base class for objects that can provide residual contributions for both local and neighbor ...
virtual void computeResidual() override
Computes the residual for the current side.
Definition: DGKernelBase.C:123
const Elem *const & _current_elem
Current element.
Definition: DGKernelBase.h:84
An interface that restricts an object to subdomains via the &#39;blocks&#39; input parameter.
const unsigned int & _current_side
Current side.
Definition: DGKernelBase.h:96
This interface is designed for DGKernel, InternalSideUserObject, InterfaceUserObject, where material properties on a side of both its primary side (face) and its secondary side (neighbor) all required.
const Elem *const & _neighbor_elem
The neighboring element.
Definition: DGKernelBase.h:90
std::vector< BoundaryID > getBoundaryIDs(const Elem *const elem, const unsigned short int side) const
Returns a vector of boundary IDs for the requested element on the requested side. ...
Definition: MooseMesh.C:2719
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
virtual const MooseVariableBase & variable() const =0
Returns the variable that this object operates on.
virtual void computeElemNeighJacobian(Moose::DGJacobianType type)=0
Computes the element/neighbor-element/neighbor Jacobian.
virtual void prepareFaceShapes(unsigned int var, const THREAD_ID tid)=0
static Threads::spin_mutex _jacoby_vars_mutex
Definition: DGKernelBase.h:141
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)
This method takes a space delimited list of parameter names and adds them to the specified group name...
virtual void computeOffDiagElemNeighJacobian(Moose::DGJacobianType type, const MooseVariableFEBase &jvar)=0
Computes the element-element off-diagonal Jacobian.