www.mooseframework.org
ADNodalBC.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 "ADNodalBC.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  : NodalBCBase(parameters),
27  MooseVariableInterface<T>(this,
28  true,
29  "variable",
31  std::is_same<T, Real>::value ? Moose::VarFieldType::VAR_FIELD_STANDARD
33  _var(*this->mooseVariable()),
34  _current_node(_var.node()),
35  _u(_var.template adNodalValue<compute_stage>())
36 {
38 }
39 
40 template <typename T>
41 T &
42 conversionHelper(T & value, const unsigned int &)
43 {
44  return value;
45 }
46 
47 template <typename T>
48 T &
49 conversionHelper(libMesh::VectorValue<T> & value, const unsigned int & i)
50 {
51  return value(i);
52 }
53 
54 template <typename T, ComputeStage compute_stage>
55 void
57 {
58  const std::vector<dof_id_type> & dof_indices = _var.dofIndices();
59 
60  auto residual = computeQpResidual();
61 
62  for (auto tag_id : _vector_tags)
63  if (_sys.hasVector(tag_id))
64  for (size_t i = 0; i < dof_indices.size(); ++i)
65  _sys.getVector(tag_id).set(dof_indices[i], conversionHelper(residual, i));
66 }
67 
68 template <>
69 void
71 {
72 }
73 
74 template <>
75 void
77 {
78 }
79 
80 template <typename T, ComputeStage compute_stage>
81 void
83 {
84  auto ad_offset = _var.number() * _sys.getMaxVarNDofsPerNode();
85  auto residual = computeQpResidual();
86  const std::vector<dof_id_type> & cached_rows = _var.dofIndices();
87 
88  // Cache the user's computeQpJacobian() value for later use.
89  for (auto tag : _matrix_tags)
90  if (_sys.hasMatrix(tag))
91  for (size_t i = 0; i < cached_rows.size(); ++i)
92  _fe_problem.assembly(0).cacheJacobianContribution(
93  cached_rows[i],
94  cached_rows[i],
95  conversionHelper(residual, i).derivatives()[ad_offset + i],
96  tag);
97 }
98 
99 template <>
100 void
102 {
103 }
104 template <>
105 void
107 {
108 }
109 
110 template <typename T, ComputeStage compute_stage>
111 void
113 {
114  if (jvar == _var.number())
115  computeJacobian();
116  else
117  {
118  auto ad_offset = jvar * _sys.getMaxVarNDofsPerNode();
119  auto residual = computeQpResidual();
120  const std::vector<dof_id_type> & cached_rows = _var.dofIndices();
121 
122  // Note: this only works for Lagrange variables...
123  dof_id_type cached_col = _current_node->dof_number(_sys.number(), jvar, 0);
124 
125  // Cache the user's computeQpJacobian() value for later use.
126  for (auto tag : _matrix_tags)
127  if (_sys.hasMatrix(tag))
128  for (size_t i = 0; i < cached_rows.size(); ++i)
129  _fe_problem.assembly(0).cacheJacobianContribution(
130  cached_rows[i],
131  cached_col,
132  conversionHelper(residual, i).derivatives()[ad_offset + i],
133  tag);
134  }
135 }
136 
137 template <>
138 void
140 {
141 }
142 template <>
143 void
145 {
146 }
147 
148 template class ADNodalBCTempl<Real, RESIDUAL>;
149 template class ADNodalBCTempl<Real, JACOBIAN>;
VarFieldType
Definition: MooseTypes.h:488
void computeJacobian() override
Definition: ADNodalBC.C:82
void computeOffDiagJacobian(unsigned int jvar) override
Definition: ADNodalBC.C:112
Base class for deriving any automatic differentiation boundary condition of a integrated type...
Definition: ADNodalBC.h:19
ADNodalBCTempl(const InputParameters &parameters)
Definition: ADNodalBC.C:25
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
void addMooseVariableDependency(MooseVariableFEBase *var)
Call this function to add the passed in MooseVariableFEBase as a variable that this object depends on...
void computeResidual() override
Definition: ADNodalBC.C:56
MooseVariableFE< T > * mooseVariable() const
Get the variable that this object is using.
defineADBaseValidParams(ADNodalBC, NodalBCBase,)
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:481
Base class for deriving any boundary condition that works at nodes.
Definition: NodalBCBase.h:32
Interface for objects that need to get values of MooseVariables.
Definition: Moose.h:112
T & conversionHelper(T &value, const unsigned int &)
Definition: ADNodalBC.C:42