Line data Source code
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 "FunctionParserUtils.h"
11 :
12 : // MOOSE includes
13 : #include "InputParameters.h"
14 : #include "MooseEnum.h"
15 :
16 : template <bool is_ad>
17 : InputParameters
18 374692 : FunctionParserUtils<is_ad>::validParams()
19 : {
20 374692 : InputParameters params = emptyInputParameters();
21 :
22 1124076 : params.addParam<bool>(
23 : "enable_jit",
24 : #ifdef LIBMESH_HAVE_FPARSER_JIT
25 749384 : true,
26 : #else
27 : false,
28 : #endif
29 : "Enable just-in-time compilation of function expressions for faster evaluation");
30 1124076 : params.addParam<bool>(
31 749384 : "enable_ad_cache", true, "Enable caching of function derivatives for faster startup time");
32 1124076 : params.addParam<bool>(
33 749384 : "enable_auto_optimize", true, "Enable automatic immediate optimization of derivatives");
34 1124076 : params.addParam<bool>(
35 749384 : "disable_fpoptimizer", false, "Disable the function parser algebraic optimizer");
36 374692 : MooseEnum evalerror("nan nan_warning error exception", "nan");
37 374692 : params.addParam<MooseEnum>("evalerror_behavior",
38 : evalerror,
39 : "What to do if evaluation error occurs. Options are to pass a nan, "
40 : "pass a nan with a warning, throw a error, or throw an exception");
41 :
42 374692 : params.addParamNamesToGroup(
43 : "enable_jit enable_ad_cache enable_auto_optimize disable_fpoptimizer evalerror_behavior",
44 : "Parsed expression advanced");
45 374692 : params.addParam<Real>("epsilon", 0, "Fuzzy comparison tolerance");
46 749384 : return params;
47 374692 : }
48 :
49 : template <bool is_ad>
50 : const char * FunctionParserUtils<is_ad>::_eval_error_msg[] = {
51 : "Unknown",
52 : "Division by zero",
53 : "Square root of a negative value",
54 : "Logarithm of negative value",
55 : "Trigonometric error (asin or acos of illegal value)",
56 : "Maximum recursion level reached"};
57 :
58 : template <bool is_ad>
59 10435 : FunctionParserUtils<is_ad>::FunctionParserUtils(const InputParameters & parameters)
60 10435 : : _enable_jit(parameters.isParamValid("enable_jit") && parameters.get<bool>("enable_jit")),
61 10435 : _enable_ad_cache(parameters.get<bool>("enable_ad_cache")),
62 10435 : _disable_fpoptimizer(parameters.get<bool>("disable_fpoptimizer")),
63 10435 : _enable_auto_optimize(parameters.get<bool>("enable_auto_optimize") && !_disable_fpoptimizer),
64 10435 : _evalerror_behavior(parameters.get<MooseEnum>("evalerror_behavior").getEnum<FailureMethod>()),
65 10435 : _quiet_nan(std::numeric_limits<Real>::quiet_NaN()),
66 20870 : _epsilon(parameters.get<Real>("epsilon"))
67 : {
68 : #ifndef LIBMESH_HAVE_FPARSER_JIT
69 : if (_enable_jit)
70 : {
71 : mooseWarning("Tried to enable FParser JIT but libmesh does not have it compiled in.");
72 : _enable_jit = false;
73 : }
74 : #endif
75 10435 : }
76 :
77 : template <bool is_ad>
78 : void
79 11758 : FunctionParserUtils<is_ad>::setParserFeatureFlags(SymFunctionPtr & parser) const
80 : {
81 11758 : parser->SetADFlags(SymFunction::ADCacheDerivatives, _enable_ad_cache);
82 11758 : parser->SetADFlags(SymFunction::ADAutoOptimize, _enable_auto_optimize);
83 11758 : }
84 :
85 : template <bool is_ad>
86 : GenericReal<is_ad>
87 22567498 : FunctionParserUtils<is_ad>::evaluate(SymFunctionPtr & parser, const std::string & name)
88 : {
89 22567498 : return evaluate(parser, _func_params, name);
90 : }
91 :
92 : template <bool is_ad>
93 : GenericReal<is_ad>
94 22569406 : FunctionParserUtils<is_ad>::evaluate(SymFunctionPtr & parser,
95 : const std::vector<GenericReal<is_ad>> & func_params,
96 : const std::string & name)
97 : {
98 : // null pointer is a shortcut for vanishing derivatives, see functionsOptimize()
99 22569406 : if (parser == NULL)
100 0 : return 0.0;
101 :
102 : // set desired epsilon
103 22569406 : auto tmp_eps = parser->epsilon();
104 22569406 : parser->setEpsilon(_epsilon);
105 :
106 : // evaluate expression
107 22569406 : auto result = parser->Eval(func_params.data());
108 :
109 : // restore epsilon
110 22569406 : parser->setEpsilon(tmp_eps);
111 :
112 : // fetch fparser evaluation error (set to unknown if the JIT result is nan)
113 22569406 : int error_code = _enable_jit ? (std::isnan(result) ? -1 : 0) : parser->EvalError();
114 :
115 : // no error
116 22569406 : if (error_code == 0)
117 22569289 : return result;
118 :
119 : // hard fail or return not a number
120 117 : switch (_evalerror_behavior)
121 : {
122 42 : case FailureMethod::nan:
123 42 : return _quiet_nan;
124 :
125 42 : case FailureMethod::nan_warning:
126 42 : mooseWarning("In ",
127 : name,
128 : ": Parsed function evaluation encountered an error: ",
129 42 : _eval_error_msg[(error_code < 0 || error_code > 5) ? 0 : error_code]);
130 42 : return _quiet_nan;
131 :
132 12 : case FailureMethod::error:
133 12 : mooseError("In ",
134 : name,
135 : ": Parsed function evaluation encountered an error: ",
136 12 : _eval_error_msg[(error_code < 0 || error_code > 5) ? 0 : error_code]);
137 :
138 21 : case FailureMethod::exception:
139 21 : mooseException("In ",
140 : name,
141 : ": Parsed function evaluation encountered an error: ",
142 : _eval_error_msg[(error_code < 0 || error_code > 5) ? 0 : error_code],
143 : "\n Cutting timestep");
144 : }
145 :
146 0 : return _quiet_nan;
147 5390687 : }
148 :
149 : template <bool is_ad>
150 : void
151 11332 : FunctionParserUtils<is_ad>::addFParserConstants(
152 : SymFunctionPtr & parser,
153 : const std::vector<std::string> & constant_names,
154 : const std::vector<std::string> & constant_expressions) const
155 : {
156 : // check constant vectors
157 11332 : unsigned int nconst = constant_expressions.size();
158 11332 : if (nconst != constant_names.size())
159 0 : mooseError("The parameter vectors constant_names (size " +
160 : std::to_string(constant_names.size()) + ") and constant_expressions (size " +
161 : std::to_string(nconst) + ") must have equal length.");
162 :
163 : // previously evaluated constant_expressions may be used in following constant_expressions
164 11332 : std::vector<Real> constant_values(nconst);
165 :
166 13575 : for (unsigned int i = 0; i < nconst; ++i)
167 : {
168 : // no need to use dual numbers for the constant expressions
169 2243 : auto expression = std::make_shared<FunctionParserADBase<Real>>();
170 :
171 : // add previously evaluated constants
172 3638 : for (unsigned int j = 0; j < i; ++j)
173 1395 : if (!expression->AddConstant(constant_names[j], constant_values[j]))
174 0 : mooseError("Invalid constant name: ", constant_names[j], " and value ", constant_values[j]);
175 :
176 : // build the temporary constant expression function
177 2243 : if (expression->Parse(constant_expressions[i], "") >= 0)
178 0 : mooseError("Invalid constant expression\n",
179 : constant_expressions[i],
180 : "\n in parsed function object.\n",
181 0 : expression->ErrorMsg());
182 :
183 2243 : constant_values[i] = expression->Eval(NULL);
184 :
185 2243 : if (!parser->AddConstant(constant_names[i], constant_values[i]))
186 0 : mooseError("Invalid constant name in parsed function object");
187 : }
188 11332 : }
189 :
190 : template <>
191 : void
192 2256 : FunctionParserUtils<false>::functionsOptimize(SymFunctionPtr & parsed_function)
193 : {
194 : // set desired epsilon for optimization!
195 2256 : auto tmp_eps = parsed_function->epsilon();
196 2256 : parsed_function->setEpsilon(_epsilon);
197 :
198 : // base function
199 2256 : if (!_disable_fpoptimizer)
200 2028 : parsed_function->Optimize();
201 2256 : if (_enable_jit && !parsed_function->JITCompile())
202 0 : mooseInfo("Failed to JIT compile expression, falling back to byte code interpretation.");
203 :
204 2256 : parsed_function->setEpsilon(tmp_eps);
205 2256 : }
206 :
207 : template <>
208 : void
209 1605 : FunctionParserUtils<true>::functionsOptimize(SymFunctionPtr & parsed_function)
210 : {
211 : // set desired epsilon for optimization!
212 1605 : auto tmp_eps = parsed_function->epsilon();
213 1605 : parsed_function->setEpsilon(_epsilon);
214 :
215 : // base function
216 1605 : if (!_disable_fpoptimizer)
217 1491 : parsed_function->Optimize();
218 1605 : if (!_enable_jit || !parsed_function->JITCompile())
219 0 : mooseError("AD parsed objects require JIT compilation to be enabled and working.");
220 :
221 1605 : parsed_function->setEpsilon(tmp_eps);
222 1605 : }
223 :
224 : template <bool is_ad>
225 : void
226 6224 : FunctionParserUtils<is_ad>::parsedFunctionSetup(
227 : SymFunctionPtr & function,
228 : const std::string & expression,
229 : const std::string & variables,
230 : const std::vector<std::string> & constant_names,
231 : const std::vector<std::string> & constant_expressions,
232 : const libMesh::Parallel::Communicator & comm) const
233 : {
234 : // set FParser internal feature flags
235 6224 : setParserFeatureFlags(function);
236 :
237 : // add the constant expressions
238 6224 : addFParserConstants(function, constant_names, constant_expressions);
239 :
240 : // parse function
241 6224 : if (function->Parse(expression, variables) >= 0)
242 8 : mooseError("Invalid function\n", expression, "\nError:\n", function->ErrorMsg());
243 :
244 : // optimize
245 6216 : if (!_disable_fpoptimizer)
246 6216 : function->Optimize();
247 :
248 : // just-in-time compile
249 6216 : if (_enable_jit)
250 : {
251 : // let rank 0 do the JIT compilation first
252 6216 : if (comm.rank() != 0)
253 1177 : comm.barrier();
254 :
255 6216 : function->JITCompile();
256 :
257 : // wait for ranks > 0 to catch up
258 6216 : if (comm.rank() == 0)
259 5039 : comm.barrier();
260 : }
261 6216 : }
262 :
263 : // explicit instantiation
264 : template class FunctionParserUtils<false>;
265 : template class FunctionParserUtils<true>;
|
