www.mooseframework.org
Function.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 "Function.h"
11 
12 using namespace Moose;
13 
16 {
18  params += SetupInterface::validParams();
19 
20  params.registerBase("Function");
21 
22  return params;
23 }
24 
26  : MooseObject(parameters),
27  SetupInterface(this),
28  TransientInterface(this),
30  UserObjectInterface(this),
31  Restartable(this, "Functions"),
32  MeshChangedInterface(parameters),
33  ScalarCoupleable(this),
35 {
36 }
37 
39 
40 Real
41 Function::value(Real /*t*/, const Point & /*p*/) const
42 {
43  mooseError("value method not implemented");
44  return 0.0;
45 }
46 
47 ADReal
48 Function::value(const ADReal & t, const ADPoint & p) const
49 {
50  const auto rt = MetaPhysicL::raw_value(t);
51  const auto rp = MetaPhysicL::raw_value(p);
52  const auto grad = gradient(rt, rp);
53  ADReal ret = value(rt, rp);
54  ret.derivatives() = grad(0) * p(0).derivatives()
55 #if LIBMESH_DIM > 1
56  + grad(1) * p(1).derivatives()
57 #endif
58 #if LIBMESH_DIM > 2
59  + grad(2) * p(2).derivatives()
60 #endif
61  + timeDerivative(rt, rp) * t.derivatives();
62  return ret;
63 }
64 
66 Function::value(const ChainedReal & t) const
67 {
68  static const Point p;
69  return {value(t.value(), p), timeDerivative(t.value(), p) * t.derivatives()};
70 }
71 
73 Function::gradient(Real /*t*/, const Point & /*p*/) const
74 {
75  mooseError("gradient method not implemented");
76  return RealGradient(0, 0, 0);
77 }
78 
79 Real
80 Function::timeDerivative(Real /*t*/, const Point & /*p*/) const
81 {
82  mooseError("timeDerivative method not implemented");
83  return 0;
84 }
85 
87 Function::vectorValue(Real /*t*/, const Point & /*p*/) const
88 {
89  mooseError("vectorValue method not implemented");
90  return RealVectorValue(0, 0, 0);
91 }
92 
94 Function::curl(Real /*t*/, const Point & /*p*/) const
95 {
96  mooseError("curl method not implemented");
97  return RealVectorValue(0, 0, 0);
98 }
99 
100 Real
101 Function::div(Real /*t*/, const Point & /*p*/) const
102 {
103  mooseError("div method not implemented");
104  return 0;
105 }
106 
107 Real
109 {
110  mooseError("Integral method not implemented for function ", name());
111  return 0;
112 }
113 
114 Real
116 {
117  mooseError("Average method not implemented for function ", name());
118  return 0;
119 }
120 
121 template <typename R>
122 typename Function::ValueType
123 Function::evaluateHelper(const R & r, const Moose::StateArg & state) const
124 {
125  return value(_ti_feproblem.getTimeFromStateArg(state), r.getPoint());
126 }
127 
128 typename Function::ValueType
129 Function::evaluate(const ElemArg & elem_arg, const Moose::StateArg & state) const
130 {
131  return evaluateHelper(elem_arg, state);
132 }
133 
134 typename Function::ValueType
135 Function::evaluate(const FaceArg & face, const Moose::StateArg & state) const
136 {
137  if (face.face_side && face.fi->neighborPtr() &&
138  (face.fi->elem().subdomain_id() != face.fi->neighbor().subdomain_id()))
139  {
140  // Some users like to put discontinuities in their functions at subdomain changes in which case
141  // in order to always get the proper discontinuous effect we should evaluate ever so slightly
142  // off the face. Consider evaluation of: if(x < 0, -1, 1) if the face centroid is right at x ==
143  // 0 for example. The user likely doesn't want you to return 1 if they've asked for a 0-
144  // evaluation
145  //
146  // I can't quite tell but I think the tolerance for comparing coordinates (x, y, z, t) in
147  // fparser is ~1e-9 so we need to use something larger than that. The comparison is absolute
148  static constexpr Real offset_tolerance = 1e-8;
149  auto offset = offset_tolerance * face.fi->normal();
150  if (face.face_side == face.fi->elemPtr())
151  offset *= -1;
152  return value(_ti_feproblem.getTimeFromStateArg(state), face.getPoint() + offset);
153  }
154  else
155  return value(_ti_feproblem.getTimeFromStateArg(state), face.getPoint());
156 }
157 
158 typename Function::ValueType
159 Function::evaluate(const ElemQpArg & elem_qp, const Moose::StateArg & state) const
160 {
161  return evaluateHelper(elem_qp, state);
162 }
163 
164 typename Function::ValueType
165 Function::evaluate(const ElemSideQpArg & elem_side_qp, const Moose::StateArg & state) const
166 {
167  return evaluateHelper(elem_side_qp, state);
168 }
169 
170 typename Function::ValueType
171 Function::evaluate(const ElemPointArg & elem_point_arg, const Moose::StateArg & state) const
172 {
173  return evaluateHelper(elem_point_arg, state);
174 }
175 
176 typename Function::ValueType
177 Function::evaluate(const NodeArg & node_arg, const Moose::StateArg & state) const
178 {
179  return evaluateHelper(node_arg, state);
180 }
181 
182 template <typename R>
183 typename Function::GradientType
184 Function::evaluateGradientHelper(const R & r, const Moose::StateArg & state) const
185 {
186  return gradient(_ti_feproblem.getTimeFromStateArg(state), r.getPoint());
187 }
188 
189 typename Function::GradientType
190 Function::evaluateGradient(const ElemArg & elem_arg, const Moose::StateArg & state) const
191 {
192  return evaluateGradientHelper(elem_arg, state);
193 }
194 
195 typename Function::GradientType
196 Function::evaluateGradient(const FaceArg & face, const Moose::StateArg & state) const
197 {
198  return evaluateGradientHelper(face, state);
199 }
200 
201 typename Function::GradientType
202 Function::evaluateGradient(const ElemQpArg & elem_qp, const Moose::StateArg & state) const
203 {
204  return evaluateGradientHelper(elem_qp, state);
205 }
206 
207 typename Function::GradientType
208 Function::evaluateGradient(const ElemSideQpArg & elem_side_qp, const Moose::StateArg & state) const
209 {
210  return evaluateGradientHelper(elem_side_qp, state);
211 }
212 
213 typename Function::GradientType
214 Function::evaluateGradient(const ElemPointArg & elem_point_arg, const Moose::StateArg & state) const
215 {
216  return evaluateGradientHelper(elem_point_arg, state);
217 }
218 
219 typename Function::GradientType
220 Function::evaluateGradient(const NodeArg & node_arg, const Moose::StateArg & state) const
221 {
222  return evaluateGradientHelper(node_arg, state);
223 }
224 
225 template <typename R>
226 typename Function::DotType
227 Function::evaluateDotHelper(const R & r, const Moose::StateArg & state) const
228 {
229  return timeDerivative(_ti_feproblem.getTimeFromStateArg(state), r.getPoint());
230 }
231 
232 typename Function::DotType
233 Function::evaluateDot(const ElemArg & elem_arg, const Moose::StateArg & state) const
234 {
235  return evaluateDotHelper(elem_arg, state);
236 }
237 
238 typename Function::DotType
239 Function::evaluateDot(const FaceArg & face, const Moose::StateArg & state) const
240 {
241  return evaluateDotHelper(face, state);
242 }
243 
244 typename Function::DotType
245 Function::evaluateDot(const ElemQpArg & elem_qp, const Moose::StateArg & state) const
246 {
247  return evaluateDotHelper(elem_qp, state);
248 }
249 
250 typename Function::DotType
251 Function::evaluateDot(const ElemSideQpArg & elem_side_qp, const Moose::StateArg & state) const
252 {
253  return evaluateDotHelper(elem_side_qp, state);
254 }
255 
256 typename Function::DotType
257 Function::evaluateDot(const ElemPointArg & elem_point_arg, const Moose::StateArg & state) const
258 {
259  return evaluateDotHelper(elem_point_arg, state);
260 }
261 
262 typename Function::DotType
263 Function::evaluateDot(const NodeArg & node_arg, const Moose::StateArg & state) const
264 {
265  return evaluateDotHelper(node_arg, state);
266 }
267 
268 void
270 {
272 }
273 
274 void
276 {
278 }
279 
280 void
282 {
284 }
285 
286 void
288 {
290 }
std::string name(const ElemQuality q)
A class for creating restricted objects.
Definition: Restartable.h:28
FEProblemBase & _ti_feproblem
Base class template for functor objects.
Definition: MooseFunctor.h:137
virtual Real div(Real t, const Point &p) const
Override this to evaluate the divergence of the vector function at a point (t,x,y,z), by default this returns zero, you must override it.
Definition: Function.C:101
DotType evaluateDotHelper(const R &r, const Moose::StateArg &state) const
Definition: Function.C:227
virtual RealVectorValue curl(Real t, const Point &p) const
Override this to evaluate the curl of the vector function at a point (t,x,y,z), by default this retur...
Definition: Function.C:94
virtual Real timeDerivative(Real t, const Point &p) const
Get the time derivative of the function.
Definition: Function.C:80
const Elem & elem() const
Definition: FaceInfo.h:80
DualNumber< Real, Real > ChainedReal
Definition: ChainedReal.h:30
auto raw_value(const Eigen::Map< T > &in)
Definition: ADReal.h:73
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
A structure that is used to evaluate Moose functors at an arbitrary physical point contained within a...
virtual const std::string & name() const
Get the name of the class.
Definition: MooseBase.h:56
typename FunctorReturnType< Real, FunctorEvaluationKind::Gradient >::type GradientType
This rigmarole makes it so that a user can create functors that return containers (std::vector...
Definition: MooseFunctor.h:149
Function(const InputParameters &parameters)
Definition: Function.C:25
virtual Real average() const
Definition: Function.C:115
void registerBase(const std::string &value)
This method must be called from every base "Moose System" to create linkage with the Action System...
void residualSetup() override
Gets called just before the residual is computed and before this object is asked to do its job...
Definition: Function.C:275
Interface for objects that needs transient capabilities.
const Elem * neighborPtr() const
Definition: FaceInfo.h:83
Interface for notifications that the mesh has changed.
A structure defining a "face" evaluation calling argument for Moose functors.
Every object that can be built by the factory should be derived from this class.
Definition: MooseObject.h:33
Real getTimeFromStateArg(const Moose::StateArg &state) const
Returns the time associated with the requested state.
const FaceInfo * fi
a face information object which defines our location in space
const Elem * face_side
A member that can be used to indicate whether there is a sidedness to this face.
const Elem & neighbor() const
Definition: FaceInfo.h:205
A structure that is used to evaluate Moose functors logically at an element/cell center.
DualReal ADReal
Definition: ADRealForward.h:14
Argument for requesting functor evaluation at a quadrature point location in an element.
Interface for objects that need to use UserObjects.
const Point & normal() const
Returns the unit normal vector for the face oriented outward from the face&#39;s elem element...
Definition: FaceInfo.h:67
void timestepSetup() override
Gets called at the beginning of the timestep before this object is asked to do its job...
Definition: Function.C:269
static InputParameters validParams()
const Elem * elemPtr() const
Definition: FaceInfo.h:81
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
GradientType evaluateGradientHelper(const R &r, const Moose::StateArg &state) const
Definition: Function.C:184
ValueType evaluate(const ElemArg &elem, const Moose::StateArg &state) const override final
Evaluate the functor with a given element.
Definition: Function.C:129
Class for containing MooseEnum item information.
Definition: MooseEnumItem.h:18
void customSetup(const ExecFlagType &exec_type) override
Gets called in FEProblemBase::execute() for execute flags other than initial, timestep_begin, nonlinear, linear and subdomain.
Definition: Function.C:287
virtual RealGradient gradient(Real t, const Point &p) const
Function objects can optionally provide a gradient at a point.
Definition: Function.C:73
DotType evaluateDot(const ElemArg &elem, const Moose::StateArg &state) const override final
Evaluate the functor time derivative with a given element.
Definition: Function.C:233
Interface for objects that needs scalar coupling capabilities.
virtual RealVectorValue vectorValue(Real t, const Point &p) const
Override this to evaluate the vector function at a point (t,x,y,z), by default this returns a zero ve...
Definition: Function.C:87
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
virtual Real integral() const
Definition: Function.C:108
State argument for evaluating functors.
virtual ~Function()
Function destructor.
Definition: Function.C:38
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
static InputParameters validParams()
Definition: MooseObject.C:24
virtual Real value(Real t, const Point &p) const
Override this to evaluate the scalar function at point (t,x,y,z), by default this returns zero...
Definition: Function.C:41
Point getPoint() const
GradientType evaluateGradient(const ElemArg &elem, const Moose::StateArg &state) const override final
Evaluate the functor gradient with a given element.
Definition: Function.C:190
static InputParameters validParams()
Class constructor.
Definition: Function.C:15
Argument for requesting functor evaluation at quadrature point locations on an element side...
ValueType evaluateHelper(const R &r, const Moose::StateArg &state) const
Definition: Function.C:123
Interface class for classes which interact with Postprocessors.
void jacobianSetup() override
Gets called just before the Jacobian is computed and before this object is asked to do its job...
Definition: Function.C:281