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 
17 template <>
20 {
25  params += validParams<RandomInterface>();
27 
28  params.addRequiredParam<NonlinearVariableName>(
29  "variable", "The name of the variable that this boundary condition applies to");
30 
31  params.addParam<std::vector<AuxVariableName>>(
32  "save_in",
33  "The name of auxiliary variables to save this BC's residual contributions to. "
34  "Everything about that variable must match everything about this variable (the "
35  "type, what blocks it's on, etc.)");
36 
37  params.addParam<std::vector<AuxVariableName>>(
38  "diag_save_in",
39  "The name of auxiliary variables to save this BC's diagonal jacobian "
40  "contributions to. Everything about that variable must match everything "
41  "about this variable (the type, what blocks it's on, etc.)");
42 
43  params.addParam<bool>("use_displaced_mesh",
44  false,
45  "Whether or not this object should use the "
46  "displaced mesh for computation. Note that "
47  "in the case this is true but no "
48  "displacements are provided in the Mesh block "
49  "the undisplaced mesh will still be used.");
50  params.addParamNamesToGroup("use_displaced_mesh", "Advanced");
51 
52  params.declareControllable("enable");
53 
54  params.registerBase("NodalKernel");
55 
56  return params;
57 }
58 
60  : MooseObject(parameters),
61  BlockRestrictable(this),
62  BoundaryRestrictable(this, true), // true for applying to nodesets
63  SetupInterface(this),
64  FunctionInterface(this),
65  UserObjectInterface(this),
66  TransientInterface(this),
69  Restartable(this, "BCs"),
70  MeshChangedInterface(parameters),
71  RandomInterface(parameters,
72  *parameters.getCheckedPointerParam<FEProblemBase *>("_fe_problem_base"),
73  parameters.get<THREAD_ID>("_tid"),
74  true),
76  MooseVariableInterface<Real>(this,
77  true,
78  "variable",
81  TaggingInterface(this),
82  _subproblem(*getCheckedPointerParam<SubProblem *>("_subproblem")),
83  _fe_problem(*getCheckedPointerParam<FEProblemBase *>("_fe_problem_base")),
84  _sys(*getCheckedPointerParam<SystemBase *>("_sys")),
85  _tid(parameters.get<THREAD_ID>("_tid")),
86  _assembly(_subproblem.assembly(_tid)),
87  _var(*mooseVariable()),
88  _mesh(_subproblem.mesh()),
89  _current_node(_var.node()),
90  _u(_var.dofValues()),
91  _save_in_strings(parameters.get<std::vector<AuxVariableName>>("save_in")),
92  _diag_save_in_strings(parameters.get<std::vector<AuxVariableName>>("diag_save_in"))
93 
94 {
95  _save_in.resize(_save_in_strings.size());
96  _diag_save_in.resize(_diag_save_in_strings.size());
97 
98  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
99  {
101 
102  if (var->feType() != _var.feType())
103  paramError(
104  "save_in",
105  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
107 
108  _save_in[i] = var;
111  }
112 
113  _has_save_in = _save_in.size() > 0;
114 
115  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
116  {
118 
119  if (var->feType() != _var.feType())
120  paramError(
121  "diag_save_in",
122  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
124 
125  _diag_save_in[i] = var;
128  }
129 
130  _has_diag_save_in = _diag_save_in.size() > 0;
131 }
132 
135 {
136  return _var;
137 }
138 
139 SubProblem &
141 {
142  return _subproblem;
143 }
144 
145 void
147 {
148  if (_var.isNodalDefined())
149  {
150  const dof_id_type & dof_idx = _var.nodalDofIndex();
151  _qp = 0;
152  Real res = computeQpResidual();
154 
155  if (_has_save_in)
156  {
157  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
158  for (const auto & var : _save_in)
159  var->sys().solution().add(var->nodalDofIndex(), res);
160  }
161  }
162 }
163 
164 void
166 {
167  if (_var.isNodalDefined())
168  {
169  _qp = 0;
170  Real cached_val = computeQpJacobian();
171  dof_id_type cached_row = _var.nodalDofIndex();
172 
173  _assembly.cacheJacobianContribution(cached_row, cached_row, cached_val, _matrix_tags);
174 
175  if (_has_diag_save_in)
176  {
177  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
178  for (const auto & var : _diag_save_in)
179  var->sys().solution().add(var->nodalDofIndex(), cached_val);
180  }
181  }
182 }
183 
184 void
186 {
187  if (jvar == _var.number())
188  computeJacobian();
189  else
190  {
191  _qp = 0;
192  Real cached_val = computeQpOffDiagJacobian(jvar);
193  dof_id_type cached_row = _var.nodalDofIndex();
194  // Note: this only works for Lagrange variables...
195  dof_id_type cached_col = _current_node->dof_number(_sys.number(), jvar, 0);
196 
197  _assembly.cacheJacobianContribution(cached_row, cached_col, cached_val, _matrix_tags);
198  }
199 }
200 
201 Real
203 {
204  return 0.;
205 }
206 
207 Real
209 {
210  return 0.;
211 }
Interface for objects that need parallel consistent random numbers without patterns over the course o...
VarFieldType
Definition: MooseTypes.h:488
MooseVariable & variable()
Gets the variable this is active on.
Definition: NodalKernel.C:134
std::vector< AuxVariableName > _save_in_strings
Definition: NodalKernel.h:159
A class for creating restricted objects.
Definition: Restartable.h:29
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.
std::string incompatVarMsg(MooseVariableFEBase &var1, MooseVariableFEBase &var2)
Builds and returns a string of the form:
Definition: MooseError.C:22
InputParameters validParams< BlockRestrictable >()
unsigned int number() const
Get variable number coming from libMesh.
SubProblem & subProblem()
Get a reference to the subproblem.
Definition: NodalKernel.C:140
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:208
Assembly & _assembly
Reference to assembly.
Definition: NodalKernel.h:139
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()
Compute the Jacobian at one node.
Definition: NodalKernel.C:165
virtual void computeResidual()
Compute the residual at the current node.
Definition: NodalKernel.C:146
Base class for a system (of equations)
Definition: SystemBase.h:92
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
SubProblem & _subproblem
Reference to SubProblem.
Definition: NodalKernel.h:127
const FEType & feType() const
Get the type of finite element object.
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...
void addMooseVariableDependency(MooseVariableFEBase *var)
Call this function to add the passed in MooseVariableFEBase as a variable that this object depends on...
MooseVariable & _var
variable this works on
Definition: NodalKernel.h:142
Interface for objects that needs transient capabilities.
std::vector< AuxVariableName > _diag_save_in_strings
Definition: NodalKernel.h:164
virtual void computeOffDiagJacobian(unsigned int jvar)
Compute the off-diagonal Jacobian at one node.
Definition: NodalKernel.C:185
Interface for notifications that the mesh has changed.
InputParameters validParams< RandomInterface >()
Every object that can be built by the factory should be derived from this class.
Definition: MooseObject.h:42
void cacheJacobianContribution(numeric_index_type i, numeric_index_type j, Real value, TagID tag=0)
Caches the Jacobian entry &#39;value&#39;, to eventually be added/set in the (i,j) location of the matrix...
Definition: Assembly.C:3430
std::set< TagID > _matrix_tags
The matrices this Kernel will contribute to.
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
virtual Real computeQpJacobian()
The user can override this function to compute the "on-diagonal" Jacobian contribution.
Definition: NodalKernel.C:202
std::set< TagID > _vector_tags
The vectors this Kernel will contribute to.
Interface for objects that need to use UserObjects.
virtual unsigned int number() const
Gets the number of this system.
Definition: SystemBase.C:926
SystemBase & _sys
Reference to SystemBase.
Definition: NodalKernel.h:133
InputParameters validParams< MooseObject >()
Definition: MooseObject.C:25
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: NodalKernel.h:162
const Node *const & _current_node
current node being processed
Definition: NodalKernel.h:148
void paramError(const std::string &param, Args... args)
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
Definition: MooseObject.h:108
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: NodalKernel.h:157
bool isNodalDefined() const
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:59
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:178
InputParameters validParams< NodalKernel >()
Definition: NodalKernel.C:19
InputParameters validParams< BoundaryRestrictable >()
InputParameters validParams< TaggingInterface >()
Intermediate base class that ties together all the interfaces for getting MooseVariableFEBases with t...
Interface for objects that need to get values of MooseVariables.
An interface that restricts an object to subdomains via the &#39;blocks&#39; input parameter.
unsigned int _qp
Quadrature point index.
Definition: NodalKernel.h:151
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
void cacheResidualContribution(dof_id_type dof, Real value, TagID tag_id)
Cache individual residual contributions.
Definition: Assembly.C:2825
std::vector< MooseVariableFEBase * > _save_in
Definition: NodalKernel.h:158
Definition: Moose.h:112
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: NodalKernel.h:163
SystemBase & sys()
Get the system this variable is part of.
Interface for objects that need to use functions.
InputParameters validParams< TransientInterface >()
THREAD_ID _tid
Thread id.
Definition: NodalKernel.h:136
virtual Real computeQpResidual()=0
The user can override this function to compute the residual at a node.
Interface class for classes which interact with Postprocessors.
NodalKernel(const InputParameters &parameters)
Class constructor.
Definition: NodalKernel.C:59
unsigned int THREAD_ID
Definition: MooseTypes.h:161