www.mooseframework.org
InterfaceKernel.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 "InterfaceKernel.h"
16 
17 // MOOSE includes
18 #include "Assembly.h"
19 #include "MooseVariable.h"
20 #include "SystemBase.h"
21 
22 #include "libmesh/quadrature.h"
23 
24 template <>
27 {
32  params.addRequiredParam<NonlinearVariableName>(
33  "variable", "The name of the variable that this boundary condition applies to");
34  params.addParam<bool>("use_displaced_mesh",
35  false,
36  "Whether or not this object should use the "
37  "displaced mesh for computation. Note that in "
38  "the case this is true but no displacements "
39  "are provided in the Mesh block the "
40  "undisplaced mesh will still be used.");
41  params.addParamNamesToGroup("use_displaced_mesh", "Advanced");
42 
43  params.declareControllable("enable");
44  params.addRequiredCoupledVar("neighbor_var", "The variable on the other side of the interface.");
45  params.set<std::string>("_moose_base") = "InterfaceKernel";
46  params.addParam<std::vector<AuxVariableName>>(
47  "save_in",
48  "The name of auxiliary variables to save this Kernel's residual contributions to. "
49  " Everything about that variable must match everything about this variable (the "
50  "type, what blocks it's on, etc.)");
51  params.addParam<std::vector<AuxVariableName>>(
52  "diag_save_in",
53  "The name of auxiliary variables to save this Kernel's diagonal Jacobian "
54  "contributions to. Everything about that variable must match everything "
55  "about this variable (the type, what blocks it's on, etc.)");
56 
57  MultiMooseEnum save_in_var_side("m s");
58  params.addParam<MultiMooseEnum>(
59  "save_in_var_side",
60  save_in_var_side,
61  "This parameter must exist if save_in variables are specified and must have the same length "
62  "as save_in. This vector specifies whether the corresponding aux_var should save-in "
63  "residual contributions from the master ('m') or slave side ('s').");
64  params.addParam<MultiMooseEnum>(
65  "diag_save_in_var_side",
66  save_in_var_side,
67  "This parameter must exist if diag_save_in variables are specified and must have the same "
68  "length as diag_save_in. This vector specifies whether the corresponding aux_var should "
69  "save-in jacobian contributions from the master ('m') or slave side ('s').");
70 
71  return params;
72 }
73 
74 // Static mutex definitions
75 Threads::spin_mutex InterfaceKernel::_resid_vars_mutex;
76 Threads::spin_mutex InterfaceKernel::_jacoby_vars_mutex;
77 
79  : MooseObject(parameters),
80  BoundaryRestrictable(this, false), // false for _not_ nodal
81  SetupInterface(this),
82  TransientInterface(this),
83  FunctionInterface(this),
84  UserObjectInterface(this),
86  Restartable(parameters, "InterfaceKernels"),
87  ZeroInterface(parameters),
88  MeshChangedInterface(parameters),
89  TwoMaterialPropertyInterface(this, boundaryIDs()),
90  _subproblem(*parameters.get<SubProblem *>("_subproblem")),
91  _sys(*parameters.get<SystemBase *>("_sys")),
92  _tid(parameters.get<THREAD_ID>("_tid")),
93  _assembly(_subproblem.assembly(_tid)),
94  _var(_sys.getVariable(_tid, parameters.get<NonlinearVariableName>("variable"))),
95  _mesh(_subproblem.mesh()),
96  _current_elem(_assembly.elem()),
97  _current_elem_volume(_assembly.elemVolume()),
98  _neighbor_elem(_assembly.neighbor()),
99  _current_side(_assembly.side()),
100  _current_side_elem(_assembly.sideElem()),
101  _current_side_volume(_assembly.sideElemVolume()),
102  _coord_sys(_assembly.coordSystem()),
103  _q_point(_assembly.qPointsFace()),
104  _qrule(_assembly.qRuleFace()),
105  _JxW(_assembly.JxWFace()),
106  _coord(_assembly.coordTransformation()),
107  _u(_is_implicit ? _var.sln() : _var.slnOld()),
108  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld()),
109  _phi(_assembly.phiFace()),
110  _grad_phi(_assembly.gradPhiFace()),
111  _test(_var.phiFace()),
112  _grad_test(_var.gradPhiFace()),
113  _normals(_var.normals()),
114  _neighbor_var(*getVar("neighbor_var", 0)),
115  _neighbor_value(_neighbor_var.slnNeighbor()),
116  _grad_neighbor_value(_neighbor_var.gradSlnNeighbor()),
117  _phi_neighbor(_assembly.phiFaceNeighbor()),
118  _grad_phi_neighbor(_assembly.gradPhiFaceNeighbor()),
119  _test_neighbor(_neighbor_var.phiFaceNeighbor()),
120  _grad_test_neighbor(_neighbor_var.gradPhiFaceNeighbor()),
121  _save_in_var_side(parameters.get<MultiMooseEnum>("save_in_var_side")),
122  _save_in_strings(parameters.get<std::vector<AuxVariableName>>("save_in")),
123  _diag_save_in_var_side(parameters.get<MultiMooseEnum>("diag_save_in_var_side")),
124  _diag_save_in_strings(parameters.get<std::vector<AuxVariableName>>("diag_save_in"))
125 
126 {
127  if (!parameters.isParamValid("boundary"))
128  mooseError(
129  "In order to use an interface kernel, you must specify a boundary where it will live.");
130 
131  if (parameters.isParamSetByUser("save_in"))
132  {
133  if (_save_in_strings.size() != _save_in_var_side.size())
134  mooseError("save_in and save_in_var_side must be the same length");
135  else
136  {
137  for (unsigned i = 0; i < _save_in_strings.size(); ++i)
138  {
140 
142  mooseError("Trying to use solution variable " + _save_in_strings[i] +
143  " as a save_in variable in " + name());
144 
145  if (_save_in_var_side[i] == "m")
146  {
147  if (var->feType() != _var.feType())
148  mooseError(
149  "Error in " + name() +
150  ". There is a mismatch between the fe_type of the save-in Auxiliary variable "
151  "and the fe_type of the the master side nonlinear "
152  "variable this interface kernel object is acting on.");
153  _master_save_in_residual_variables.push_back(var);
154  }
155  else
156  {
157  if (var->feType() != _neighbor_var.feType())
158  mooseError(
159  "Error in " + name() +
160  ". There is a mismatch between the fe_type of the save-in Auxiliary variable "
161  "and the fe_type of the the slave side nonlinear "
162  "variable this interface kernel object is acting on.");
163  _slave_save_in_residual_variables.push_back(var);
164  }
165 
168  }
169  }
170  }
171 
174 
175  if (parameters.isParamSetByUser("diag_save_in"))
176  {
178  mooseError("diag_save_in and diag_save_in_var_side must be the same length");
179  else
180  {
181  for (unsigned i = 0; i < _diag_save_in_strings.size(); ++i)
182  {
184 
186  mooseError("Trying to use solution variable " + _diag_save_in_strings[i] +
187  " as a save_in variable in " + name());
188 
189  if (_diag_save_in_var_side[i] == "m")
190  {
191  if (var->feType() != _var.feType())
192  mooseError(
193  "Error in " + name() +
194  ". There is a mismatch between the fe_type of the save-in Auxiliary variable "
195  "and the fe_type of the the master side nonlinear "
196  "variable this interface kernel object is acting on.");
197  _master_save_in_jacobian_variables.push_back(var);
198  }
199  else
200  {
201  if (var->feType() != _neighbor_var.feType())
202  mooseError(
203  "Error in " + name() +
204  ". There is a mismatch between the fe_type of the save-in Auxiliary variable "
205  "and the fe_type of the the slave side nonlinear "
206  "variable this interface kernel object is acting on.");
207  _slave_save_in_jacobian_variables.push_back(var);
208  }
209 
212  }
213  }
214  }
215 
218 }
219 
220 void
222 {
223  bool is_elem;
224  if (type == Moose::Element)
225  is_elem = true;
226  else
227  is_elem = false;
228 
229  const VariableTestValue & test_space = is_elem ? _test : _test_neighbor;
230  DenseVector<Number> & re = is_elem ? _assembly.residualBlock(_var.number())
232 
233  _local_re.resize(re.size());
234  _local_re.zero();
235 
236  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
237  for (_i = 0; _i < test_space.size(); _i++)
238  _local_re(_i) += _JxW[_qp] * _coord[_qp] * computeQpResidual(type);
239 
240  re += _local_re;
241 
242  if (_has_master_residuals_saved_in && is_elem)
243  {
244  Threads::spin_mutex::scoped_lock lock(_resid_vars_mutex);
245  for (const auto & var : _master_save_in_residual_variables)
246  {
247  var->sys().solution().add_vector(_local_re, var->dofIndices());
248  }
249  }
250  else if (_has_slave_residuals_saved_in && !is_elem)
251  {
252  Threads::spin_mutex::scoped_lock lock(_resid_vars_mutex);
253  for (const auto & var : _slave_save_in_residual_variables)
254  var->sys().solution().add_vector(_local_re, var->dofIndicesNeighbor());
255  }
256 }
257 
258 void
260 {
261  // Compute the residual for this element
263 
264  // Compute the residual for the neighbor
266 }
267 
268 void
270 {
271  const VariableTestValue & test_space =
273  const VariableTestValue & loc_phi =
275  DenseMatrix<Number> & Kxx =
276  type == Moose::ElementElement
278  : type == Moose::ElementNeighbor
281  : type == Moose::NeighborElement
287 
288  _local_kxx.resize(Kxx.m(), Kxx.n());
289  _local_kxx.zero();
290 
291  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
292  for (_i = 0; _i < test_space.size(); _i++)
293  for (_j = 0; _j < loc_phi.size(); _j++)
294  _local_kxx(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpJacobian(type);
295 
296  Kxx += _local_kxx;
297 
299  {
300  auto rows = _local_kxx.m();
301  DenseVector<Number> diag(rows);
302  for (decltype(rows) i = 0; i < rows; i++)
303  diag(i) = _local_kxx(i, i);
304 
305  Threads::spin_mutex::scoped_lock lock(_jacoby_vars_mutex);
306  for (const auto & var : _master_save_in_jacobian_variables)
307  var->sys().solution().add_vector(diag, var->dofIndices());
308  }
310  {
311  auto rows = _local_kxx.m();
312  DenseVector<Number> diag(rows);
313  for (decltype(rows) i = 0; i < rows; i++)
314  diag(i) = _local_kxx(i, i);
315 
316  Threads::spin_mutex::scoped_lock lock(_jacoby_vars_mutex);
317  for (const auto & var : _slave_save_in_jacobian_variables)
318  var->sys().solution().add_vector(diag, var->dofIndicesNeighbor());
319  }
320 }
321 
322 void
324 {
327 }
328 
329 void
331 {
332  const VariableTestValue & test_space =
334  const VariableTestValue & loc_phi =
336  DenseMatrix<Number> & Kxx =
337  type == Moose::ElementElement
339  : type == Moose::ElementNeighbor
341  : type == Moose::NeighborElement
346 
347  // Prevent calling of Jacobian computation if jvar doesn't lie in the current block
348  if ((Kxx.m() == test_space.size()) && (Kxx.n() == loc_phi.size()))
349  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
350  for (_i = 0; _i < test_space.size(); _i++)
351  for (_j = 0; _j < loc_phi.size(); _j++)
352  Kxx(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(type, jvar);
353 }
354 
355 void
357 {
358  bool is_jvar_not_interface_var = true;
359  if (jvar == _var.number())
360  {
362  is_jvar_not_interface_var = false;
363  }
364  if (jvar == _neighbor_var.number())
365  {
367  is_jvar_not_interface_var = false;
368  }
369 
370  if (is_jvar_not_interface_var)
371  {
374  }
375 }
376 
377 void
379 {
380  bool is_jvar_not_interface_var = true;
381  if (jvar == _var.number())
382  {
384  is_jvar_not_interface_var = false;
385  }
386  if (jvar == _neighbor_var.number())
387  {
389  is_jvar_not_interface_var = false;
390  }
391 
392  if (is_jvar_not_interface_var)
393  {
396  }
397 }
398 
399 Real
401 {
402  return 0.;
403 }
404 
407 {
408  return _var;
409 }
410 
411 const MooseVariable &
413 {
414  return _neighbor_var;
415 }
416 
417 SubProblem &
419 {
420  return _subproblem;
421 }
422 
423 const Real &
425 {
426  return _assembly.neighborVolume();
427 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
const VariablePhiValue & _phi
shape function
InterfaceKernel(const InputParameters &parameters)
MultiMooseEnum _save_in_var_side
MultiMooseEnum specifying whether residual save-in aux variables correspond to master or slave side...
A class for creating restricted objects.
Definition: Restartable.h:31
bool isParamSetByUser(const std::string &name) const
Method returns true if the parameter was by the user.
virtual void computeOffDiagElemNeighJacobian(Moose::DGJacobianType type, unsigned int jvar)
Using the passed DGJacobian type, selects the correct test function and trial function spaces and jac...
bool _has_slave_residuals_saved_in
Whether there are slave residual aux variables.
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...
const MooseArray< Real > & _coord
Coordinate transformation value; relevant in axisymmetric simulations for example.
Class for stuff related to variables.
Definition: MooseVariable.h:43
MultiMooseEnum _diag_save_in_var_side
MultiMooseEnum specifying whether jacobian save-in aux variables correspond to master or slave side...
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. ...
virtual Real computeQpOffDiagJacobian(Moose::DGJacobianType type, unsigned int jvar)
compute off-diagonal jacobian components at quadrature points
MooseVariable & _neighbor_var
Coupled neighbor variable.
DenseVector< Number > _local_re
Holds residual entries as they are accumulated by this InterfaceKernel.
SubProblem & subProblem()
Return a reference to the subproblem.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int _qp
Current quadrature point.
void addMooseVariableDependency(MooseVariable *var)
Call this function to add the passed in MooseVariable as a variable that this object depends on...
/class BoundaryRestrictable /brief Provides functionality for limiting the object to certain boundary...
virtual Real computeQpResidual(Moose::DGResidualType type)=0
Compute residuals at quadrature points.
const Real & neighborVolume()
Returns the reference to the current neighbor volume.
Definition: Assembly.C:353
const MooseArray< Real > & _JxW
Elemtn Jacobian/quadrature weight.
unsigned int _j
DGResidualType
Definition: MooseTypes.h:184
Base class for a system (of equations)
Definition: SystemBase.h:91
virtual void computeElemNeighJacobian(Moose::DGJacobianType type)
Using the passed DGJacobian type, selects the correct test function and trial function spaces and jac...
std::vector< MooseVariable * > _master_save_in_residual_variables
The aux variables to save the master residual contributions to.
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...
DenseMatrix< Number > _local_kxx
Holds residual entries as they are accumulated by this InterfaceKernel.
static Threads::spin_mutex _resid_vars_mutex
Mutex that prevents multiple threads from saving into the residual aux_var at the same time...
virtual void computeElemNeighResidual(Moose::DGResidualType type)
Using the passed DGResidual type, selects the correct test function space and residual block...
unsigned int _i
Index for test and trial functions.
const VariablePhiValue & _phi_neighbor
Side neighbor shape function.
bool _has_slave_jacobians_saved_in
Whether there are slave jacobian aux variables.
Interface for objects that needs transient capabilities.
bool isParamValid(const std::string &name) const
This method returns parameters that have been initialized in one fashion or another, i.e.
std::vector< AuxVariableName > _diag_save_in_strings
The names of the aux variables that will be used to save-in jacobians (includes both master and slave...
Interface for notifications that the mesh has changed.
virtual void computeResidual()
Computes the residual for the current side.
virtual void computeNeighborOffDiagJacobian(unsigned int jvar)
Selects the correct Jacobian type and routine to call for the slave variable jacobian.
Every object that can be built by the factory should be derived from this class.
Definition: MooseObject.h:36
const Real & getNeighborElemVolume()
The volume of the current neighbor.
InputParameters validParams< MeshChangedInterface >()
SystemBase & _sys
Reference to the nonlinear system.
const MooseVariable & neighborVariable() const
The neighbor variable number that this interface kernel operates on.
Interface for objects that need to use UserObjects.
std::vector< MooseVariable * > _master_save_in_jacobian_variables
The aux variables to save the diagonal Jacobian contributions of the master variables to...
virtual bool hasVariable(const std::string &var_name)
Query a system for a variable.
Definition: SystemBase.C:556
const VariableTestValue & _test
Side shape function.
InputParameters validParams< InterfaceKernel >()
InputParameters validParams< MooseObject >()
Definition: MooseObject.C:22
DGJacobianType
Definition: MooseTypes.h:190
MooseVariable & variable() const
The master variable that this interface kernel operates on.
const VariableTestValue & _test_neighbor
Side neighbor test function.
MooseVariable & _var
The master side MooseVariable.
MatType type
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpJacobian(Moose::DGJacobianType type)=0
Compute jacobians at quadrature points.
bool _has_master_jacobians_saved_in
Whether there are master jacobian aux variables.
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 >()
Assembly & _assembly
Problem assembly.
QBase *& _qrule
Quadrature rule.
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
unsigned int size() const
Return the number of items in the MultiMooseEnum.
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
THREAD_ID _tid
The thread ID.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
bool _has_master_residuals_saved_in
Whether there are master residual aux variables.
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:887
static Threads::spin_mutex _jacoby_vars_mutex
Mutex that prevents multiple threads from saving into the jacobian aux_var at the same time...
std::vector< AuxVariableName > _save_in_strings
The names of the aux variables that will be used to save-in residuals (includes both master and slave...
virtual void computeElementOffDiagJacobian(unsigned int jvar)
Selects the correct Jacobian type and routine to call for the master variable jacobian.
SystemBase & sys()
Get the system this variable is part of.
std::vector< MooseVariable * > _slave_save_in_jacobian_variables
The aux variables to save the diagonal Jacobian contributions of the slave variables to...
std::vector< MooseVariable * > _slave_save_in_residual_variables
The aux variables to save the slave contributions to.
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.
SubProblem & _subproblem
Reference to the controlling finite element problem.
InputParameters validParams< TransientInterface >()
unsigned int THREAD_ID
Definition: MooseTypes.h:79
virtual void computeJacobian()
Computes the jacobian for the current side.