www.mooseframework.org
InterfaceKernel.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 "InterfaceKernel.h"
11 
12 // MOOSE includes
13 #include "Assembly.h"
14 #include "MooseVariableFE.h"
15 #include "SystemBase.h"
16 
17 #include "libmesh/quadrature.h"
18 
19 template <>
22 {
24  params.registerBase("InterfaceKernel");
25  return params;
26 }
27 
28 template <>
31 {
33  params.registerBase("VectorInterfaceKernel");
34  return params;
35 }
36 
37 template <typename T>
39  : InterfaceKernelBase(parameters),
41  false,
43  std::is_same<T, Real>::value
46  _var(*this->mooseVariable()),
47  _normals(_assembly.normals()),
48  _u(_is_implicit ? _var.sln() : _var.slnOld()),
49  _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld()),
50  _phi(_assembly.phiFace(_var)),
51  _grad_phi(_assembly.gradPhiFace(_var)),
52  _test(_var.phiFace()),
53  _grad_test(_var.gradPhiFace()),
54  _neighbor_var(*getVarHelper<T>("neighbor_var", 0)),
55  _neighbor_value(_is_implicit ? _neighbor_var.slnNeighbor() : _neighbor_var.slnOldNeighbor()),
56  _grad_neighbor_value(_neighbor_var.gradSlnNeighbor()),
57  _phi_neighbor(_assembly.phiFaceNeighbor(_neighbor_var)),
58  _grad_phi_neighbor(_assembly.gradPhiFaceNeighbor(_neighbor_var)),
59  _test_neighbor(_neighbor_var.phiFaceNeighbor()),
60  _grad_test_neighbor(_neighbor_var.gradPhiFaceNeighbor())
61 
62 {
64 
65  if (!parameters.isParamValid("boundary"))
66  mooseError(
67  "In order to use an interface kernel, you must specify a boundary where it will live.");
68 
69  if (parameters.isParamSetByUser("save_in"))
70  {
71  if (_save_in_strings.size() != _save_in_var_side.size())
72  mooseError("save_in and save_in_var_side must be the same length");
73  else
74  {
75  for (unsigned i = 0; i < _save_in_strings.size(); ++i)
76  {
78 
80  mooseError("Trying to use solution variable " + _save_in_strings[i] +
81  " as a save_in variable in " + name());
82 
83  if (_save_in_var_side[i] == "m")
84  {
85  if (var->feType() != _var.feType())
86  mooseError(
87  "Error in " + name() +
88  ". There is a mismatch between the fe_type of the save-in Auxiliary variable "
89  "and the fe_type of the the master side nonlinear "
90  "variable this interface kernel object is acting on.");
92  }
93  else
94  {
95  if (var->feType() != _neighbor_var.feType())
96  mooseError(
97  "Error in " + name() +
98  ". There is a mismatch between the fe_type of the save-in Auxiliary variable "
99  "and the fe_type of the the slave side nonlinear "
100  "variable this interface kernel object is acting on.");
101  _slave_save_in_residual_variables.push_back(var);
102  }
103 
106  }
107  }
108  }
109 
112 
113  if (parameters.isParamSetByUser("diag_save_in"))
114  {
116  mooseError("diag_save_in and diag_save_in_var_side must be the same length");
117  else
118  {
119  for (unsigned i = 0; i < _diag_save_in_strings.size(); ++i)
120  {
122 
124  mooseError("Trying to use solution variable " + _diag_save_in_strings[i] +
125  " as a save_in variable in " + name());
126 
127  if (_diag_save_in_var_side[i] == "m")
128  {
129  if (var->feType() != _var.feType())
130  mooseError(
131  "Error in " + name() +
132  ". There is a mismatch between the fe_type of the save-in Auxiliary variable "
133  "and the fe_type of the the master side nonlinear "
134  "variable this interface kernel object is acting on.");
135  _master_save_in_jacobian_variables.push_back(var);
136  }
137  else
138  {
139  if (var->feType() != _neighbor_var.feType())
140  mooseError(
141  "Error in " + name() +
142  ". There is a mismatch between the fe_type of the save-in Auxiliary variable "
143  "and the fe_type of the the slave side nonlinear "
144  "variable this interface kernel object is acting on.");
145  _slave_save_in_jacobian_variables.push_back(var);
146  }
147 
150  }
151  }
152  }
153 
156 }
157 
158 template <typename T>
159 void
161 {
162  bool is_elem;
163  if (type == Moose::Element)
164  is_elem = true;
165  else
166  is_elem = false;
167 
168  const TemplateVariableTestValue & test_space = is_elem ? _test : _test_neighbor;
169 
170  if (is_elem)
171  prepareVectorTag(_assembly, _var.number());
172  else
173  prepareVectorTagNeighbor(_assembly, _neighbor_var.number());
174 
175  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
176  for (_i = 0; _i < test_space.size(); _i++)
177  _local_re(_i) += _JxW[_qp] * _coord[_qp] * computeQpResidual(type);
178 
179  accumulateTaggedLocalResidual();
180 
181  if (_has_master_residuals_saved_in && is_elem)
182  {
183  Threads::spin_mutex::scoped_lock lock(_resid_vars_mutex);
184  for (const auto & var : _master_save_in_residual_variables)
185  {
186  var->sys().solution().add_vector(_local_re, var->dofIndices());
187  }
188  }
189  else if (_has_slave_residuals_saved_in && !is_elem)
190  {
191  Threads::spin_mutex::scoped_lock lock(_resid_vars_mutex);
192  for (const auto & var : _slave_save_in_residual_variables)
193  var->sys().solution().add_vector(_local_re, var->dofIndicesNeighbor());
194  }
195 }
196 
197 template <typename T>
198 void
200 {
201  // Compute the residual for this element
202  computeElemNeighResidual(Moose::Element);
203 
204  // Compute the residual for the neighbor
205  computeElemNeighResidual(Moose::Neighbor);
206 }
207 
208 template <typename T>
209 void
211 {
212  const TemplateVariableTestValue & test_space =
213  (type == Moose::ElementElement || type == Moose::ElementNeighbor) ? _test : _test_neighbor;
214  const TemplateVariableTestValue & loc_phi =
215  (type == Moose::ElementElement || type == Moose::NeighborElement) ? _phi : _phi_neighbor;
216 
217  unsigned int ivar, jvar;
218 
219  switch (type)
220  {
222  ivar = jvar = _var.number();
223  break;
225  ivar = _var.number(), jvar = _neighbor_var.number();
226  break;
228  ivar = _neighbor_var.number(), jvar = _var.number();
229  break;
231  ivar = _neighbor_var.number(), jvar = _neighbor_var.number();
232  break;
233  default:
234  mooseError("Unknown DGJacobianType ", type);
235  }
236 
238  prepareMatrixTag(_assembly, ivar, jvar);
239  else
240  prepareMatrixTagNeighbor(_assembly, ivar, jvar, type);
241 
242  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
243  for (_i = 0; _i < test_space.size(); _i++)
244  for (_j = 0; _j < loc_phi.size(); _j++)
245  _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpJacobian(type);
246 
247  accumulateTaggedLocalMatrix();
248 
249  if (_has_master_jacobians_saved_in && type == Moose::ElementElement)
250  {
251  auto rows = _local_ke.m();
252  DenseVector<Number> diag(rows);
253  for (decltype(rows) i = 0; i < rows; i++)
254  diag(i) = _local_ke(i, i);
255 
256  Threads::spin_mutex::scoped_lock lock(_jacoby_vars_mutex);
257  for (const auto & var : _master_save_in_jacobian_variables)
258  var->sys().solution().add_vector(diag, var->dofIndices());
259  }
260  else if (_has_slave_jacobians_saved_in && type == Moose::NeighborNeighbor)
261  {
262  auto rows = _local_ke.m();
263  DenseVector<Number> diag(rows);
264  for (decltype(rows) i = 0; i < rows; i++)
265  diag(i) = _local_ke(i, i);
266 
267  Threads::spin_mutex::scoped_lock lock(_jacoby_vars_mutex);
268  for (const auto & var : _slave_save_in_jacobian_variables)
269  var->sys().solution().add_vector(diag, var->dofIndicesNeighbor());
270  }
271 }
272 
273 template <typename T>
274 void
276 {
277  computeElemNeighJacobian(Moose::ElementElement);
278  computeElemNeighJacobian(Moose::NeighborNeighbor);
279 }
280 
281 template <typename T>
282 void
284  unsigned int jvar)
285 {
286  const TemplateVariableTestValue & test_space =
287  (type == Moose::ElementElement || type == Moose::ElementNeighbor) ? _test : _test_neighbor;
288  const TemplateVariableTestValue & loc_phi =
289  (type == Moose::ElementElement || type == Moose::NeighborElement) ? _phi : _phi_neighbor;
290 
291  unsigned int ivar;
292 
294  ivar = _var.number();
295  else
296  ivar = _neighbor_var.number();
297 
299  prepareMatrixTag(_assembly, ivar, jvar);
300  else
301  prepareMatrixTagNeighbor(_assembly, ivar, jvar, type);
302 
303  // Prevent calling of Jacobian computation if jvar doesn't lie in the current block
304  if ((_local_ke.m() == test_space.size()) && (_local_ke.n() == loc_phi.size()))
305  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
306  for (_i = 0; _i < test_space.size(); _i++)
307  for (_j = 0; _j < loc_phi.size(); _j++)
308  _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(type, jvar);
309 
310  accumulateTaggedLocalMatrix();
311 }
312 
313 template <typename T>
314 void
316 {
317  bool is_jvar_not_interface_var = true;
318  if (jvar == _var.number())
319  {
320  computeElemNeighJacobian(Moose::ElementElement);
321  is_jvar_not_interface_var = false;
322  }
323  if (jvar == _neighbor_var.number())
324  {
325  computeElemNeighJacobian(Moose::ElementNeighbor);
326  is_jvar_not_interface_var = false;
327  }
328 
329  if (is_jvar_not_interface_var)
330  {
331  computeOffDiagElemNeighJacobian(Moose::ElementElement, jvar);
332  computeOffDiagElemNeighJacobian(Moose::ElementNeighbor, jvar);
333  }
334 }
335 
336 template <typename T>
337 void
339 {
340  bool is_jvar_not_interface_var = true;
341  if (jvar == _var.number())
342  {
343  computeElemNeighJacobian(Moose::NeighborElement);
344  is_jvar_not_interface_var = false;
345  }
346  if (jvar == _neighbor_var.number())
347  {
348  computeElemNeighJacobian(Moose::NeighborNeighbor);
349  is_jvar_not_interface_var = false;
350  }
351 
352  if (is_jvar_not_interface_var)
353  {
354  computeOffDiagElemNeighJacobian(Moose::NeighborElement, jvar);
355  computeOffDiagElemNeighJacobian(Moose::NeighborNeighbor, jvar);
356  }
357 }
358 
359 // Explicitly instantiates the two versions of the InterfaceKernelTempl class
360 template class InterfaceKernelTempl<Real>;
virtual void computeElemNeighJacobian(Moose::DGJacobianType type) override
Using the passed DGJacobian type, selects the correct test function and trial function spaces and jac...
VarFieldType
Definition: MooseTypes.h:488
THREAD_ID _tid
The thread ID.
SystemBase & _sys
Reference to the nonlinear system.
virtual void computeElemNeighResidual(Moose::DGResidualType type) override
Using the passed DGResidual type, selects the correct test function space and residual block...
MultiMooseEnum _diag_save_in_var_side
MultiMooseEnum specifying whether jacobian save-in aux variables correspond to master or slave side...
virtual MooseVariable & getStandardVariable(THREAD_ID tid, const std::string &var_name)=0
Returns the variable reference for requested MooseVariable which may be in any system.
MooseVariableFE< T > & _neighbor_var
Coupled neighbor variable.
std::vector< MooseVariableFEBase * > _slave_save_in_jacobian_variables
The aux variables to save the diagonal Jacobian contributions of the slave variables to...
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:207
InterfaceKernel and VectorInterfaceKernel is responsible for interfacing physics across subdomains...
InterfaceKernelTempl(const InputParameters &parameters)
unsigned int size() const
Return the number of active items in the MultiMooseEnum.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
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...
bool _has_master_residuals_saved_in
Whether there are master residual aux variables.
InputParameters validParams< InterfaceKernelBase >()
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
DGResidualType
Definition: MooseTypes.h:509
MultiMooseEnum _save_in_var_side
MultiMooseEnum specifying whether residual save-in aux variables correspond to master or slave side...
const FEType & feType() const
Get the type of finite element object.
void addMooseVariableDependency(MooseVariableFEBase *var)
Call this function to add the passed in MooseVariableFEBase as a variable that this object depends on...
void registerBase(const std::string &value)
This method must be called from every base "Moose System" to create linkage with the Action System...
MooseVariableFE< T > * mooseVariable() const
Get the variable that this object is using.
Enhances MooseVariableInterface interface provide values from neighbor elements.
virtual void computeElementOffDiagJacobian(unsigned int jvar) override
Selects the correct Jacobian type and routine to call for the master variable jacobian.
virtual void computeOffDiagElemNeighJacobian(Moose::DGJacobianType type, unsigned int jvar) override
Using the passed DGJacobian type, selects the correct test function and trial function spaces and jac...
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseObject.h:65
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
virtual void computeJacobian() override
Computes the jacobian for the current side.
std::vector< MooseVariableFEBase * > _slave_save_in_residual_variables
The aux variables to save the slave contributions to.
bool _has_slave_residuals_saved_in
Whether there are slave residual aux variables.
InputParameters validParams< InterfaceKernel >()
virtual bool hasVariable(const std::string &var_name) const
Query a system for a variable.
Definition: SystemBase.C:668
InputParameters validParams< VectorInterfaceKernel >()
SubProblem & _subproblem
Reference to the controlling finite element problem.
DGJacobianType
Definition: MooseTypes.h:515
bool isParamSetByUser(const std::string &name) const
Method returns true if the parameter was by the user.
MatType type
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:178
bool _has_master_jacobians_saved_in
Whether there are master jacobian aux variables.
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...
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:59
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:184
InterfaceKernelBase is the base class for all InterfaceKernel type classes.
Definition: Moose.h:112
bool _has_slave_jacobians_saved_in
Whether there are slave jacobian aux variables.
std::vector< MooseVariableFEBase * > _master_save_in_jacobian_variables
The aux variables to save the diagonal Jacobian contributions of the master variables to...
MooseVariableFE< T > & _var
The master side MooseVariable.
SystemBase & sys()
Get the system this variable is part of.
std::vector< MooseVariableFEBase * > _master_save_in_residual_variables
The aux variables to save the master residual contributions to.
virtual void computeResidual() override
Computes the residual for the current side.
virtual void computeNeighborOffDiagJacobian(unsigned int jvar) override
Selects the correct Jacobian type and routine to call for the slave variable jacobian.
bool isParamValid(const std::string &name) const
This method returns parameters that have been initialized in one fashion or another, i.e.