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 376688 : FunctionParserUtils<is_ad>::validParams()
19 : {
20 376688 : InputParameters params = emptyInputParameters();
21 :
22 1130064 : params.addParam<bool>(
23 : "enable_jit",
24 : #ifdef LIBMESH_HAVE_FPARSER_JIT
25 753376 : true,
26 : #else
27 : false,
28 : #endif
29 : "Enable just-in-time compilation of function expressions for faster evaluation");
30 1130064 : params.addParam<bool>(
31 753376 : "enable_ad_cache", true, "Enable caching of function derivatives for faster startup time");
32 1130064 : params.addParam<bool>(
33 753376 : "enable_auto_optimize", true, "Enable automatic immediate optimization of derivatives");
34 1130064 : params.addParam<bool>(
35 753376 : "disable_fpoptimizer", false, "Disable the function parser algebraic optimizer");
36 376688 : MooseEnum evalerror("nan nan_warning error exception", "nan");
37 376688 : 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 376688 : params.addParamNamesToGroup(
43 : "enable_jit enable_ad_cache enable_auto_optimize disable_fpoptimizer evalerror_behavior",
44 : "Parsed expression advanced");
45 376688 : params.addParam<Real>("epsilon", 0, "Fuzzy comparison tolerance");
46 753376 : return params;
47 376688 : }
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 11535 : FunctionParserUtils<is_ad>::FunctionParserUtils(const InputParameters & parameters)
60 11535 : : _enable_jit(parameters.isParamValid("enable_jit") && parameters.get<bool>("enable_jit")),
61 11535 : _enable_ad_cache(parameters.get<bool>("enable_ad_cache")),
62 11535 : _disable_fpoptimizer(parameters.get<bool>("disable_fpoptimizer")),
63 11535 : _enable_auto_optimize(parameters.get<bool>("enable_auto_optimize") && !_disable_fpoptimizer),
64 11535 : _evalerror_behavior(parameters.get<MooseEnum>("evalerror_behavior").getEnum<FailureMethod>()),
65 11535 : _quiet_nan(std::numeric_limits<Real>::quiet_NaN()),
66 23070 : _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 11535 : }
76 :
77 : template <bool is_ad>
78 : void
79 12959 : FunctionParserUtils<is_ad>::setParserFeatureFlags(SymFunctionPtr & parser) const
80 : {
81 12959 : parser->SetADFlags(SymFunction::ADCacheDerivatives, _enable_ad_cache);
82 12959 : parser->SetADFlags(SymFunction::ADAutoOptimize, _enable_auto_optimize);
83 12959 : }
84 :
85 : template <bool is_ad>
86 : GenericReal<is_ad>
87 25964942 : FunctionParserUtils<is_ad>::evaluate(SymFunctionPtr & parser, const std::string & name)
88 : {
89 25964942 : return evaluate(parser, _func_params, name);
90 : }
91 :
92 : template <bool is_ad>
93 : GenericReal<is_ad>
94 25967038 : 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 25967038 : if (parser == NULL)
100 0 : return 0.0;
101 :
102 : // set desired epsilon
103 25967038 : auto tmp_eps = parser->epsilon();
104 25967038 : parser->setEpsilon(_epsilon);
105 :
106 : // evaluate expression
107 25967038 : auto result = parser->Eval(func_params.data());
108 :
109 : // restore epsilon
110 25967038 : parser->setEpsilon(tmp_eps);
111 :
112 : // fetch fparser evaluation error (set to unknown if the JIT result is nan)
113 25967038 : int error_code = _enable_jit ? (std::isnan(result) ? -1 : 0) : parser->EvalError();
114 :
115 : // no error
116 25967038 : if (error_code == 0)
117 25966915 : return result;
118 :
119 : // hard fail or return not a number
120 123 : switch (_evalerror_behavior)
121 : {
122 48 : case FailureMethod::nan:
123 48 : 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 6358887 : }
148 :
149 : template <bool is_ad>
150 : void
151 12492 : 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 12492 : unsigned int nconst = constant_expressions.size();
158 12492 : 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 12492 : std::vector<Real> constant_values(nconst);
165 :
166 15223 : for (unsigned int i = 0; i < nconst; ++i)
167 : {
168 : // no need to use dual numbers for the constant expressions
169 2731 : auto expression = std::make_shared<FunctionParserADBase<Real>>();
170 :
171 : // add previously evaluated constants
172 4402 : for (unsigned int j = 0; j < i; ++j)
173 1671 : 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 2731 : 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 2731 : constant_values[i] = expression->Eval(NULL);
184 :
185 2731 : if (!parser->AddConstant(constant_names[i], constant_values[i]))
186 0 : mooseError("Invalid constant name in parsed function object");
187 : }
188 12492 : }
189 :
190 : template <>
191 : void
192 2418 : FunctionParserUtils<false>::functionsOptimize(SymFunctionPtr & parsed_function)
193 : {
194 : // set desired epsilon for optimization!
195 2418 : auto tmp_eps = parsed_function->epsilon();
196 2418 : parsed_function->setEpsilon(_epsilon);
197 :
198 : // base function
199 2418 : if (!_disable_fpoptimizer)
200 2184 : parsed_function->Optimize();
201 2418 : if (_enable_jit && !parsed_function->JITCompile())
202 0 : mooseInfo("Failed to JIT compile expression, falling back to byte code interpretation.");
203 :
204 2418 : parsed_function->setEpsilon(tmp_eps);
205 2418 : }
206 :
207 : template <>
208 : void
209 1719 : FunctionParserUtils<true>::functionsOptimize(SymFunctionPtr & parsed_function)
210 : {
211 : // set desired epsilon for optimization!
212 1719 : auto tmp_eps = parsed_function->epsilon();
213 1719 : parsed_function->setEpsilon(_epsilon);
214 :
215 : // base function
216 1719 : if (!_disable_fpoptimizer)
217 1602 : parsed_function->Optimize();
218 1719 : if (!_enable_jit || !parsed_function->JITCompile())
219 0 : mooseError("AD parsed objects require JIT compilation to be enabled and working.");
220 :
221 1719 : parsed_function->setEpsilon(tmp_eps);
222 1719 : }
223 :
224 : template <bool is_ad>
225 : void
226 7021 : 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 7021 : setParserFeatureFlags(function);
236 :
237 : // add the constant expressions
238 7021 : addFParserConstants(function, constant_names, constant_expressions);
239 :
240 : // parse function
241 7021 : if (function->Parse(expression, variables) >= 0)
242 8 : mooseError("Invalid function\n", expression, "\nError:\n", function->ErrorMsg());
243 :
244 : // optimize
245 7013 : if (!_disable_fpoptimizer)
246 7013 : function->Optimize();
247 :
248 : // just-in-time compile
249 7013 : if (_enable_jit)
250 : {
251 : // let rank 0 do the JIT compilation first
252 7013 : if (comm.rank() != 0)
253 1242 : comm.barrier();
254 :
255 7013 : function->JITCompile();
256 :
257 : // wait for ranks > 0 to catch up
258 7013 : if (comm.rank() == 0)
259 5771 : comm.barrier();
260 : }
261 7013 : }
262 :
263 : // explicit instantiation
264 : template class FunctionParserUtils<false>;
265 : template class FunctionParserUtils<true>;
|
