www.mooseframework.org
NodalKernel.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 "NodalKernel.h"
16 #include "Problem.h"
17 #include "SubProblem.h"
18 #include "SystemBase.h"
19 #include "MooseVariable.h"
20 #include "Assembly.h"
21 
22 template <>
25 {
30  params += validParams<RandomInterface>();
31 
32  params.addRequiredParam<NonlinearVariableName>(
33  "variable", "The name of the variable that this boundary condition applies to");
34 
35  params.addParam<std::vector<AuxVariableName>>(
36  "save_in",
37  "The name of auxiliary variables to save this BC's residual contributions to. "
38  "Everything about that variable must match everything about this variable (the "
39  "type, what blocks it's on, etc.)");
40 
41  params.addParam<std::vector<AuxVariableName>>(
42  "diag_save_in",
43  "The name of auxiliary variables to save this BC's diagonal jacobian "
44  "contributions to. Everything about that variable must match everything "
45  "about this variable (the type, what blocks it's on, etc.)");
46 
47  params.addParam<bool>("use_displaced_mesh",
48  false,
49  "Whether or not this object should use the "
50  "displaced mesh for computation. Note that "
51  "in the case this is true but no "
52  "displacements are provided in the Mesh block "
53  "the undisplaced mesh will still be used.");
54  params.addParamNamesToGroup("use_displaced_mesh", "Advanced");
55 
56  params.declareControllable("enable");
57 
58  params.registerBase("NodalKernel");
59 
60  return params;
61 }
62 
64  : MooseObject(parameters),
65  BlockRestrictable(this),
66  BoundaryRestrictable(this, true), // true for applying to nodesets
67  SetupInterface(this),
68  FunctionInterface(this),
69  UserObjectInterface(this),
70  TransientInterface(this),
73  Restartable(parameters, "BCs"),
74  ZeroInterface(parameters),
75  MeshChangedInterface(parameters),
76  RandomInterface(parameters,
77  *parameters.get<FEProblemBase *>("_fe_problem_base"),
78  parameters.get<THREAD_ID>("_tid"),
79  true),
81  _subproblem(*parameters.get<SubProblem *>("_subproblem")),
82  _fe_problem(*parameters.get<FEProblemBase *>("_fe_problem_base")),
83  _sys(*parameters.get<SystemBase *>("_sys")),
84  _tid(parameters.get<THREAD_ID>("_tid")),
85  _assembly(_subproblem.assembly(_tid)),
86  _var(_sys.getVariable(_tid, parameters.get<NonlinearVariableName>("variable"))),
87  _mesh(_subproblem.mesh()),
88  _current_node(_var.node()),
89  _u(_var.nodalSln()),
90  _u_dot(_var.nodalSlnDot()),
91  _du_dot_du(_var.nodalSlnDuDotDu()),
92  _save_in_strings(parameters.get<std::vector<AuxVariableName>>("save_in")),
93  _diag_save_in_strings(parameters.get<std::vector<AuxVariableName>>("diag_save_in"))
94 
95 {
96  _save_in.resize(_save_in_strings.size());
97  _diag_save_in.resize(_diag_save_in_strings.size());
98 
99  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
100  {
102 
103  if (var->feType() != _var.feType())
104  mooseError("Error in " + name() + ". When saving residual values in an Auxiliary variable "
105  "the AuxVariable must be the same type as the nonlinear "
106  "variable the object is acting on.");
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  mooseError("Error in " + name() + ". When saving diagonal Jacobian values in an Auxiliary "
121  "variable the AuxVariable must be the same type as the "
122  "nonlinear variable the object is acting on.");
123 
124  _diag_save_in[i] = var;
127  }
128 
129  _has_diag_save_in = _diag_save_in.size() > 0;
130 }
131 
134 {
135  return _var;
136 }
137 
138 SubProblem &
140 {
141  return _subproblem;
142 }
143 
144 void
146 {
147  if (_var.isNodalDefined())
148  {
149  dof_id_type & dof_idx = _var.nodalDofIndex();
150  _qp = 0;
151  Real res = computeQpResidual();
153 
154  if (_has_save_in)
155  {
156  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
157  for (const auto & var : _save_in)
158  var->sys().solution().add(var->nodalDofIndex(), res);
159  }
160  }
161 }
162 
163 void
165 {
166  if (_var.isNodalDefined())
167  {
168  _qp = 0;
169  Real cached_val = computeQpJacobian();
170  dof_id_type cached_row = _var.nodalDofIndex();
171 
172  _assembly.cacheJacobianContribution(cached_row, cached_row, cached_val);
173 
174  if (_has_diag_save_in)
175  {
176  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
177  for (const auto & var : _diag_save_in)
178  var->sys().solution().add(var->nodalDofIndex(), cached_val);
179  }
180  }
181 }
182 
183 void
185 {
186  if (jvar == _var.number())
187  computeJacobian();
188  else
189  {
190  _qp = 0;
191  Real cached_val = computeQpOffDiagJacobian(jvar);
192  dof_id_type cached_row = _var.nodalDofIndex();
193  // Note: this only works for Lagrange variables...
194  dof_id_type cached_col = _current_node->dof_number(_sys.number(), jvar, 0);
195 
196  _assembly.cacheJacobianContribution(cached_row, cached_col, cached_val);
197  }
198 }
199 
200 Real
202 {
203  return 0.;
204 }
205 
206 Real
208 {
209  return 0.;
210 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
Interface for objects that need parallel consistent random numbers without patterns over the course o...
MooseVariable & variable()
Gets the variable this is active on.
Definition: NodalKernel.C:133
std::vector< AuxVariableName > _save_in_strings
Definition: NodalKernel.h:166
A class for creating restricted objects.
Definition: Restartable.h:31
const FEType & feType() const
Get the type of finite element object.
Class for stuff related to variables.
Definition: MooseVariable.h:43
InputParameters validParams< BlockRestrictable >()
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. ...
void cacheResidualContribution(dof_id_type dof, Real value, Moose::KernelType type)
Cache individual residual contributions.
Definition: Assembly.C:1375
SubProblem & subProblem()
Get a reference to the subproblem.
Definition: NodalKernel.C:139
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:207
const Node *& _current_node
current node being processed
Definition: NodalKernel.h:149
std::vector< MooseVariable * > _diag_save_in
Definition: NodalKernel.h:170
Assembly & _assembly
Reference to assembly.
Definition: NodalKernel.h:140
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
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 void computeJacobian()
Compute the Jacobian at one node.
Definition: NodalKernel.C:164
virtual void computeResidual()
Compute the residual at the current node.
Definition: NodalKernel.C:145
Base class for a system (of equations)
Definition: SystemBase.h:91
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
SubProblem & _subproblem
Reference to SubProblem.
Definition: NodalKernel.h:128
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...
MooseVariable & _var
variable this works on
Definition: NodalKernel.h:143
dof_id_type & nodalDofIndex()
Interface for objects that needs transient capabilities.
std::vector< AuxVariableName > _diag_save_in_strings
Definition: NodalKernel.h:171
virtual void computeOffDiagJacobian(unsigned int jvar)
Compute the off-diagonal Jacobian at one node.
Definition: NodalKernel.C:184
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:36
virtual Real computeQpJacobian()
The user can override this function to compute the "on-diagonal" Jacobian contribution.
Definition: NodalKernel.C:201
Interface for objects that need to use UserObjects.
SystemBase & _sys
Reference to SystemBase.
Definition: NodalKernel.h:134
void cacheJacobianContribution(numeric_index_type i, numeric_index_type j, Real value)
Caches the Jacobian entry &#39;value&#39;, to eventually be added/set in the (i,j) location of the matrix...
Definition: Assembly.C:1903
InputParameters validParams< MooseObject >()
Definition: MooseObject.C:22
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: NodalKernel.h:169
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: NodalKernel.h:164
virtual unsigned int number()
Gets the number of this system.
Definition: SystemBase.C:604
unsigned int number() const
Get variable number coming from libMesh.
std::vector< MooseVariable * > _save_in
Definition: NodalKernel.h:165
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:53
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< NodalKernel >()
Definition: NodalKernel.C:24
InputParameters validParams< BoundaryRestrictable >()
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
An interface that restricts an object to subdomains via the &#39;blocks&#39; input parameter.
unsigned int _qp
Quadrature point index.
Definition: NodalKernel.h:152
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
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
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:137
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:63
unsigned int THREAD_ID
Definition: MooseTypes.h:79
bool isNodalDefined()