www.mooseframework.org
NodalBC.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 "NodalBC.h"
11 
12 #include "Assembly.h"
13 #include "MooseVariableFE.h"
14 #include "SystemBase.h"
15 #include "NonlinearSystemBase.h"
16 
19 {
21 
22  return params;
23 }
24 
25 NodalBC::NodalBC(const InputParameters & parameters)
26  : NodalBCBase(parameters),
28  true,
29  "variable",
32  _var(*mooseVariable()),
33  _current_node(_var.node()),
34  _u(_var.dofValues())
35 {
37 
38  _save_in.resize(_save_in_strings.size());
39  _diag_save_in.resize(_diag_save_in_strings.size());
40 
41  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
42  {
44 
45  if (var->feType() != _var.feType())
46  paramError(
47  "save_in",
48  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
50 
51  _save_in[i] = var;
54  }
55 
56  _has_save_in = _save_in.size() > 0;
57 
58  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
59  {
61 
62  if (var->feType() != _var.feType())
63  paramError(
64  "diag_save_in",
65  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
67 
68  _diag_save_in[i] = var;
71  }
72 
73  _has_diag_save_in = _diag_save_in.size() > 0;
74 }
75 
76 void
78 {
79  if (_var.isNodalDefined())
80  {
81  const Real res = computeQpResidual();
82  setResidual(_sys, res, _var);
83 
84  if (_has_save_in)
85  for (unsigned int i = 0; i < _save_in.size(); i++)
86  _save_in[i]->sys().solution().set(_save_in[i]->nodalDofIndex(), res);
87  }
88 }
89 
90 void
92 {
93  // We call the user's computeQpJacobian() function and store the
94  // results in the _assembly object. We can't store them directly in
95  // the element stiffness matrix, as they will only be inserted after
96  // all the assembly is done.
97  if (_var.isNodalDefined())
98  {
99  const Real cached_val = computeQpJacobian();
100  const dof_id_type cached_row = _var.nodalDofIndex();
101 
102  // Cache the user's computeQpJacobian() value for later use.
104  cached_val,
105  cached_row,
106  cached_row,
107  /*scaling_factor=*/1);
108 
109  if (_has_diag_save_in)
110  for (unsigned int i = 0; i < _diag_save_in.size(); i++)
111  _diag_save_in[i]->sys().solution().set(_diag_save_in[i]->nodalDofIndex(), cached_val);
112  }
113 }
114 
115 void
116 NodalBC::computeOffDiagJacobian(const unsigned int jvar_num)
117 {
118  if (jvar_num == _var.number())
119  computeJacobian();
120  else
121  {
122  if (!_var.isNodalDefined())
123  return;
124 
125  const Real cached_val = computeQpOffDiagJacobian(jvar_num);
126 
127  if (cached_val == 0.)
128  // there's no reason to cache this if it's zero, and it can even lead to new nonzero
129  // allocations
130  return;
131 
132  const dof_id_type cached_row = _var.nodalDofIndex();
133  // Note: this only works for Lagrange variables...
134  const dof_id_type cached_col = _current_node->dof_number(_sys.number(), jvar_num, 0);
135 
136  // Cache the user's computeQpJacobian() value for later use.
138  cached_val,
139  cached_row,
140  cached_col,
141  /*scaling_factor=*/1);
142  }
143 }
144 
145 Real
147 {
148  return 1.;
149 }
150 
151 Real
152 NodalBC::computeQpOffDiagJacobian(unsigned int /*jvar*/)
153 {
154  return 0.;
155 }
156 
157 void
159 {
160  computeResidual();
161 
162  for (const auto & [ivariable, jvariable] : _fe_problem.couplingEntries(_tid, _sys.number()))
163  {
164  const unsigned int ivar = ivariable->number();
165  const unsigned int jvar = jvariable->number();
166 
167  if (ivar != _var.number())
168  continue;
169 
170  if (_is_implicit)
172  }
173 
175 }
VarFieldType
Definition: MooseTypes.h:634
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes this object&#39;s contribution to off-diagonal blocks of the system Jacobian matrix...
Definition: NodalBC.C:116
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase * > > & couplingEntries(const THREAD_ID tid, const unsigned int nl_sys_num)
Class for stuff related to variables.
Definition: Adaptivity.h:31
std::string incompatVarMsg(MooseVariableFieldBase &var1, MooseVariableFieldBase &var2)
Builds and returns a string of the form:
Definition: MooseError.C:23
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: NodalBCBase.h:41
virtual void computeResidualAndJacobian() override
Compute this object&#39;s contribution to the residual and Jacobian simultaneously.
Definition: NodalBC.C:158
std::vector< AuxVariableName > _save_in_strings
Definition: NodalBCBase.h:43
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpResidual()=0
const Node *const & _current_node
current node being processed
Definition: NodalBC.h:41
NodalBC(const InputParameters &parameters)
Definition: NodalBC.C:25
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
THREAD_ID _tid
The thread ID for this kernel.
const FEType & feType() const
Get the type of finite element object.
virtual void computeResidual() override
Compute this object&#39;s contribution to the residual.
Definition: NodalBC.C:77
SystemBase & _sys
Reference to the EquationSystem object.
MooseVariableFE< Real > * mooseVariable() const
MooseVariable & _var
Definition: NodalBC.h:38
virtual Real computeQpJacobian()
The user can override this function to compute the "on-diagonal" Jacobian contribution for this Nodal...
Definition: NodalBC.C:146
SubProblem & _subproblem
Reference to this kernel&#39;s SubProblem.
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:627
FEProblemBase & _fe_problem
Reference to this kernel&#39;s FEProblemBase.
std::vector< MooseVariableFEBase * > _save_in
Definition: NodalBCBase.h:42
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
Base class for deriving any boundary condition that works at nodes.
Definition: NodalBCBase.h:26
virtual Real computeQpOffDiagJacobian(unsigned int jvar)
This is the virtual that derived classes should override for computing an off-diagonal jacobian compo...
Definition: NodalBC.C:152
std::vector< AuxVariableName > _diag_save_in_strings
Definition: NodalBCBase.h:48
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: NodalBCBase.h:46
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
Interface for objects that need to get values of MooseVariables.
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 * > _diag_save_in
Definition: NodalBCBase.h:47
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
static InputParameters validParams()
Definition: NodalBC.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.
Assembly & assembly(const THREAD_ID tid, const unsigned int nl_sys_num) override
void setResidual(SystemBase &sys, const T &residual, MooseVariableFE< T > &var)
Set residual using the variables&#39; insertion API.
SystemBase & sys()
Get the system this variable is part of.
virtual void computeJacobian() override
Compute this object&#39;s contribution to the diagonal Jacobian entries.
Definition: NodalBC.C:91
bool _is_implicit
If the object is using implicit or explicit form.
uint8_t dof_id_type
static InputParameters validParams()
Definition: NodalBCBase.C:13