libMesh
parsed_fem_function.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_PARSED_FEM_FUNCTION_H
21 #define LIBMESH_PARSED_FEM_FUNCTION_H
22 
23 #include "libmesh/libmesh_config.h"
24 
25 // Local Includes
26 #include "libmesh/fem_function_base.h"
27 #include "libmesh/point.h"
28 #include "libmesh/system.h"
29 
30 #ifdef LIBMESH_HAVE_FPARSER
31 // FParser includes
32 #include "libmesh/fparser.hh"
33 #endif
34 
35 // C++ includes
36 #include <cctype>
37 #include <iomanip>
38 #include <sstream>
39 #include <string>
40 #include <vector>
41 
42 namespace libMesh
43 {
44 
54 template <typename Output=Number>
55 class ParsedFEMFunction : public FEMFunctionBase<Output>
56 {
57 public:
58 
62  explicit
63  ParsedFEMFunction (const System & sys,
64  const std::string & expression,
65  const std::vector<std::string> * additional_vars=libmesh_nullptr,
66  const std::vector<Output> * initial_vals=libmesh_nullptr);
67 
71  virtual ~ParsedFEMFunction () {}
72 
76  void reparse (const std::string & expression);
77 
78  virtual void init_context (const FEMContext & c) libmesh_override;
79 
80  virtual UniquePtr<FEMFunctionBase<Output>> clone () const libmesh_override;
81 
82  virtual Output operator() (const FEMContext & c,
83  const Point & p,
84  const Real time = 0.) libmesh_override;
85 
86  void operator() (const FEMContext & c,
87  const Point & p,
88  const Real time,
89  DenseVector<Output> & output) libmesh_override;
90 
91  virtual Output component(const FEMContext & c,
92  unsigned int i,
93  const Point & p,
94  Real time=0.) libmesh_override;
95 
96  const std::string & expression() { return _expression; }
97 
106  Output get_inline_value(const std::string & inline_var_name) const;
107 
118  void set_inline_value(const std::string & inline_var_name,
119  Output newval);
120 
121 protected:
122  // Helper function for reparsing minor changes to expression
123  void partial_reparse (const std::string & expression);
124 
125  // Helper function for parsing out variable names
126  std::size_t find_name (const std::string & varname,
127  const std::string & expr) const;
128 
129  // Helper function for evaluating function arguments
130  void eval_args(const FEMContext & c,
131  const Point & p,
132  const Real time);
133 
134  // Evaluate the ith FunctionParser and check the result
135 #ifdef LIBMESH_HAVE_FPARSER
136  inline Output eval(FunctionParserBase<Output> & parser,
137  const std::string & libmesh_dbg_var(function_name),
138  unsigned int libmesh_dbg_var(component_idx)) const;
139 #else // LIBMESH_HAVE_FPARSER
140  inline Output eval(char & libmesh_dbg_var(parser),
141  const std::string & libmesh_dbg_var(function_name),
142  unsigned int libmesh_dbg_var(component_idx)) const;
143 #endif
144 
145 private:
146  const System & _sys;
147  std::string _expression;
148  std::vector<std::string> _subexpressions;
149  unsigned int _n_vars,
154 #ifdef LIBMESH_HAVE_FPARSER
155  std::vector<FunctionParserBase<Output>> parsers;
156 #else
157  std::vector<char> parsers;
158 #endif
159  std::vector<Output> _spacetime;
160 
161  // Flags for which variables need to be computed
162 
163  // _need_var[v] is true iff value(v) is needed
164  std::vector<bool> _need_var;
165 
166  // _need_var_grad[v*LIBMESH_DIM+i] is true iff grad(v,i) is needed
167  std::vector<bool> _need_var_grad;
168 
169 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
170  // _need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+i*LIBMESH_DIM+j] is true
171  // iff grad(v,i,j) is needed
172  std::vector<bool> _need_var_hess;
173 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
174 
175  // Variables/values that can be parsed and handled by the function parser
176  std::string variables;
177  std::vector<std::string> _additional_vars;
178  std::vector<Output> _initial_vals;
179 };
180 
181 
182 /*----------------------- Inline functions ----------------------------------*/
183 
184 template <typename Output>
185 inline
187  const std::string & expression,
188  const std::vector<std::string> * additional_vars,
189  const std::vector<Output> * initial_vals) :
190  _sys(sys),
191  _expression (), // overridden by parse()
192  _n_vars(sys.n_vars()),
196  _requested_normals(false),
197  _need_var(_n_vars, false),
198  _need_var_grad(_n_vars*LIBMESH_DIM, false),
199 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
200  _need_var_hess(_n_vars*LIBMESH_DIM*LIBMESH_DIM, false),
201 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
202  _additional_vars (additional_vars ? *additional_vars : std::vector<std::string>()),
203  _initial_vals (initial_vals ? *initial_vals : std::vector<Output>())
204 {
205  this->reparse(expression);
206 }
207 
208 
209 template <typename Output>
210 inline
211 void
213 {
214  variables = "x";
215 #if LIBMESH_DIM > 1
216  variables += ",y";
217 #endif
218 #if LIBMESH_DIM > 2
219  variables += ",z";
220 #endif
221  variables += ",t";
222 
223  for (unsigned int v=0; v != _n_vars; ++v)
224  {
225  const std::string & varname = _sys.variable_name(v);
226  std::size_t varname_i = find_name(varname, expression);
227 
228  // If we didn't find our variable name then let's go to the
229  // next.
230  if (varname_i == std::string::npos)
231  continue;
232 
233  _need_var[v] = true;
234  variables += ',';
235  variables += varname;
237  }
238 
239  for (unsigned int v=0; v != _n_vars; ++v)
240  {
241  const std::string & varname = _sys.variable_name(v);
242 
243  for (unsigned int d=0; d != LIBMESH_DIM; ++d)
244  {
245  std::string gradname = std::string("grad_") +
246  "xyz"[d] + '_' + varname;
247  std::size_t gradname_i = find_name(gradname, expression);
248 
249  // If we didn't find that gradient component of our
250  // variable name then let's go to the next.
251  if (gradname_i == std::string::npos)
252  continue;
253 
254  _need_var_grad[v*LIBMESH_DIM+d] = true;
255  variables += ',';
256  variables += gradname;
258  }
259  }
260 
261 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
262  for (unsigned int v=0; v != _n_vars; ++v)
263  {
264  const std::string & varname = _sys.variable_name(v);
265 
266  for (unsigned int d1=0; d1 != LIBMESH_DIM; ++d1)
267  for (unsigned int d2=0; d2 != LIBMESH_DIM; ++d2)
268  {
269  std::string hessname = std::string("hess_") +
270  "xyz"[d1] + "xyz"[d2] + '_' + varname;
271  std::size_t hessname_i = find_name(hessname, expression);
272 
273  // If we didn't find that hessian component of our
274  // variable name then let's go to the next.
275  if (hessname_i == std::string::npos)
276  continue;
277 
278  _need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+d1*LIBMESH_DIM+d2]
279  = true;
280  variables += ',';
281  variables += hessname;
283  }
284  }
285 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
286 
287  {
288  std::size_t nx_i = find_name("n_x", expression);
289  std::size_t ny_i = find_name("n_y", expression);
290  std::size_t nz_i = find_name("n_z", expression);
291 
292  // If we found any requests for normal components then we'll
293  // compute normals
294  if (nx_i != std::string::npos ||
295  ny_i != std::string::npos ||
296  nz_i != std::string::npos)
297  {
298  _requested_normals = true;
299  variables += ',';
300  variables += "n_x";
301  if (LIBMESH_DIM > 1)
302  {
303  variables += ',';
304  variables += "n_y";
305  }
306  if (LIBMESH_DIM > 2)
307  {
308  variables += ',';
309  variables += "n_z";
310  }
311  }
312  }
313 
314  _spacetime.resize
315  (LIBMESH_DIM + 1 + _n_requested_vars +
317  (_requested_normals ? LIBMESH_DIM : 0) +
318  _additional_vars.size());
319 
320  // If additional vars were passed, append them to the string
321  // that we send to the function parser. Also add them to the
322  // end of our spacetime vector
323  unsigned int offset = LIBMESH_DIM + 1 + _n_requested_vars +
325 
326  for (std::size_t i=0; i < _additional_vars.size(); ++i)
327  {
328  variables += "," + _additional_vars[i];
329  // Initialize extra variables to the vector passed in or zero
330  // Note: The initial_vals vector can be shorter than the additional_vars vector
331  _spacetime[offset + i] =
332  (i < _initial_vals.size()) ? _initial_vals[i] : 0;
333  }
334 
335  this->partial_reparse(expression);
336 }
337 
338 template <typename Output>
339 inline
340 void
342 {
343  for (unsigned int v=0; v != _n_vars; ++v)
344  {
345  FEBase * elem_fe;
346  c.get_element_fe(v, elem_fe);
347  if (_n_requested_vars)
348  elem_fe->get_phi();
350  elem_fe->get_dphi();
351 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
353  elem_fe->get_d2phi();
354 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
355  }
356 
357  if (_requested_normals)
358  {
359  FEBase * side_fe;
360  c.get_side_fe(0, side_fe);
361 
362  side_fe->get_normals();
363 
364  // FIXME: this is a hack to support normals at quadrature
365  // points; we don't support normals elsewhere.
366  side_fe->get_xyz();
367  }
368 }
369 
370 template <typename Output>
371 inline
374 {
377 }
378 
379 template <typename Output>
380 inline
381 Output
383  const Point & p,
384  const Real time)
385 {
386  eval_args(c, p, time);
387 
388  return eval(parsers[0], "f", 0);
389 }
390 
391 
392 
393 template <typename Output>
394 inline
395 void
397  const Point & p,
398  const Real time,
399  DenseVector<Output> & output)
400 {
401  eval_args(c, p, time);
402 
403  unsigned int size = output.size();
404 
405  libmesh_assert_equal_to (size, parsers.size());
406 
407  for (unsigned int i=0; i != size; ++i)
408  output(i) = eval(parsers[i], "f", i);
409 }
410 
411 
412 template <typename Output>
413 inline
414 Output
416  unsigned int i,
417  const Point & p,
418  Real time)
419 {
420  eval_args(c, p, time);
421 
422  libmesh_assert_less (i, parsers.size());
423  return eval(parsers[i], "f", i);
424 }
425 
426 template <typename Output>
427 inline
428 Output
429 ParsedFEMFunction<Output>::get_inline_value(const std::string & inline_var_name) const
430 {
431  libmesh_assert_greater (_subexpressions.size(), 0);
432 
433 #ifndef NDEBUG
434  bool found_var_name = false;
435 #endif
436  Output old_var_value(0.);
437 
438  for (std::size_t s=0; s != _subexpressions.size(); ++s)
439  {
440  const std::string & subexpression = _subexpressions[s];
441  const std::size_t varname_i =
442  find_name(inline_var_name, subexpression);
443  if (varname_i == std::string::npos)
444  continue;
445 
446  const std::size_t assignment_i =
447  subexpression.find(":", varname_i+1);
448 
449  libmesh_assert_not_equal_to(assignment_i, std::string::npos);
450 
451  libmesh_assert_equal_to(subexpression[assignment_i+1], '=');
452  for (unsigned int i = varname_i+1; i != assignment_i; ++i)
453  libmesh_assert_equal_to(subexpression[i], ' ');
454 
455  std::size_t end_assignment_i =
456  subexpression.find(";", assignment_i+1);
457 
458  libmesh_assert_not_equal_to(end_assignment_i, std::string::npos);
459 
460  std::string new_subexpression =
461  subexpression.substr(0, end_assignment_i+1) +
462  inline_var_name;
463 
464 #ifdef LIBMESH_HAVE_FPARSER
465  // Parse and evaluate the new subexpression.
466  // Add the same constants as we used originally.
467  FunctionParserBase<Output> fp;
468  fp.AddConstant("NaN", std::numeric_limits<Real>::quiet_NaN());
469  fp.AddConstant("pi", std::acos(Real(-1)));
470  fp.AddConstant("e", std::exp(Real(1)));
471  if (fp.Parse(new_subexpression, variables) != -1) // -1 for success
472  libmesh_error_msg
473  ("ERROR: FunctionParser is unable to parse modified expression: "
474  << new_subexpression << '\n' << fp.ErrorMsg());
475 
476  Output new_var_value = this->eval(fp, new_subexpression, 0);
477 #ifdef NDEBUG
478  return new_var_value;
479 #else
480  if (found_var_name)
481  {
482  libmesh_assert_equal_to(old_var_value, new_var_value);
483  }
484  else
485  {
486  old_var_value = new_var_value;
487  found_var_name = true;
488  }
489 #endif
490 
491 #else
492  libmesh_error_msg("ERROR: This functionality requires fparser!");
493 #endif
494  }
495 
496  libmesh_assert(found_var_name);
497  return old_var_value;
498 }
499 
500 template <typename Output>
501 inline
502 void
503 ParsedFEMFunction<Output>::set_inline_value (const std::string & inline_var_name,
504  Output newval)
505 {
506  libmesh_assert_greater (_subexpressions.size(), 0);
507 
508 #ifndef NDEBUG
509  bool found_var_name = false;
510 #endif
511  for (std::size_t s=0; s != _subexpressions.size(); ++s)
512  {
513  const std::string & subexpression = _subexpressions[s];
514  const std::size_t varname_i =
515  find_name(inline_var_name, subexpression);
516  if (varname_i == std::string::npos)
517  continue;
518 
519 #ifndef NDEBUG
520  found_var_name = true;
521 #endif
522 
523  const std::size_t assignment_i =
524  subexpression.find(":", varname_i+1);
525 
526  libmesh_assert_not_equal_to(assignment_i, std::string::npos);
527 
528  libmesh_assert_equal_to(subexpression[assignment_i+1], '=');
529  for (unsigned int i = varname_i+1; i != assignment_i; ++i)
530  libmesh_assert_equal_to(subexpression[i], ' ');
531 
532  std::size_t end_assignment_i =
533  subexpression.find(";", assignment_i+1);
534 
535  libmesh_assert_not_equal_to(end_assignment_i, std::string::npos);
536 
537  std::ostringstream new_subexpression;
538  new_subexpression << subexpression.substr(0, assignment_i+2)
539  << std::setprecision(std::numeric_limits<Output>::digits10+2)
540 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
541  << '(' << newval.real() << '+'
542  << newval.imag() << 'i' << ')'
543 #else
544  << newval
545 #endif
546  << subexpression.substr(end_assignment_i,
547  std::string::npos);
548  _subexpressions[s] = new_subexpression.str();
549  }
550 
551  libmesh_assert(found_var_name);
552 
553  std::string new_expression;
554 
555  for (std::size_t s=0; s != _subexpressions.size(); ++s)
556  {
557  new_expression += '{';
558  new_expression += _subexpressions[s];
559  new_expression += '}';
560  }
561 
562  this->partial_reparse(new_expression);
563 }
564 
565 
566 template <typename Output>
567 inline
568 void
570 {
572  _subexpressions.clear();
573  parsers.clear();
574 
575  size_t nextstart = 0, end = 0;
576 
577  while (end != std::string::npos)
578  {
579  // If we're past the end of the string, we can't make any more
580  // subparsers
581  if (nextstart >= expression.size())
582  break;
583 
584  // If we're at the start of a brace delimited section, then we
585  // parse just that section:
586  if (expression[nextstart] == '{')
587  {
588  nextstart++;
589  end = expression.find('}', nextstart);
590  }
591  // otherwise we parse the whole thing
592  else
593  end = std::string::npos;
594 
595  // We either want the whole end of the string (end == npos) or
596  // a substring in the middle.
597  _subexpressions.push_back
598  (expression.substr(nextstart, (end == std::string::npos) ?
599  std::string::npos : end - nextstart));
600 
601  // fparser can crash on empty expressions
602  if (_subexpressions.back().empty())
603  libmesh_error_msg("ERROR: FunctionParser is unable to parse empty expression.\n");
604 
605 
606 #ifdef LIBMESH_HAVE_FPARSER
607  // Parse (and optimize if possible) the subexpression.
608  // Add some basic constants, to Real precision.
609  FunctionParserBase<Output> fp;
610  fp.AddConstant("NaN", std::numeric_limits<Real>::quiet_NaN());
611  fp.AddConstant("pi", std::acos(Real(-1)));
612  fp.AddConstant("e", std::exp(Real(1)));
613  if (fp.Parse(_subexpressions.back(), variables) != -1) // -1 for success
614  libmesh_error_msg
615  ("ERROR: FunctionParser is unable to parse expression: "
616  << _subexpressions.back() << '\n' << fp.ErrorMsg());
617  fp.Optimize();
618  parsers.push_back(fp);
619 #else
620  libmesh_error_msg("ERROR: This functionality requires fparser!");
621 #endif
622 
623  // If at end, use nextstart=maxSize. Else start at next
624  // character.
625  nextstart = (end == std::string::npos) ?
626  std::string::npos : end + 1;
627  }
628 }
629 
630 template <typename Output>
631 inline
632 std::size_t
633 ParsedFEMFunction<Output>::find_name (const std::string & varname,
634  const std::string & expr) const
635 {
636  const std::size_t namesize = varname.size();
637  std::size_t varname_i = expr.find(varname);
638 
639  while ((varname_i != std::string::npos) &&
640  (((varname_i > 0) &&
641  (std::isalnum(expr[varname_i-1]) ||
642  (expr[varname_i-1] == '_'))) ||
643  ((varname_i+namesize < expr.size()) &&
644  (std::isalnum(expr[varname_i+namesize]) ||
645  (expr[varname_i+namesize] == '_')))))
646  {
647  varname_i = expr.find(varname, varname_i+1);
648  }
649 
650  return varname_i;
651 }
652 
653 
654 
655 // Helper function for evaluating function arguments
656 template <typename Output>
657 inline
658 void
660  const Point & p,
661  const Real time)
662 {
663  _spacetime[0] = p(0);
664 #if LIBMESH_DIM > 1
665  _spacetime[1] = p(1);
666 #endif
667 #if LIBMESH_DIM > 2
668  _spacetime[2] = p(2);
669 #endif
670  _spacetime[LIBMESH_DIM] = time;
671 
672  unsigned int request_index = 0;
673  for (unsigned int v=0; v != _n_vars; ++v)
674  {
675  if (!_need_var[v])
676  continue;
677 
678  c.point_value(v, p, _spacetime[LIBMESH_DIM+1+request_index]);
679  request_index++;
680  }
681 
683  for (unsigned int v=0; v != _n_vars; ++v)
684  {
685  if (!_need_var_grad[v*LIBMESH_DIM]
686 #if LIBMESH_DIM > 1
687  && !_need_var_grad[v*LIBMESH_DIM+1]
688 #if LIBMESH_DIM > 2
689  && !_need_var_grad[v*LIBMESH_DIM+2]
690 #endif
691 #endif
692  )
693  continue;
694 
695  Gradient g;
696  c.point_gradient(v, p, g);
697 
698  for (unsigned int d=0; d != LIBMESH_DIM; ++d)
699  {
700  if (!_need_var_grad[v*LIBMESH_DIM+d])
701  continue;
702 
703  _spacetime[LIBMESH_DIM+1+request_index] = g(d);
704  request_index++;
705  }
706  }
707 
708 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
710  for (unsigned int v=0; v != _n_vars; ++v)
711  {
712  if (!_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM]
713 #if LIBMESH_DIM > 1
714  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+1]
715  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+2]
716  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+3]
717 #if LIBMESH_DIM > 2
718  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+4]
719  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+5]
720  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+6]
721  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+7]
722  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+8]
723 #endif
724 #endif
725  )
726  continue;
727 
728  Tensor h;
729  c.point_hessian(v, p, h);
730 
731  for (unsigned int d1=0; d1 != LIBMESH_DIM; ++d1)
732  for (unsigned int d2=0; d2 != LIBMESH_DIM; ++d2)
733  {
734  if (!_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+d1*LIBMESH_DIM+d2])
735  continue;
736 
737  _spacetime[LIBMESH_DIM+1+request_index] = h(d1,d2);
738  request_index++;
739  }
740  }
741 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
742 
743  if (_requested_normals)
744  {
745  FEBase * side_fe;
746  c.get_side_fe(0, side_fe);
747 
748  const std::vector<Point> & normals = side_fe->get_normals();
749 
750  const std::vector<Point> & xyz = side_fe->get_xyz();
751 
752  libmesh_assert_equal_to(normals.size(), xyz.size());
753 
754  // We currently only support normals at quadrature points!
755 #ifndef NDEBUG
756  bool at_quadrature_point = false;
757 #endif
758  for (std::size_t qp = 0; qp != normals.size(); ++qp)
759  {
760  if (p == xyz[qp])
761  {
762  const Point & n = normals[qp];
763  for (unsigned int d=0; d != LIBMESH_DIM; ++d)
764  {
765  _spacetime[LIBMESH_DIM+1+request_index] = n(d);
766  request_index++;
767  }
768 #ifndef NDEBUG
769  at_quadrature_point = true;
770 #endif
771  break;
772  }
773  }
774 
775  libmesh_assert(at_quadrature_point);
776  }
777 
778  // The remaining locations in _spacetime are currently set at construction
779  // but could potentially be made dynamic
780 }
781 
782 
783 // Evaluate the ith FunctionParser and check the result
784 #ifdef LIBMESH_HAVE_FPARSER
785 template <typename Output>
786 inline
787 Output
788 ParsedFEMFunction<Output>::eval (FunctionParserBase<Output> & parser,
789  const std::string & libmesh_dbg_var(function_name),
790  unsigned int libmesh_dbg_var(component_idx)) const
791 {
792 #ifndef NDEBUG
793  Output result = parser.Eval(&_spacetime[0]);
794  int error_code = parser.EvalError();
795  if (error_code)
796  {
797  libMesh::err << "ERROR: FunctionParser is unable to evaluate component "
798  << component_idx
799  << " of expression '"
800  << function_name
801  << "' with arguments:\n";
802  for (std::size_t j=0; j<_spacetime.size(); ++j)
803  libMesh::err << '\t' << _spacetime[j] << '\n';
804  libMesh::err << '\n';
805 
806  // Currently no API to report error messages, we'll do it manually
807  std::string error_message = "Reason: ";
808 
809  switch (error_code)
810  {
811  case 1:
812  error_message += "Division by zero";
813  break;
814  case 2:
815  error_message += "Square Root error (negative value)";
816  break;
817  case 3:
818  error_message += "Log error (negative value)";
819  break;
820  case 4:
821  error_message += "Trigonometric error (asin or acos of illegal value)";
822  break;
823  case 5:
824  error_message += "Maximum recursion level reached";
825  break;
826  default:
827  error_message += "Unknown";
828  break;
829  }
830  libmesh_error_msg(error_message);
831  }
832 
833  return result;
834 #else
835  return parser.Eval(&_spacetime[0]);
836 #endif
837 }
838 #else // LIBMESH_HAVE_FPARSER
839 template <typename Output>
840 inline
841 Output
842 ParsedFEMFunction<Output>::eval (char & /*parser*/,
843  const std::string & /*function_name*/,
844  unsigned int /*component_idx*/) const
845 {
846  libmesh_error_msg("ERROR: This functionality requires fparser!");
847  return Output(0);
848 }
849 #endif
850 
851 
852 } // namespace libMesh
853 
854 #endif // LIBMESH_PARSED_FEM_FUNCTION_H
virtual UniquePtr< FEMFunctionBase< Output > > clone() const libmesh_override
OStreamProxy err
void set_inline_value(const std::string &inline_var_name, Output newval)
Changes the value of an inline variable.
std::vector< bool > _need_var_grad
std::vector< Output > _initial_vals
ParsedFEMFunction(const System &sys, const std::string &expression, const std::vector< std::string > *additional_vars=libmesh_nullptr, const std::vector< Output > *initial_vals=libmesh_nullptr)
Constructor.
const std::string & variable_name(const unsigned int i) const
Definition: system.h:2134
void partial_reparse(const std::string &expression)
Number point_value(unsigned int var, const Point &p) const
Definition: fem_context.C:788
ImplicitSystem & sys
virtual Output component(const FEMContext &c, unsigned int i, const Point &p, Real time=0.) libmesh_override
std::vector< std::string > _additional_vars
const class libmesh_nullptr_t libmesh_nullptr
void get_element_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for interior finite element object for variable var for the largest dimension in the mesh...
Definition: fem_context.h:262
IterBase * end
Also have a polymorphic pointer to the end object, this prevents iterating past the end...
std::vector< std::string > _subexpressions
std::vector< Output > _spacetime
const unsigned int n_vars
Definition: tecplot_io.C:68
The libMesh namespace provides an interface to certain functionality in the library.
virtual unsigned int size() const libmesh_override
Definition: dense_vector.h:87
libmesh_assert(j)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
const std::vector< Point > & get_normals() const
Definition: fe_abstract.h:377
void reparse(const std::string &expression)
Re-parse with new expression.
Output get_inline_value(const std::string &inline_var_name) const
virtual Output operator()(const FEMContext &c, const Point &p, const Real time=0.) libmesh_override
ParsedFEMFunction provides support for FParser-based parsed functions in FEMSystem.
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
const std::vector< std::vector< OutputShape > > & get_phi() const
Definition: fe_base.h:208
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:61
Tensor point_hessian(unsigned int var, const Point &p) const
Definition: fem_context.C:886
PetscErrorCode Vec Mat libmesh_dbg_var(j)
std::vector< FunctionParserBase< Output > > parsers
std::vector< bool > _need_var_hess
void eval_args(const FEMContext &c, const Point &p, const Real time)
const std::vector< std::vector< OutputGradient > > & get_dphi() const
Definition: fe_base.h:216
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const std::vector< std::vector< OutputTensor > > & get_d2phi() const
Definition: fe_base.h:290
Output eval(FunctionParserBase< Output > &parser, const std::string &libmesh_dbg_var(function_name), unsigned int libmesh_dbg_var(component_idx)) const
std::size_t find_name(const std::string &varname, const std::string &expr) const
Gradient point_gradient(unsigned int var, const Point &p) const
Definition: fem_context.C:834
virtual ~ParsedFEMFunction()
Destructor.
Defines a dense vector for use in Finite Element-type computations.
virtual void init_context(const FEMContext &c) libmesh_override
Prepares a context object for use.
void get_side_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for edge/face (2D/3D) finite element object for variable var for the largest dimension in th...
Definition: fem_context.h:299
FEMFunctionBase is a base class from which users can derive in order to define "function-like" object...
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
const std::vector< Point > & get_xyz() const
Definition: fe_abstract.h:230
This class forms the foundation from which generic finite elements may be derived.
const std::string & expression()
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.