www.mooseframework.org
ADIntegratedBC.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 "ADIntegratedBC.h"
11 
12 // MOOSE includes
13 #include "Assembly.h"
14 #include "SubProblem.h"
15 #include "SystemBase.h"
16 #include "MooseVariableFE.h"
17 #include "MooseVariableScalar.h"
18 
19 #include "libmesh/quadrature.h"
20 
23 
24 template <typename T, ComputeStage compute_stage>
26  : IntegratedBCBase(parameters),
27  MooseVariableInterface<T>(this,
28  false,
29  "variable",
31  std::is_same<T, Real>::value ? Moose::VarFieldType::VAR_FIELD_STANDARD
33  _var(*this->mooseVariable()),
34  _normals(_assembly.adNormals<compute_stage>()),
35  _ad_q_points(_assembly.adQPointsFace<compute_stage>()),
36  _test(_var.phiFace()),
37  _grad_test(_var.template adGradPhiFace<compute_stage>()),
38  _u(_var.template adSln<compute_stage>()),
39  _grad_u(_var.template adGradSln<compute_stage>()),
40  _ad_JxW(_assembly.adJxWFace<compute_stage>()),
41  _ad_coord(_assembly.template adCoordTransformation<compute_stage>())
42 {
44 
45  _save_in.resize(_save_in_strings.size());
46  _diag_save_in.resize(_diag_save_in_strings.size());
47 
48  for (unsigned int i = 0; i < _save_in_strings.size(); i++)
49  {
51 
52  if (var->feType() != _var.feType())
53  paramError(
54  "save_in",
55  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
57  _save_in[i] = var;
60  }
61 
62  _has_save_in = _save_in.size() > 0;
63 
64  for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
65  {
67 
68  if (var->feType() != _var.feType())
69  paramError(
70  "diag_save_in",
71  "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
73 
74  _diag_save_in[i] = var;
77  }
78 
79  _has_diag_save_in = _diag_save_in.size() > 0;
80 }
81 
82 template <typename T, ComputeStage compute_stage>
83 void
85 {
86  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
87  _local_re.resize(re.size());
88  _local_re.zero();
89 
90  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
91  for (_i = 0; _i < _test.size(); _i++)
92  _local_re(_i) += _ad_JxW[_qp] * _ad_coord[_qp] * computeQpResidual();
93 
94  re += _local_re;
95 
96  if (_has_save_in)
97  {
98  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
99  for (unsigned int i = 0; i < _save_in.size(); i++)
100  _save_in[i]->sys().solution().add_vector(_local_re, _save_in[i]->dofIndices());
101  }
102 }
103 
104 template <>
105 void
107 {
108 }
109 
110 template <>
111 void
113 {
114 }
115 
116 template <typename T, ComputeStage compute_stage>
117 void
119 {
120  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
121  _local_ke.resize(ke.m(), ke.n());
122  _local_ke.zero();
123 
124  size_t ad_offset = _var.number() * _sys.getMaxVarNDofsPerElem();
125 
126  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
127  for (_i = 0; _i < _test.size(); _i++)
128  {
129  DualReal residual = computeQpResidual();
130  for (_j = 0; _j < _var.phiSize(); ++_j)
131  _local_ke(_i, _j) +=
132  (_ad_JxW[_qp] * _ad_coord[_qp] * residual).derivatives()[ad_offset + _j];
133  }
134 
135  ke += _local_ke;
136 
137  if (_has_diag_save_in)
138  {
139  unsigned int rows = ke.m();
140  DenseVector<Number> diag(rows);
141  for (unsigned int i = 0; i < rows; i++)
142  diag(i) = _local_ke(i, i);
143 
144  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
145  for (unsigned int i = 0; i < _diag_save_in.size(); i++)
146  _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
147  }
148 }
149 
150 template <>
151 void
153 {
154 }
155 template <>
156 void
158 {
159 }
160 
161 template <typename T, ComputeStage compute_stage>
162 void
164 {
165  auto jvar_num = jvar.number();
166 
167  if (jvar_num == _var.number())
168  computeJacobian();
169  else
170  {
171  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar_num);
172  if (jvar.phiFaceSize() != ke.n())
173  return;
174 
175  size_t ad_offset = jvar_num * _sys.getMaxVarNDofsPerElem();
176 
177  for (_i = 0; _i < _test.size(); _i++)
178  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
179  {
180  DualReal residual = _ad_JxW[_qp] * _coord[_qp] * computeQpResidual();
181 
182  for (_j = 0; _j < jvar.phiFaceSize(); _j++)
183  ke(_i, _j) += residual.derivatives()[ad_offset + _j];
184  }
185  }
186 }
187 
188 template <>
189 void
191 {
192 }
193 template <>
194 void
196 {
197 }
198 
199 template <typename T, ComputeStage compute_stage>
200 void
202 {
203 }
204 
VarFieldType
Definition: MooseTypes.h:488
SubProblem & _subproblem
Reference to SubProblem.
std::vector< MooseVariableFEBase * > _diag_save_in
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
std::vector< MooseVariableFEBase * > _save_in
unsigned int number() const
Get variable number coming from libMesh.
void computeJacobianBlock(MooseVariableFEBase &jvar) override
Computes d-ivar-residual / d-jvar...
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
DualNumber< Real, NumberArray< AD_MAX_DOFS_PER_ELEM, Real > > DualReal
Definition: DualReal.h:29
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...
MooseVariableFE< T > * mooseVariable() const
Get the variable that this object is using.
std::vector< AuxVariableName > _save_in_strings
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
Base class for deriving any boundary condition of a integrated type.
void computeResidual() override
bool _has_save_in
The aux variables to save the residual contributions to.
virtual size_t phiFaceSize() const =0
Return phiFace size.
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
std::vector< AuxVariableName > _diag_save_in_strings
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
Interface for objects that need to get values of MooseVariables.
Base class for deriving any boundary condition of a integrated type.
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
MooseVariableFE< T > & _var
The variable that this IntegratedBC operates on.
Definition: Moose.h:112
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
SystemBase & sys()
Get the system this variable is part of.
defineADBaseValidParams(ADIntegratedBC, IntegratedBCBase,)
void computeJacobian() override
THREAD_ID _tid
Thread id.
void computeJacobianBlockScalar(unsigned int jvar) override
Computes jacobian block with respect to a scalar variable.
ADIntegratedBCTempl(const InputParameters &parameters)