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 
17 template <>
20 {
22  params += validParams<RandomInterface>();
23 
24  return params;
25 }
26 
27 NodalBC::NodalBC(const InputParameters & parameters)
28  : NodalBCBase(parameters),
29  MooseVariableInterface<Real>(this,
30  true,
31  "variable",
34  _var(*mooseVariable()),
35  _current_node(_var.node()),
36  _u(_var.dofValues())
37 {
39 
40  _save_in.resize(_save_in_strings.size());
41  _diag_save_in.resize(_diag_save_in_strings.size());
42 
43  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
44  {
46 
47  if (var->feType() != _var.feType())
48  paramError(
49  "save_in",
50  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
52 
53  _save_in[i] = var;
56  }
57 
58  _has_save_in = _save_in.size() > 0;
59 
60  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
61  {
63 
64  if (var->feType() != _var.feType())
65  paramError(
66  "diag_save_in",
67  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
69 
70  _diag_save_in[i] = var;
73  }
74 
75  _has_diag_save_in = _diag_save_in.size() > 0;
76 }
77 
78 void
80 {
81  if (_var.isNodalDefined())
82  {
83  const dof_id_type & dof_idx = _var.nodalDofIndex();
84  _qp = 0;
85  Real res = 0;
86  res = computeQpResidual();
87 
88  for (auto tag_id : _vector_tags)
89  if (_sys.hasVector(tag_id))
90  _sys.getVector(tag_id).set(dof_idx, res);
91 
92  if (_has_save_in)
93  {
94  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
95  for (unsigned int i = 0; i < _save_in.size(); i++)
96  _save_in[i]->sys().solution().set(_save_in[i]->nodalDofIndex(), res);
97  }
98  }
99 }
100 
101 void
103 {
104  // We call the user's computeQpJacobian() function and store the
105  // results in the _assembly object. We can't store them directly in
106  // the element stiffness matrix, as they will only be inserted after
107  // all the assembly is done.
108  if (_var.isNodalDefined())
109  {
110  _qp = 0;
111  Real cached_val = 0.;
112  cached_val = computeQpJacobian();
113 
114  dof_id_type cached_row = _var.nodalDofIndex();
115 
116  // Cache the user's computeQpJacobian() value for later use.
117  for (auto tag : _matrix_tags)
118  if (_sys.hasMatrix(tag))
119  _fe_problem.assembly(0).cacheJacobianContribution(cached_row, cached_row, cached_val, tag);
120 
121  if (_has_diag_save_in)
122  {
123  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
124  for (unsigned int i = 0; i < _diag_save_in.size(); i++)
125  _diag_save_in[i]->sys().solution().set(_diag_save_in[i]->nodalDofIndex(), cached_val);
126  }
127  }
128 }
129 
130 void
132 {
133  if (jvar == _var.number())
134  computeJacobian();
135  else
136  {
137  _qp = 0;
138  Real cached_val = 0.0;
139  cached_val = computeQpOffDiagJacobian(jvar);
140 
141  dof_id_type cached_row = _var.nodalDofIndex();
142  // Note: this only works for Lagrange variables...
143  dof_id_type cached_col = _current_node->dof_number(_sys.number(), jvar, 0);
144 
145  // Cache the user's computeQpJacobian() value for later use.
146  for (auto tag : _matrix_tags)
147  if (_sys.hasMatrix(tag))
148  _fe_problem.assembly(0).cacheJacobianContribution(cached_row, cached_col, cached_val, tag);
149  }
150 }
151 
152 Real
154 {
155  return 1.;
156 }
157 
158 Real
159 NodalBC::computeQpOffDiagJacobian(unsigned int /*jvar*/)
160 {
161  return 0.;
162 }
VarFieldType
Definition: MooseTypes.h:488
virtual bool hasMatrix(TagID tag) const
Check if the tagged matrix exists in the system.
Definition: SystemBase.C:805
virtual void computeOffDiagJacobian(unsigned int jvar) override
Definition: NodalBC.C:131
SubProblem & _subproblem
Reference to SubProblem.
bool hasVector(const std::string &tag_name) const
Check if the named vector exists in the system.
Definition: SystemBase.C:701
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
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: NodalBCBase.h:44
std::vector< AuxVariableName > _save_in_strings
Definition: NodalBCBase.h:46
FEProblemBase & _fe_problem
Reference to FEProblemBase.
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:44
Assembly & assembly(THREAD_ID tid) override
NodalBC(const InputParameters &parameters)
Definition: NodalBC.C:27
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
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...
virtual void computeResidual() override
Definition: NodalBC.C:79
MooseVariableFE< Real > * mooseVariable() const
Get the variable that this object is using.
MooseVariable & _var
Definition: NodalBC.h:41
InputParameters validParams< NodalBC >()
Definition: NodalBC.C:19
virtual Real computeQpJacobian()
The user can override this function to compute the "on-diagonal" Jacobian contribution for this Nodal...
Definition: NodalBC.C:153
InputParameters validParams< RandomInterface >()
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
std::vector< MooseVariableFEBase * > _save_in
Definition: NodalBCBase.h:45
std::set< TagID > _vector_tags
The vectors this Kernel will contribute to.
virtual unsigned int number() const
Gets the number of this system.
Definition: SystemBase.C:926
Base class for deriving any boundary condition that works at nodes.
Definition: NodalBCBase.h:32
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
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:159
bool isNodalDefined() const
std::vector< AuxVariableName > _diag_save_in_strings
Definition: NodalBCBase.h:51
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: NodalBCBase.h:49
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
SystemBase & _sys
Reference to SystemBase.
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:184
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: NodalBCBase.h:50
Definition: Moose.h:112
SystemBase & sys()
Get the system this variable is part of.
virtual void computeJacobian() override
Definition: NodalBC.C:102
virtual NumericVector< Number > & getVector(const std::string &name)
Get a raw NumericVector.
Definition: SystemBase.C:751
unsigned int _qp
Quadrature point index.
Definition: NodalBC.h:47
THREAD_ID _tid
Thread id.
InputParameters validParams< NodalBCBase >()
Definition: NodalBCBase.C:14