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 468484 : FunctionParserUtils<is_ad>::validParams()
19 : {
20 468484 : InputParameters params = emptyInputParameters();
21 :
22 1405452 : params.addParam<bool>(
23 : "enable_jit",
24 : #ifdef LIBMESH_HAVE_FPARSER_JIT
25 936968 : true,
26 : #else
27 : false,
28 : #endif
29 : "Enable just-in-time compilation of function expressions for faster evaluation");
30 1405452 : params.addParam<bool>(
31 936968 : "enable_ad_cache", true, "Enable caching of function derivatives for faster startup time");
32 1405452 : params.addParam<bool>(
33 936968 : "enable_auto_optimize", true, "Enable automatic immediate optimization of derivatives");
34 1405452 : params.addParam<bool>(
35 936968 : "disable_fpoptimizer", false, "Disable the function parser algebraic optimizer");
36 1873936 : MooseEnum evalerror("nan nan_warning error exception", "nan");
37 1873936 : 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 1873936 : params.addParamNamesToGroup(
43 : "enable_jit enable_ad_cache enable_auto_optimize disable_fpoptimizer evalerror_behavior",
44 : "Parsed expression advanced");
45 1405452 : params.addParam<Real>("epsilon", 0, "Fuzzy comparison tolerance");
46 936968 : return params;
47 468484 : }
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 12834 : FunctionParserUtils<is_ad>::FunctionParserUtils(const InputParameters & parameters)
60 12834 : : _enable_jit(parameters.isParamValid("enable_jit") && parameters.get<bool>("enable_jit")),
61 12834 : _enable_ad_cache(parameters.get<bool>("enable_ad_cache")),
62 12834 : _disable_fpoptimizer(parameters.get<bool>("disable_fpoptimizer")),
63 12834 : _enable_auto_optimize(parameters.get<bool>("enable_auto_optimize") && !_disable_fpoptimizer),
64 12834 : _evalerror_behavior(parameters.get<MooseEnum>("evalerror_behavior").getEnum<FailureMethod>()),
65 12834 : _quiet_nan(std::numeric_limits<Real>::quiet_NaN()),
66 25668 : _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 12834 : }
76 :
77 : template <bool is_ad>
78 : void
79 14098 : FunctionParserUtils<is_ad>::setParserFeatureFlags(SymFunctionPtr & parser) const
80 : {
81 14098 : parser->SetADFlags(SymFunction::ADCacheDerivatives, _enable_ad_cache);
82 14098 : parser->SetADFlags(SymFunction::ADAutoOptimize, _enable_auto_optimize);
83 14098 : }
84 :
85 : template <bool is_ad>
86 : GenericReal<is_ad>
87 28717760 : FunctionParserUtils<is_ad>::evaluate(SymFunctionPtr & parser, const std::string & name)
88 : {
89 28717760 : return evaluate(parser, _func_params, name);
90 : }
91 :
92 : template <bool is_ad>
93 : GenericReal<is_ad>
94 28731995 : FunctionParserUtils<is_ad>::evaluate(SymFunctionPtr & parser,
95 : const std::vector<GenericReal<is_ad>> & func_params,
96 : const std::string & name)
97 : {
98 : using std::isnan;
99 :
100 : // null pointer is a shortcut for vanishing derivatives, see functionsOptimize()
101 28731995 : if (parser == NULL)
102 0 : return 0.0;
103 :
104 : // set desired epsilon
105 28731995 : auto tmp_eps = parser->epsilon();
106 28731995 : parser->setEpsilon(_epsilon);
107 :
108 : // evaluate expression
109 28731995 : auto result = parser->Eval(func_params.data());
110 :
111 : // restore epsilon
112 28731995 : parser->setEpsilon(tmp_eps);
113 :
114 : // fetch fparser evaluation error (set to unknown if the JIT result is nan)
115 28731995 : int error_code = _enable_jit ? (isnan(result) ? -1 : 0) : parser->EvalError();
116 :
117 : // no error
118 28731995 : if (error_code == 0)
119 28731872 : return result;
120 :
121 : // hard fail or return not a number
122 123 : switch (_evalerror_behavior)
123 : {
124 48 : case FailureMethod::nan:
125 48 : return _quiet_nan;
126 :
127 42 : case FailureMethod::nan_warning:
128 42 : mooseWarning("In ",
129 : name,
130 : ": Parsed function evaluation encountered an error: ",
131 42 : _eval_error_msg[(error_code < 0 || error_code > 5) ? 0 : error_code]);
132 42 : return _quiet_nan;
133 :
134 12 : case FailureMethod::error:
135 12 : mooseError("In ",
136 : name,
137 : ": Parsed function evaluation encountered an error: ",
138 12 : _eval_error_msg[(error_code < 0 || error_code > 5) ? 0 : error_code]);
139 :
140 21 : case FailureMethod::exception:
141 21 : mooseException("In ",
142 : name,
143 : ": Parsed function evaluation encountered an error: ",
144 : _eval_error_msg[(error_code < 0 || error_code > 5) ? 0 : error_code],
145 : "\n Cutting timestep");
146 : }
147 :
148 0 : return _quiet_nan;
149 6629591 : }
150 :
151 : template <bool is_ad>
152 : void
153 13556 : FunctionParserUtils<is_ad>::addFParserConstants(
154 : SymFunctionPtr & parser,
155 : const std::vector<std::string> & constant_names,
156 : const std::vector<std::string> & constant_expressions) const
157 : {
158 : // check constant vectors
159 13556 : unsigned int nconst = constant_expressions.size();
160 13556 : if (nconst != constant_names.size())
161 0 : mooseError("The parameter vectors constant_names (size " +
162 : std::to_string(constant_names.size()) + ") and constant_expressions (size " +
163 : std::to_string(nconst) + ") must have equal length.");
164 :
165 : // previously evaluated constant_expressions may be used in following constant_expressions
166 13556 : std::vector<Real> constant_values(nconst);
167 :
168 16548 : for (unsigned int i = 0; i < nconst; ++i)
169 : {
170 : // no need to use dual numbers for the constant expressions
171 2992 : auto expression = std::make_shared<FunctionParserADBase<Real>>();
172 :
173 : // add previously evaluated constants
174 4795 : for (unsigned int j = 0; j < i; ++j)
175 1803 : if (!expression->AddConstant(constant_names[j], constant_values[j]))
176 0 : mooseError("Invalid constant name: ", constant_names[j], " and value ", constant_values[j]);
177 :
178 : // build the temporary constant expression function
179 8976 : if (expression->Parse(constant_expressions[i], "") >= 0)
180 0 : mooseError("Invalid constant expression\n",
181 : constant_expressions[i],
182 : "\n in parsed function object.\n",
183 0 : expression->ErrorMsg());
184 :
185 2992 : constant_values[i] = expression->Eval(NULL);
186 :
187 2992 : if (!parser->AddConstant(constant_names[i], constant_values[i]))
188 0 : mooseError("Invalid constant name in parsed function object");
189 : }
190 13556 : }
191 :
192 : template <>
193 : void
194 2502 : FunctionParserUtils<false>::functionsOptimize(SymFunctionPtr & parsed_function)
195 : {
196 : // set desired epsilon for optimization!
197 2502 : auto tmp_eps = parsed_function->epsilon();
198 2502 : parsed_function->setEpsilon(_epsilon);
199 :
200 : // base function
201 2502 : if (!_disable_fpoptimizer)
202 2268 : parsed_function->Optimize();
203 2502 : if (_enable_jit && !parsed_function->JITCompile())
204 0 : mooseInfo("Failed to JIT compile expression, falling back to byte code interpretation.");
205 :
206 2502 : parsed_function->setEpsilon(tmp_eps);
207 2502 : }
208 :
209 : template <>
210 : void
211 1929 : FunctionParserUtils<true>::functionsOptimize(SymFunctionPtr & parsed_function)
212 : {
213 : // set desired epsilon for optimization!
214 1929 : auto tmp_eps = parsed_function->epsilon();
215 1929 : parsed_function->setEpsilon(_epsilon);
216 :
217 : // base function
218 1929 : if (!_disable_fpoptimizer)
219 1812 : parsed_function->Optimize();
220 1929 : if (!_enable_jit || !parsed_function->JITCompile())
221 0 : mooseError("AD parsed objects require JIT compilation to be enabled and working.");
222 :
223 1929 : parsed_function->setEpsilon(tmp_eps);
224 1929 : }
225 :
226 : template <bool is_ad>
227 : void
228 7778 : FunctionParserUtils<is_ad>::parsedFunctionSetup(
229 : SymFunctionPtr & function,
230 : const std::string & expression,
231 : const std::string & variables,
232 : const std::vector<std::string> & constant_names,
233 : const std::vector<std::string> & constant_expressions,
234 : const libMesh::Parallel::Communicator & comm) const
235 : {
236 : // set FParser internal feature flags
237 7778 : setParserFeatureFlags(function);
238 :
239 : // add the constant expressions
240 7778 : addFParserConstants(function, constant_names, constant_expressions);
241 :
242 : // parse function
243 7778 : if (function->Parse(expression, variables) >= 0)
244 8 : mooseError("Invalid function\n", expression, "\nError:\n", function->ErrorMsg());
245 :
246 : // optimize
247 7770 : if (!_disable_fpoptimizer)
248 7770 : function->Optimize();
249 :
250 : // just-in-time compile
251 7770 : if (_enable_jit)
252 : {
253 : // let rank 0 do the JIT compilation first
254 7770 : if (comm.rank() != 0)
255 1380 : comm.barrier();
256 :
257 7770 : function->JITCompile();
258 :
259 : // wait for ranks > 0 to catch up
260 7770 : if (comm.rank() == 0)
261 6390 : comm.barrier();
262 : }
263 7770 : }
264 :
265 : // explicit instantiation
266 : template class FunctionParserUtils<false>;
267 : template class FunctionParserUtils<true>;
|
