www.mooseframework.org
NodalBC.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 "NodalBC.h"
16 
17 #include "Assembly.h"
18 #include "MooseVariable.h"
19 #include "SystemBase.h"
20 
21 template <>
24 {
26  params += validParams<RandomInterface>();
27  params.addParam<std::vector<AuxVariableName>>(
28  "save_in",
29  "The name of auxiliary variables to save this BC's residual contributions to. "
30  "Everything about that variable must match everything about this variable (the "
31  "type, what blocks it's on, etc.)");
32  params.addParam<std::vector<AuxVariableName>>(
33  "diag_save_in",
34  "The name of auxiliary variables to save this BC's diagonal jacobian "
35  "contributions to. Everything about that variable must match everything "
36  "about this variable (the type, what blocks it's on, etc.)");
37 
38  return params;
39 }
40 
41 NodalBC::NodalBC(const InputParameters & parameters)
42  : BoundaryCondition(parameters, true), // true is for being Nodal
43  RandomInterface(parameters, _fe_problem, _tid, true),
45  _current_node(_var.node()),
46  _u(_var.nodalSln()),
47  _save_in_strings(parameters.get<std::vector<AuxVariableName>>("save_in")),
48  _diag_save_in_strings(parameters.get<std::vector<AuxVariableName>>("diag_save_in")),
49  _is_eigen(false)
50 {
51  _save_in.resize(_save_in_strings.size());
52  _diag_save_in.resize(_diag_save_in_strings.size());
53 
54  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
55  {
57 
58  if (var->feType() != _var.feType())
59  mooseError("Error in " + name() + ". When saving residual values in an Auxiliary variable "
60  "the AuxVariable must be the same type as the nonlinear "
61  "variable the object is acting on.");
62 
63  _save_in[i] = var;
66  }
67 
68  _has_save_in = _save_in.size() > 0;
69 
70  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
71  {
73 
74  if (var->feType() != _var.feType())
75  mooseError("Error in " + name() + ". When saving diagonal Jacobian values in an Auxiliary "
76  "variable the AuxVariable must be the same type as the "
77  "nonlinear variable the object is acting on.");
78 
79  _diag_save_in[i] = var;
82  }
83 
84  _has_diag_save_in = _diag_save_in.size() > 0;
85 }
86 
87 void
88 NodalBC::computeResidual(NumericVector<Number> & residual)
89 {
90  if (_var.isNodalDefined())
91  {
92  dof_id_type & dof_idx = _var.nodalDofIndex();
93  _qp = 0;
94  Real res = 0;
95 
96  if (!_is_eigen)
97  res = computeQpResidual();
98 
99  residual.set(dof_idx, res);
100 
101  if (_has_save_in)
102  {
103  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
104  for (unsigned int i = 0; i < _save_in.size(); i++)
105  _save_in[i]->sys().solution().set(_save_in[i]->nodalDofIndex(), res);
106  }
107  }
108 }
109 
110 void
112 {
113  // We call the user's computeQpJacobian() function and store the
114  // results in the _assembly object. We can't store them directly in
115  // the element stiffness matrix, as they will only be inserted after
116  // all the assembly is done.
117  if (_var.isNodalDefined())
118  {
119  _qp = 0;
120  Real cached_val = computeQpJacobian();
121  dof_id_type cached_row = _var.nodalDofIndex();
122 
123  // Cache the user's computeQpJacobian() value for later use.
124  _fe_problem.assembly(0).cacheJacobianContribution(cached_row, cached_row, cached_val);
125 
126  if (_has_diag_save_in)
127  {
128  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
129  for (unsigned int i = 0; i < _diag_save_in.size(); i++)
130  _diag_save_in[i]->sys().solution().set(_diag_save_in[i]->nodalDofIndex(), cached_val);
131  }
132  }
133 }
134 
135 void
137 {
138  if (jvar == _var.number())
139  computeJacobian();
140  else
141  {
142  _qp = 0;
143  Real cached_val = computeQpOffDiagJacobian(jvar);
144  dof_id_type cached_row = _var.nodalDofIndex();
145  // Note: this only works for Lagrange variables...
146  dof_id_type cached_col = _current_node->dof_number(_sys.number(), jvar, 0);
147 
148  // Cache the user's computeQpJacobian() value for later use.
149  _fe_problem.assembly(0).cacheJacobianContribution(cached_row, cached_col, cached_val);
150  }
151 }
152 
153 Real
155 {
156  return 1.;
157 }
158 
159 Real
160 NodalBC::computeQpOffDiagJacobian(unsigned int /*jvar*/)
161 {
162  return 0.;
163 }
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...
SubProblem & _subproblem
Reference to SubProblem.
const FEType & feType() const
Get the type of finite element object.
Class for stuff related to variables.
Definition: MooseVariable.h:43
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. ...
FEProblemBase & _fe_problem
Reference to FEProblemBase.
virtual Assembly & assembly(THREAD_ID tid) override
virtual Real computeQpResidual()=0
std::vector< AuxVariableName > _save_in_strings
Definition: NodalBC.h:63
NodalBC(const InputParameters &parameters)
Definition: NodalBC.C:41
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...
dof_id_type & nodalDofIndex()
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: NodalBC.h:66
InputParameters validParams< NodalBC >()
Definition: NodalBC.C:23
MooseVariable & _var
variable this BC works on
virtual Real computeQpJacobian()
The user can override this function to compute the "on-diagonal" Jacobian contribution for this Nodal...
Definition: NodalBC.C:154
InputParameters validParams< RandomInterface >()
virtual void computeJacobian()
Definition: NodalBC.C:111
virtual void computeOffDiagJacobian(unsigned int jvar)
Definition: NodalBC.C:136
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
std::vector< AuxVariableName > _diag_save_in_strings
Definition: NodalBC.h:68
bool _is_eigen
Indicate whether or not the boundary condition is applied to the right hand side of eigenvalue proble...
Definition: NodalBC.h:72
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:160
Base class for creating new types of boundary conditions.
std::vector< MooseVariable * > _save_in
Definition: NodalBC.h:62
virtual unsigned int number()
Gets the number of this system.
Definition: SystemBase.C:604
unsigned int number() const
Get variable number coming from libMesh.
bool _has_save_in
The aux variables to save the residual contributions to.
Definition: NodalBC.h:61
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
SystemBase & _sys
Reference to SystemBase.
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
const Node *& _current_node
current node being processed
Definition: NodalBC.h:53
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 addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
std::vector< MooseVariable * > _diag_save_in
Definition: NodalBC.h:67
SystemBase & sys()
Get the system this variable is part of.
unsigned int _qp
Quadrature point index.
Definition: NodalBC.h:56
THREAD_ID _tid
Thread id.
InputParameters validParams< BoundaryCondition >()
virtual void computeResidual(NumericVector< Number > &residual)
Definition: NodalBC.C:88
bool isNodalDefined()