www.mooseframework.org
NodalKernel.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 "NodalKernel.h"
11 #include "Problem.h"
12 #include "SubProblem.h"
13 #include "SystemBase.h"
14 #include "MooseVariableFE.h"
15 #include "Assembly.h"
16 
19 {
21 }
22 
24  : NodalKernelBase(parameters),
25  _u(_var.dofValues()),
26  _save_in_strings(parameters.get<std::vector<AuxVariableName>>("save_in")),
27  _diag_save_in_strings(parameters.get<std::vector<AuxVariableName>>("diag_save_in"))
28 
29 {
30  _save_in.resize(_save_in_strings.size());
31  _diag_save_in.resize(_diag_save_in_strings.size());
32 
33  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
34  {
36 
37  if (var->feType() != _var.feType())
38  paramError(
39  "save_in",
40  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
42 
43  _save_in[i] = var;
46  }
47 
48  _has_save_in = _save_in.size() > 0;
49 
50  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
51  {
53 
54  if (var->feType() != _var.feType())
55  paramError(
56  "diag_save_in",
57  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
59 
60  _diag_save_in[i] = var;
63  }
64 
65  _has_diag_save_in = _diag_save_in.size() > 0;
66 }
67 
68 void
70 {
71  if (_var.isNodalDefined())
72  {
73  const dof_id_type & dof_idx = _var.nodalDofIndex();
74  _qp = 0;
75  const Real res = computeQpResidual();
77  std::array<Real, 1>{{res}},
78  std::array<dof_id_type, 1>{{dof_idx}},
80 
81  if (_has_save_in)
82  for (const auto & var : _save_in)
83  var->sys().solution().add(var->nodalDofIndex(), res);
84  }
85 }
86 
87 void
89 {
90  if (_var.isNodalDefined())
91  {
92  _qp = 0;
93  const Real cached_val = computeQpJacobian();
94  const dof_id_type cached_row = _var.nodalDofIndex();
95 
96  addJacobianElement(_assembly, cached_val, cached_row, cached_row, _var.scalingFactor());
97 
99  {
100  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
101  for (const auto & var : _diag_save_in)
102  var->sys().solution().add(var->nodalDofIndex(), cached_val);
103  }
104  }
105 }
106 
107 void
108 NodalKernel::computeOffDiagJacobian(const unsigned int jvar_num)
109 {
110  const auto & jvar = getVariable(jvar_num);
111 
112  if (_var.isNodalDefined())
113  {
114  if (jvar.number() == _var.number())
115  computeJacobian();
116  else
117  {
118  _qp = 0;
119  const Real cached_val = computeQpOffDiagJacobian(jvar.number());
120  const dof_id_type cached_row = _var.nodalDofIndex();
121 
122  // Note: this only works for equal order Lagrange variables...
123  const dof_id_type cached_col = _current_node->dof_number(_sys.number(), jvar.number(), 0);
124 
125  addJacobianElement(_assembly, cached_val, cached_row, cached_col, _var.scalingFactor());
126  }
127  }
128 }
129 
130 Real
132 {
133  return 0.;
134 }
135 
136 Real
138 {
139  return 0.;
140 }
std::vector< AuxVariableName > _save_in_strings
Definition: NodalKernel.h:77
std::string incompatVarMsg(MooseVariableFieldBase &var1, MooseVariableFieldBase &var2)
Builds and returns a string of the form:
Definition: MooseError.C:23
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpOffDiagJacobian(unsigned int jvar)
This is the virtual that derived classes should override for computing an off-diagonal jacobian compo...
Definition: NodalKernel.C:137
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
void addResiduals(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided incoming residuals corresponding to the provided dof indices.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual void computeJacobian() override
Compute the Jacobian at one node.
Definition: NodalKernel.C:88
MooseVariable & _var
variable this works on
THREAD_ID _tid
The thread ID for this kernel.
const FEType & feType() const
Get the type of finite element object.
const MooseVariableFieldBase & getVariable(unsigned int jvar_num) const
Retrieve the variable object from our system associated with jvar_num.
SystemBase & _sys
Reference to the EquationSystem object.
std::vector< AuxVariableName > _diag_save_in_strings
Definition: NodalKernel.h:82
SubProblem & _subproblem
Reference to this kernel&#39;s SubProblem.
virtual void computeResidual() override
Compute the residual at the current node.
Definition: NodalKernel.C:69
virtual Real computeQpJacobian()
The user can override this function to compute the "on-diagonal" Jacobian contribution.
Definition: NodalKernel.C:131
virtual MooseVariable & getStandardVariable(const THREAD_ID tid, const std::string &var_name)=0
Returns the variable reference for requested MooseVariable which may be in any system.
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 ...
unsigned int number() const
Gets the number of this system.
Definition: SystemBase.C:1125
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: NodalKernel.h:80
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: NodalKernel.h:75
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
static InputParameters validParams()
Class constructor.
Base class for creating new types of nodal kernels.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual bool isNodalDefined() const override
Is this variable defined at nodes.
const dof_id_type & nodalDofIndex() const override
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:163
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:169
std::vector< MooseVariableFEBase * > _save_in
Definition: NodalKernel.h:76
const Node *const & _current_node
current node being processed
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: NodalKernel.h:81
static InputParameters validParams()
Class constructor.
Definition: NodalKernel.C:18
void addJacobianElement(Assembly &assembly, Real value, dof_id_type row_index, dof_id_type column_index, Real scaling_factor)
Add into a single Jacobian element.
SystemBase & sys()
Get the system this variable is part of.
virtual void computeOffDiagJacobian(unsigned int jvar) override
Compute the off-diagonal Jacobian at one node.
Definition: NodalKernel.C:108
unsigned int _qp
Quadrature point index.
virtual Real computeQpResidual()=0
The user can override this function to compute the residual at a node.
void scalingFactor(const std::vector< Real > &factor)
Set the scaling factor for this variable.
NodalKernel(const InputParameters &parameters)
Definition: NodalKernel.C:23
uint8_t dof_id_type