https://mooseframework.inl.gov
ADNodalBC.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
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 #include "NodalBCBase.h"
12 #include "ADDirichletBCBase.h"
13 
14 // MOOSE includes
15 #include "Assembly.h"
16 #include "SubProblem.h"
17 #include "SystemBase.h"
18 #include "MooseVariableFE.h"
19 #include "MooseVariableScalar.h"
20 #include "FEProblemBase.h"
21 
22 #include "libmesh/quadrature.h"
23 
24 template <typename T, typename Base>
27 {
28  return Base::validParams();
29 }
30 
31 template <>
34 {
36  // The below parameters are useful for vector Nodal BCs
37  params.addParam<bool>("set_x_comp", true, "Whether to set the x-component of the variable");
38  params.addParam<bool>("set_y_comp", true, "Whether to set the y-component of the variable");
39  params.addParam<bool>("set_z_comp", true, "Whether to set the z-component of the variable");
40  return params;
41 }
42 
43 template <>
46 {
48  // The below parameters are useful for vector Nodal BCs
49  params.addParam<bool>("set_x_comp", true, "Whether to set the x-component of the variable");
50  params.addParam<bool>("set_y_comp", true, "Whether to set the y-component of the variable");
51  params.addParam<bool>("set_z_comp", true, "Whether to set the z-component of the variable");
52  return params;
53 }
54 
55 template <typename T, typename Base>
57  : Base(parameters),
58  MooseVariableInterface<T>(this,
59  true,
60  "variable",
64  ADFunctorInterface(this),
65  _var(*this->mooseVariable()),
66  _current_node(_var.node()),
67  _u(_var.adNodalValue()),
68  _set_components(
69  {std::is_same<T, RealVectorValue>::value ? this->template getParam<bool>("set_x_comp")
70  : true,
71  std::is_same<T, RealVectorValue>::value ? this->template getParam<bool>("set_y_comp")
72  : true,
73  std::is_same<T, RealVectorValue>::value ? this->template getParam<bool>("set_z_comp")
74  : true}),
75  _undisplaced_assembly(_fe_problem.assembly(_tid, _sys.number()))
76 {
77  _subproblem.haveADObjects(true);
78 
79  addMooseVariableDependency(this->mooseVariable());
80 }
81 
82 namespace
83 {
84 const ADReal &
85 conversionHelper(const ADReal & value, const unsigned int)
86 {
87  return value;
88 }
89 
90 const ADReal &
91 conversionHelper(const libMesh::VectorValue<ADReal> & value, const unsigned int i)
92 {
93  return value(i);
94 }
95 }
96 
97 template <typename T, typename Base>
98 template <typename ADResidual>
99 void
100 ADNodalBCTempl<T, Base>::addResidual(const ADResidual & residual,
101  const std::vector<dof_id_type> & dof_indices)
102 {
103  mooseAssert(dof_indices.size() <= _set_components.size(),
104  "The number of dof indices must be less than the number of settable components");
105 
106  for (const auto i : index_range(dof_indices))
107  if (_set_components[i])
108  setResidual(_sys, raw_value(conversionHelper(residual, i)), dof_indices[i]);
109 }
110 
111 template <typename T, typename Base>
112 template <typename ADResidual>
113 void
114 ADNodalBCTempl<T, Base>::addJacobian(const ADResidual & residual,
115  const std::vector<dof_id_type> & dof_indices)
116 {
117  mooseAssert(dof_indices.size() <= _set_components.size(),
118  "The number of dof indices must be less than the number of settable components");
119 
120  for (const auto i : index_range(dof_indices))
121  if (_set_components[i])
122  // If we store into the displaced assembly for nodal bc objects the data never actually makes
123  // it into the global Jacobian
124  addJacobian(_undisplaced_assembly,
125  std::array<ADReal, 1>{{conversionHelper(residual, i)}},
126  std::array<dof_id_type, 1>{{dof_indices[i]}},
127  /*scaling_factor=*/1);
128 }
129 
130 template <typename T, typename Base>
131 void
133 {
134  const std::vector<dof_id_type> & dof_indices = _var.dofIndices();
135  if (dof_indices.empty())
136  return;
137 
138  const auto residual = computeQpResidual();
139 
140  addResidual(residual, dof_indices);
141 }
142 
143 template <typename T, typename Base>
144 void
146 {
147  const std::vector<dof_id_type> & dof_indices = _var.dofIndices();
148  if (dof_indices.empty())
149  return;
150 
151  const auto residual = computeQpResidual();
152 
153  addJacobian(residual, dof_indices);
154 }
155 
156 template <typename T, typename Base>
157 void
159 {
160  const std::vector<dof_id_type> & dof_indices = _var.dofIndices();
161  if (dof_indices.empty())
162  return;
163 
164  const auto residual = computeQpResidual();
165 
166  addResidual(residual, dof_indices);
167  addJacobian(residual, dof_indices);
168 }
169 
170 template <typename T, typename Base>
171 void
173 {
174  // Only need to do this once because AD does all the derivatives at once
175  if (jvar_num == _var.number())
176  computeJacobian();
177 }
178 
179 template <typename T, typename Base>
180 void
182 {
183  // scalar coupling will have been included in the all-at-once handling in computeOffDiagJacobian
184 }
185 
186 template class ADNodalBCTempl<Real, NodalBCBase>;
VarFieldType
Definition: MooseTypes.h:721
void addJacobian(const ADResidual &residual, const std::vector< dof_id_type > &dof_indices)
process the Jacobian into the global data structures
Definition: ADNodalBC.C:114
void computeResidualAndJacobian() override
Definition: ADNodalBC.C:158
Base class for deriving any automatic differentiation boundary condition of a integrated type...
Definition: ADNodalBC.h:21
void computeJacobian() override final
Definition: ADNodalBC.C:145
auto raw_value(const Eigen::Map< T > &in)
Definition: EigenADReal.h:73
ADNodalBCTempl(const InputParameters &parameters)
Definition: ADNodalBC.C:56
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
DualNumber< Real, DNDerivativeType, true > ADReal
Definition: ADRealForward.h:47
InputParameters validParams()
An interface for accessing Moose::Functors for systems that care about automatic differentiation, e.g.
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:714
void computeOffDiagJacobian(unsigned int jvar) override final
Definition: ADNodalBC.C:172
static InputParameters validParams()
Definition: ADNodalBC.C:26
static InputParameters validParams()
void computeResidual() override final
Definition: ADNodalBC.C:132
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void computeOffDiagJacobianScalar(unsigned int jvar) override final
Definition: ADNodalBC.C:181
Interface for objects that need to get values of MooseVariables.
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an optional parameter and a documentation string to the InputParameters object...
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
void addResidual(const ADResidual &residual, const std::vector< dof_id_type > &dof_indices)
process the residual into the global data structures
Definition: ADNodalBC.C:100
auto index_range(const T &sizable)
static InputParameters validParams()
Definition: NodalBCBase.C:13