libMesh::ParsedFEMFunction< Output > Class Template Reference

Support for using parsed functions in FEMSystem. More...

#include <parsed_fem_function.h>

Inheritance diagram for libMesh::ParsedFEMFunction< Output >:

Public Member Functions

 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)
 
virtual ~ParsedFEMFunction ()
 
void reparse (const std::string &expression)
 
virtual void init_context (const FEMContext &c) libmesh_override
 
virtual UniquePtr< FEMFunctionBase< Output > > clone () const libmesh_override
 
virtual Output operator() (const FEMContext &c, const Point &p, const Real time=0.) libmesh_override
 
void operator() (const FEMContext &c, const Point &p, const Real time, DenseVector< Output > &output) libmesh_override
 
virtual Output component (const FEMContext &c, unsigned int i, const Point &p, Real time=0.) libmesh_override
 
const std::string & expression ()
 
Output get_inline_value (const std::string &inline_var_name) const
 
void set_inline_value (const std::string &inline_var_name, Output newval)
 
void operator() (const FEMContext &, const Point &p, DenseVector< Output > &output)
 

Protected Member Functions

void partial_reparse (const std::string &expression)
 
std::size_t find_name (const std::string &varname, const std::string &expr) const
 
void eval_args (const FEMContext &c, const Point &p, const Real time)
 
Output eval (FunctionParserBase< Output > &parser, const std::string &libmesh_dbg_var(function_name), unsigned int libmesh_dbg_var(component_idx)) const
 
Output eval (char &libmesh_dbg_var(parser), const std::string &libmesh_dbg_var(function_name), unsigned int libmesh_dbg_var(component_idx)) const
 

Private Attributes

const System_sys
 
std::string _expression
 
std::vector< std::string > _subexpressions
 
unsigned int _n_vars
 
unsigned int _n_requested_vars
 
unsigned int _n_requested_grad_components
 
unsigned int _n_requested_hess_components
 
bool _requested_normals
 
std::vector< FunctionParserBase< Output > > parsers
 
std::vector< char > parsers
 
std::vector< Output > _spacetime
 
std::vector< bool > _need_var
 
std::vector< bool > _need_var_grad
 
std::vector< bool > _need_var_hess
 
std::string variables
 
std::vector< std::string > _additional_vars
 
std::vector< Output > _initial_vals
 

Detailed Description

template<typename Output = Number>
class libMesh::ParsedFEMFunction< Output >

Support for using parsed functions in FEMSystem.

ParsedFEMFunction provides support for FParser-based parsed functions in FEMSystem.

Author
Roy Stogner
Date
2014

Definition at line 54 of file parsed_fem_function.h.

Constructor & Destructor Documentation

template<typename Output>
libMesh::ParsedFEMFunction< Output >::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 
)
inlineexplicit

Constructor.

Definition at line 204 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::reparse().

Referenced by libMesh::ParsedFEMFunction< Output >::clone().

207  :
208  _sys(sys),
209  _expression (), // overridden by parse()
210  _n_vars(sys.n_vars()),
214  _requested_normals(false),
215  _need_var(_n_vars, false),
216  _need_var_grad(_n_vars*LIBMESH_DIM, false),
217 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
218  _need_var_hess(_n_vars*LIBMESH_DIM*LIBMESH_DIM, false),
219 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
220  _additional_vars (additional_vars ? *additional_vars : std::vector<std::string>()),
221  _initial_vals (initial_vals ? *initial_vals : std::vector<Output>())
222 {
223  this->reparse(expression);
224 }
std::vector< bool > _need_var_grad
std::vector< Output > _initial_vals
ImplicitSystem & sys
std::vector< std::string > _additional_vars
void reparse(const std::string &expression)
std::vector< bool > _need_var_hess
unsigned int n_vars() const
Definition: system.h:2075
const std::string & expression()
template<typename Output = Number>
virtual libMesh::ParsedFEMFunction< Output >::~ParsedFEMFunction ( )
inlinevirtual

Destructor.

Definition at line 70 of file parsed_fem_function.h.

70 {}

Member Function Documentation

template<typename Output >
UniquePtr< FEMFunctionBase< Output > > libMesh::ParsedFEMFunction< Output >::clone ( ) const
inlinevirtual

Returns a new copy of the function. The new copy should be as ``deep'' as necessary to allow independent destruction and simultaneous evaluations of the copies in different threads.

Implements libMesh::FEMFunctionBase< Output >.

Definition at line 391 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::_additional_vars, libMesh::ParsedFEMFunction< Output >::_expression, libMesh::ParsedFEMFunction< Output >::_initial_vals, libMesh::ParsedFEMFunction< Output >::_sys, and libMesh::ParsedFEMFunction< Output >::ParsedFEMFunction().

Referenced by libMesh::ParsedFEMFunction< T >::~ParsedFEMFunction().

392 {
393  return UniquePtr<FEMFunctionBase<Output> >
395 }
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)
std::vector< std::string > _additional_vars
template<typename Output >
Output libMesh::ParsedFEMFunction< Output >::component ( const FEMContext c,
unsigned int  i,
const Point p,
Real  time = 0. 
)
inlinevirtual
Returns
the vector component i at coordinate p and time time.

Reimplemented from libMesh::FEMFunctionBase< Output >.

Definition at line 433 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::eval(), libMesh::ParsedFEMFunction< Output >::eval_args(), and libMesh::ParsedFEMFunction< Output >::parsers.

Referenced by libMesh::ParsedFEMFunction< T >::~ParsedFEMFunction().

437 {
438  eval_args(c, p, time);
439 
440  libmesh_assert_less (i, parsers.size());
441  return eval(parsers[i], "f", i);
442 }
std::vector< FunctionParserBase< Output > > parsers
void eval_args(const FEMContext &c, const Point &p, const Real time)
Output eval(FunctionParserBase< Output > &parser, const std::string &libmesh_dbg_var(function_name), unsigned int libmesh_dbg_var(component_idx)) const
template<typename Output>
Output libMesh::ParsedFEMFunction< Output >::eval ( FunctionParserBase< Output > &  parser,
const std::string &  libmesh_dbg_varfunction_name,
unsigned int   libmesh_dbg_varcomponent_idx 
) const
inlineprotected

Definition at line 806 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::_spacetime, and libMesh::err.

Referenced by libMesh::ParsedFEMFunction< Output >::component(), libMesh::ParsedFEMFunction< T >::expression(), libMesh::ParsedFEMFunction< Output >::get_inline_value(), and libMesh::ParsedFEMFunction< Output >::operator()().

809 {
810 #ifndef NDEBUG
811  Output result = parser.Eval(&_spacetime[0]);
812  int error_code = parser.EvalError();
813  if (error_code)
814  {
815  libMesh::err << "ERROR: FunctionParser is unable to evaluate component "
816  << component_idx
817  << " of expression '"
818  << function_name
819  << "' with arguments:\n";
820  for (std::size_t j=0; j<_spacetime.size(); ++j)
821  libMesh::err << '\t' << _spacetime[j] << '\n';
822  libMesh::err << '\n';
823 
824  // Currently no API to report error messages, we'll do it manually
825  std::string error_message = "Reason: ";
826 
827  switch (error_code)
828  {
829  case 1:
830  error_message += "Division by zero";
831  break;
832  case 2:
833  error_message += "Square Root error (negative value)";
834  break;
835  case 3:
836  error_message += "Log error (negative value)";
837  break;
838  case 4:
839  error_message += "Trigonometric error (asin or acos of illegal value)";
840  break;
841  case 5:
842  error_message += "Maximum recursion level reached";
843  break;
844  default:
845  error_message += "Unknown";
846  break;
847  }
848  libmesh_error_msg(error_message);
849  }
850 
851  return result;
852 #else
853  return parser.Eval(&_spacetime[0]);
854 #endif
855 }
std::vector< Output > _spacetime
OStreamProxy err(std::cerr)
template<typename Output>
Output libMesh::ParsedFEMFunction< Output >::eval ( char &  libmesh_dbg_varparser,
const std::string &  libmesh_dbg_varfunction_name,
unsigned int   libmesh_dbg_varcomponent_idx 
) const
inlineprotected

Definition at line 860 of file parsed_fem_function.h.

863 {
864  libmesh_error_msg("ERROR: This functionality requires fparser!");
865  return Output(0);
866 }
template<typename Output >
void libMesh::ParsedFEMFunction< Output >::eval_args ( const FEMContext c,
const Point p,
const Real  time 
)
inlineprotected

Definition at line 677 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::_n_requested_grad_components, libMesh::ParsedFEMFunction< Output >::_n_requested_hess_components, libMesh::ParsedFEMFunction< Output >::_n_vars, libMesh::ParsedFEMFunction< Output >::_need_var, libMesh::ParsedFEMFunction< Output >::_need_var_grad, libMesh::ParsedFEMFunction< Output >::_need_var_hess, libMesh::ParsedFEMFunction< Output >::_requested_normals, libMesh::ParsedFEMFunction< Output >::_spacetime, libMesh::FEAbstract::get_normals(), libMesh::FEMContext::get_side_fe(), libMesh::FEAbstract::get_xyz(), libMesh::libmesh_assert(), libMesh::FEMContext::point_gradient(), libMesh::FEMContext::point_hessian(), and libMesh::FEMContext::point_value().

Referenced by libMesh::ParsedFEMFunction< Output >::component(), libMesh::ParsedFEMFunction< T >::expression(), and libMesh::ParsedFEMFunction< Output >::operator()().

680 {
681  _spacetime[0] = p(0);
682 #if LIBMESH_DIM > 1
683  _spacetime[1] = p(1);
684 #endif
685 #if LIBMESH_DIM > 2
686  _spacetime[2] = p(2);
687 #endif
688  _spacetime[LIBMESH_DIM] = time;
689 
690  unsigned int request_index = 0;
691  for (unsigned int v=0; v != _n_vars; ++v)
692  {
693  if (!_need_var[v])
694  continue;
695 
696  c.point_value(v, p, _spacetime[LIBMESH_DIM+1+request_index]);
697  request_index++;
698  }
699 
701  for (unsigned int v=0; v != _n_vars; ++v)
702  {
703  if (!_need_var_grad[v*LIBMESH_DIM]
704 #if LIBMESH_DIM > 1
705  && !_need_var_grad[v*LIBMESH_DIM+1]
706 #if LIBMESH_DIM > 2
707  && !_need_var_grad[v*LIBMESH_DIM+2]
708 #endif
709 #endif
710  )
711  continue;
712 
713  Gradient g;
714  c.point_gradient(v, p, g);
715 
716  for (unsigned int d=0; d != LIBMESH_DIM; ++d)
717  {
718  if (!_need_var_grad[v*LIBMESH_DIM+d])
719  continue;
720 
721  _spacetime[LIBMESH_DIM+1+request_index] = g(d);
722  request_index++;
723  }
724  }
725 
726 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
728  for (unsigned int v=0; v != _n_vars; ++v)
729  {
730  if (!_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM]
731 #if LIBMESH_DIM > 1
732  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+1]
733  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+2]
734  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+3]
735 #if LIBMESH_DIM > 2
736  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+4]
737  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+5]
738  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+6]
739  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+7]
740  && !_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+8]
741 #endif
742 #endif
743  )
744  continue;
745 
746  Tensor h;
747  c.point_hessian(v, p, h);
748 
749  for (unsigned int d1=0; d1 != LIBMESH_DIM; ++d1)
750  for (unsigned int d2=0; d2 != LIBMESH_DIM; ++d2)
751  {
752  if (!_need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+d1*LIBMESH_DIM+d2])
753  continue;
754 
755  _spacetime[LIBMESH_DIM+1+request_index] = h(d1,d2);
756  request_index++;
757  }
758  }
759 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
760 
761  if (_requested_normals)
762  {
763  FEBase * side_fe;
764  c.get_side_fe(0, side_fe);
765 
766  const std::vector<Point> & normals = side_fe->get_normals();
767 
768  const std::vector<Point> & xyz = side_fe->get_xyz();
769 
770  libmesh_assert_equal_to(normals.size(), xyz.size());
771 
772  // We currently only support normals at quadrature points!
773 #ifndef NDEBUG
774  bool at_quadrature_point = false;
775 #endif
776  for (std::size_t qp = 0; qp != normals.size(); ++qp)
777  {
778  if (p == xyz[qp])
779  {
780  const Point & n = normals[qp];
781  for (unsigned int d=0; d != LIBMESH_DIM; ++d)
782  {
783  _spacetime[LIBMESH_DIM+1+request_index] = n(d);
784  request_index++;
785  }
786 #ifndef NDEBUG
787  at_quadrature_point = true;
788 #endif
789  break;
790  }
791  }
792 
793  libmesh_assert(at_quadrature_point);
794  }
795 
796  // The remaining locations in _spacetime are currently set at construction
797  // but could potentially be made dynamic
798 }
std::vector< bool > _need_var_grad
std::vector< Output > _spacetime
libmesh_assert(j)
NumberVectorValue Gradient
FEGenericBase< Real > FEBase
std::vector< bool > _need_var_hess
NumberTensorValue Tensor
template<typename Output = Number>
const std::string& libMesh::ParsedFEMFunction< Output >::expression ( )
inline
template<typename Output >
std::size_t libMesh::ParsedFEMFunction< Output >::find_name ( const std::string &  varname,
const std::string &  expr 
) const
inlineprotected

Definition at line 651 of file parsed_fem_function.h.

Referenced by libMesh::ParsedFEMFunction< T >::expression(), libMesh::ParsedFEMFunction< Output >::get_inline_value(), libMesh::ParsedFEMFunction< Output >::reparse(), and libMesh::ParsedFEMFunction< Output >::set_inline_value().

653 {
654  const std::size_t namesize = varname.size();
655  std::size_t varname_i = expr.find(varname);
656 
657  while ((varname_i != std::string::npos) &&
658  (((varname_i > 0) &&
659  (std::isalnum(expr[varname_i-1]) ||
660  (expr[varname_i-1] == '_'))) ||
661  ((varname_i+namesize < expr.size()) &&
662  (std::isalnum(expr[varname_i+namesize]) ||
663  (expr[varname_i+namesize] == '_')))))
664  {
665  varname_i = expr.find(varname, varname_i+1);
666  }
667 
668  return varname_i;
669 }
template<typename Output >
Output libMesh::ParsedFEMFunction< Output >::get_inline_value ( const std::string &  inline_var_name) const
inline
Returns
the value of an inline variable. Will only be correct if the inline variable value is independent of input variables, if the inline variable is not redefined within any subexpression, and if the inline variable takes the same value within any subexpressions where it appears.

Definition at line 447 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::_subexpressions, libMesh::ParsedFEMFunction< Output >::eval(), libMesh::ParsedFEMFunction< Output >::find_name(), libMesh::libmesh_assert(), libMesh::Real, and libMesh::ParsedFEMFunction< Output >::variables.

Referenced by libMesh::ParsedFEMFunction< T >::expression(), and libMesh::ParsedFEMFunctionParameter< T >::get().

448 {
449  libmesh_assert_greater (_subexpressions.size(), 0);
450 
451 #ifndef NDEBUG
452  bool found_var_name = false;
453 #endif
454  Output old_var_value(0.);
455 
456  for (std::size_t s=0; s != _subexpressions.size(); ++s)
457  {
458  const std::string & subexpression = _subexpressions[s];
459  const std::size_t varname_i =
460  find_name(inline_var_name, subexpression);
461  if (varname_i == std::string::npos)
462  continue;
463 
464  const std::size_t assignment_i =
465  subexpression.find(":", varname_i+1);
466 
467  libmesh_assert_not_equal_to(assignment_i, std::string::npos);
468 
469  libmesh_assert_equal_to(subexpression[assignment_i+1], '=');
470  for (unsigned int i = varname_i+1; i != assignment_i; ++i)
471  libmesh_assert_equal_to(subexpression[i], ' ');
472 
473  std::size_t end_assignment_i =
474  subexpression.find(";", assignment_i+1);
475 
476  libmesh_assert_not_equal_to(end_assignment_i, std::string::npos);
477 
478  std::string new_subexpression =
479  subexpression.substr(0, end_assignment_i+1) +
480  inline_var_name;
481 
482 #ifdef LIBMESH_HAVE_FPARSER
483  // Parse and evaluate the new subexpression.
484  // Add the same constants as we used originally.
485  FunctionParserBase<Output> fp;
486  fp.AddConstant("NaN", std::numeric_limits<Real>::quiet_NaN());
487  fp.AddConstant("pi", std::acos(Real(-1)));
488  fp.AddConstant("e", std::exp(Real(1)));
489  if (fp.Parse(new_subexpression, variables) != -1) // -1 for success
490  libmesh_error_msg
491  ("ERROR: FunctionParser is unable to parse modified expression: "
492  << new_subexpression << '\n' << fp.ErrorMsg());
493 
494  Output new_var_value = this->eval(fp, new_subexpression, 0);
495 #ifdef NDEBUG
496  return new_var_value;
497 #else
498  if (found_var_name)
499  {
500  libmesh_assert_equal_to(old_var_value, new_var_value);
501  }
502  else
503  {
504  old_var_value = new_var_value;
505  found_var_name = true;
506  }
507 #endif
508 
509 #else
510  libmesh_error_msg("ERROR: This functionality requires fparser!");
511 #endif
512  }
513 
514  libmesh_assert(found_var_name);
515  return old_var_value;
516 }
std::vector< std::string > _subexpressions
libmesh_assert(j)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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
template<typename Output >
void libMesh::ParsedFEMFunction< Output >::init_context ( const FEMContext c)
inlinevirtual

Prepares a context object for use.

Reimplemented from libMesh::FEMFunctionBase< Output >.

Definition at line 359 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::_n_requested_grad_components, libMesh::ParsedFEMFunction< Output >::_n_requested_hess_components, libMesh::ParsedFEMFunction< Output >::_n_requested_vars, libMesh::ParsedFEMFunction< Output >::_n_vars, libMesh::ParsedFEMFunction< Output >::_requested_normals, libMesh::FEGenericBase< OutputType >::get_d2phi(), libMesh::FEGenericBase< OutputType >::get_dphi(), libMesh::FEMContext::get_element_fe(), libMesh::FEAbstract::get_normals(), libMesh::FEGenericBase< OutputType >::get_phi(), libMesh::FEMContext::get_side_fe(), and libMesh::FEAbstract::get_xyz().

Referenced by libMesh::ParsedFEMFunction< T >::~ParsedFEMFunction().

360 {
361  for (unsigned int v=0; v != _n_vars; ++v)
362  {
363  FEBase * elem_fe;
364  c.get_element_fe(v, elem_fe);
365  if (_n_requested_vars)
366  elem_fe->get_phi();
368  elem_fe->get_dphi();
369 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
371  elem_fe->get_d2phi();
372 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
373  }
374 
375  if (_requested_normals)
376  {
377  FEBase * side_fe;
378  c.get_side_fe(0, side_fe);
379 
380  side_fe->get_normals();
381 
382  // FIXME: this is a hack to support normals at quadrature
383  // points; we don't support normals elsewhere.
384  side_fe->get_xyz();
385  }
386 }
FEGenericBase< Real > FEBase
template<typename Output >
Output libMesh::ParsedFEMFunction< Output >::operator() ( const FEMContext c,
const Point p,
const Real  time = 0. 
)
inlinevirtual
Returns
the scalar value at coordinate p and time time, which defaults to zero. Purely virtual, so you have to overload it. Note that this cannot be a const method, check MeshFunction.

Implements libMesh::FEMFunctionBase< Output >.

Definition at line 400 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::eval(), libMesh::ParsedFEMFunction< Output >::eval_args(), and libMesh::ParsedFEMFunction< Output >::parsers.

Referenced by libMesh::ParsedFEMFunction< T >::~ParsedFEMFunction().

403 {
404  eval_args(c, p, time);
405 
406  return eval(parsers[0], "f", 0);
407 }
std::vector< FunctionParserBase< Output > > parsers
void eval_args(const FEMContext &c, const Point &p, const Real time)
Output eval(FunctionParserBase< Output > &parser, const std::string &libmesh_dbg_var(function_name), unsigned int libmesh_dbg_var(component_idx)) const
template<typename Output >
void libMesh::FEMFunctionBase< Output >::operator() ( const FEMContext context,
const Point p,
DenseVector< Output > &  output 
)
inlineinherited

Return function for vectors. Returns in output the values of the data at the coordinate p.

Definition at line 144 of file fem_function_base.h.

References libMesh::FEMFunctionBase< Output >::operator()().

147 {
148  // Call the time-dependent function with t=0.
149  this->operator()(context, p, 0., output);
150 }
virtual Output operator()(const FEMContext &, const Point &p, const Real time=0.)=0
template<typename Output>
void libMesh::ParsedFEMFunction< Output >::operator() ( const FEMContext c,
const Point p,
const Real  time,
DenseVector< Output > &  output 
)
inlinevirtual

Return function for vectors. Returns in output the values of the data at the coordinate p and for time time.

Implements libMesh::FEMFunctionBase< Output >.

Definition at line 414 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::eval(), libMesh::ParsedFEMFunction< Output >::eval_args(), libMesh::ParsedFEMFunction< Output >::parsers, and libMesh::DenseVector< T >::size().

418 {
419  eval_args(c, p, time);
420 
421  unsigned int size = output.size();
422 
423  libmesh_assert_equal_to (size, parsers.size());
424 
425  for (unsigned int i=0; i != size; ++i)
426  output(i) = eval(parsers[i], "f", i);
427 }
std::vector< FunctionParserBase< Output > > parsers
void eval_args(const FEMContext &c, const Point &p, const Real time)
Output eval(FunctionParserBase< Output > &parser, const std::string &libmesh_dbg_var(function_name), unsigned int libmesh_dbg_var(component_idx)) const
template<typename Output >
void libMesh::ParsedFEMFunction< Output >::partial_reparse ( const std::string &  expression)
inlineprotected

Definition at line 587 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::_expression, libMesh::ParsedFEMFunction< Output >::_subexpressions, end, libMesh::ParsedFEMFunction< Output >::expression(), libMesh::ParsedFEMFunction< Output >::parsers, libMesh::Real, and libMesh::ParsedFEMFunction< Output >::variables.

Referenced by libMesh::ParsedFEMFunction< T >::expression(), libMesh::ParsedFEMFunction< Output >::reparse(), and libMesh::ParsedFEMFunction< Output >::set_inline_value().

588 {
590  _subexpressions.clear();
591  parsers.clear();
592 
593  size_t nextstart = 0, end = 0;
594 
595  while (end != std::string::npos)
596  {
597  // If we're past the end of the string, we can't make any more
598  // subparsers
599  if (nextstart >= expression.size())
600  break;
601 
602  // If we're at the start of a brace delimited section, then we
603  // parse just that section:
604  if (expression[nextstart] == '{')
605  {
606  nextstart++;
607  end = expression.find('}', nextstart);
608  }
609  // otherwise we parse the whole thing
610  else
611  end = std::string::npos;
612 
613  // We either want the whole end of the string (end == npos) or
614  // a substring in the middle.
615  _subexpressions.push_back
616  (expression.substr(nextstart, (end == std::string::npos) ?
617  std::string::npos : end - nextstart));
618 
619  // fparser can crash on empty expressions
620  if (_subexpressions.back().empty())
621  libmesh_error_msg("ERROR: FunctionParser is unable to parse empty expression.\n");
622 
623 
624 #ifdef LIBMESH_HAVE_FPARSER
625  // Parse (and optimize if possible) the subexpression.
626  // Add some basic constants, to Real precision.
627  FunctionParserBase<Output> fp;
628  fp.AddConstant("NaN", std::numeric_limits<Real>::quiet_NaN());
629  fp.AddConstant("pi", std::acos(Real(-1)));
630  fp.AddConstant("e", std::exp(Real(1)));
631  if (fp.Parse(_subexpressions.back(), variables) != -1) // -1 for success
632  libmesh_error_msg
633  ("ERROR: FunctionParser is unable to parse expression: "
634  << _subexpressions.back() << '\n' << fp.ErrorMsg());
635  fp.Optimize();
636  parsers.push_back(fp);
637 #else
638  libmesh_error_msg("ERROR: This functionality requires fparser!");
639 #endif
640 
641  // If at end, use nextstart=maxSize. Else start at next
642  // character.
643  nextstart = (end == std::string::npos) ?
644  std::string::npos : end + 1;
645  }
646 }
IterBase * end
std::vector< std::string > _subexpressions
std::vector< FunctionParserBase< Output > > parsers
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const std::string & expression()
template<typename Output >
void libMesh::ParsedFEMFunction< Output >::reparse ( const std::string &  expression)
inline

Definition at line 230 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::_additional_vars, libMesh::ParsedFEMFunction< Output >::_initial_vals, libMesh::ParsedFEMFunction< Output >::_n_requested_grad_components, libMesh::ParsedFEMFunction< Output >::_n_requested_hess_components, libMesh::ParsedFEMFunction< Output >::_n_requested_vars, libMesh::ParsedFEMFunction< Output >::_n_vars, libMesh::ParsedFEMFunction< Output >::_need_var, libMesh::ParsedFEMFunction< Output >::_need_var_grad, libMesh::ParsedFEMFunction< Output >::_need_var_hess, libMesh::ParsedFEMFunction< Output >::_requested_normals, libMesh::ParsedFEMFunction< Output >::_spacetime, libMesh::ParsedFEMFunction< Output >::_sys, libMesh::ParsedFEMFunction< Output >::find_name(), libMesh::ParsedFEMFunction< Output >::partial_reparse(), libMesh::System::variable_name(), and libMesh::ParsedFEMFunction< Output >::variables.

Referenced by libMesh::ParsedFEMFunction< Output >::ParsedFEMFunction(), and libMesh::ParsedFEMFunction< T >::~ParsedFEMFunction().

231 {
232  variables = "x";
233 #if LIBMESH_DIM > 1
234  variables += ",y";
235 #endif
236 #if LIBMESH_DIM > 2
237  variables += ",z";
238 #endif
239  variables += ",t";
240 
241  for (unsigned int v=0; v != _n_vars; ++v)
242  {
243  const std::string & varname = _sys.variable_name(v);
244  std::size_t varname_i = find_name(varname, expression);
245 
246  // If we didn't find our variable name then let's go to the
247  // next.
248  if (varname_i == std::string::npos)
249  continue;
250 
251  _need_var[v] = true;
252  variables += ',';
253  variables += varname;
255  }
256 
257  for (unsigned int v=0; v != _n_vars; ++v)
258  {
259  const std::string & varname = _sys.variable_name(v);
260 
261  for (unsigned int d=0; d != LIBMESH_DIM; ++d)
262  {
263  std::string gradname = std::string("grad_") +
264  "xyz"[d] + '_' + varname;
265  std::size_t gradname_i = find_name(gradname, expression);
266 
267  // If we didn't find that gradient component of our
268  // variable name then let's go to the next.
269  if (gradname_i == std::string::npos)
270  continue;
271 
272  _need_var_grad[v*LIBMESH_DIM+d] = true;
273  variables += ',';
274  variables += gradname;
276  }
277  }
278 
279 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
280  for (unsigned int v=0; v != _n_vars; ++v)
281  {
282  const std::string & varname = _sys.variable_name(v);
283 
284  for (unsigned int d1=0; d1 != LIBMESH_DIM; ++d1)
285  for (unsigned int d2=0; d2 != LIBMESH_DIM; ++d2)
286  {
287  std::string hessname = std::string("hess_") +
288  "xyz"[d1] + "xyz"[d2] + '_' + varname;
289  std::size_t hessname_i = find_name(hessname, expression);
290 
291  // If we didn't find that hessian component of our
292  // variable name then let's go to the next.
293  if (hessname_i == std::string::npos)
294  continue;
295 
296  _need_var_hess[v*LIBMESH_DIM*LIBMESH_DIM+d1*LIBMESH_DIM+d2]
297  = true;
298  variables += ',';
299  variables += hessname;
301  }
302  }
303 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
304 
305  {
306  std::size_t nx_i = find_name("n_x", expression);
307  std::size_t ny_i = find_name("n_y", expression);
308  std::size_t nz_i = find_name("n_z", expression);
309 
310  // If we found any requests for normal components then we'll
311  // compute normals
312  if (nx_i != std::string::npos ||
313  ny_i != std::string::npos ||
314  nz_i != std::string::npos)
315  {
316  _requested_normals = true;
317  variables += ',';
318  variables += "n_x";
319  if (LIBMESH_DIM > 1)
320  {
321  variables += ',';
322  variables += "n_y";
323  }
324  if (LIBMESH_DIM > 2)
325  {
326  variables += ',';
327  variables += "n_z";
328  }
329  }
330  }
331 
332  _spacetime.resize
333  (LIBMESH_DIM + 1 + _n_requested_vars +
335  (_requested_normals ? LIBMESH_DIM : 0) +
336  _additional_vars.size());
337 
338  // If additional vars were passed, append them to the string
339  // that we send to the function parser. Also add them to the
340  // end of our spacetime vector
341  unsigned int offset = LIBMESH_DIM + 1 + _n_requested_vars +
343 
344  for (std::size_t i=0; i < _additional_vars.size(); ++i)
345  {
346  variables += "," + _additional_vars[i];
347  // Initialize extra variables to the vector passed in or zero
348  // Note: The initial_vals vector can be shorter than the additional_vars vector
349  _spacetime[offset + i] =
350  (i < _initial_vals.size()) ? _initial_vals[i] : 0;
351  }
352 
354 }
std::vector< bool > _need_var_grad
std::vector< Output > _initial_vals
const std::string & variable_name(const unsigned int i) const
Definition: system.h:2123
void partial_reparse(const std::string &expression)
std::vector< std::string > _additional_vars
std::vector< Output > _spacetime
std::vector< bool > _need_var_hess
std::size_t find_name(const std::string &varname, const std::string &expr) const
const std::string & expression()
template<typename Output>
void libMesh::ParsedFEMFunction< Output >::set_inline_value ( const std::string &  inline_var_name,
Output  newval 
)
inline

Changes the value of an inline variable. Forever after the variable value will take the given constant, independent of input variables, in every subexpression where it is already defined. Currently only works if the inline variable is not redefined within any one subexpression.

Definition at line 521 of file parsed_fem_function.h.

References libMesh::ParsedFEMFunction< Output >::_subexpressions, libMesh::ParsedFEMFunction< Output >::find_name(), libMesh::libmesh_assert(), and libMesh::ParsedFEMFunction< Output >::partial_reparse().

Referenced by libMesh::ParsedFEMFunction< T >::expression(), and libMesh::ParsedFEMFunctionParameter< T >::set().

523 {
524  libmesh_assert_greater (_subexpressions.size(), 0);
525 
526 #ifndef NDEBUG
527  bool found_var_name = false;
528 #endif
529  for (std::size_t s=0; s != _subexpressions.size(); ++s)
530  {
531  const std::string & subexpression = _subexpressions[s];
532  const std::size_t varname_i =
533  find_name(inline_var_name, subexpression);
534  if (varname_i == std::string::npos)
535  continue;
536 
537 #ifndef NDEBUG
538  found_var_name = true;
539 #endif
540 
541  const std::size_t assignment_i =
542  subexpression.find(":", varname_i+1);
543 
544  libmesh_assert_not_equal_to(assignment_i, std::string::npos);
545 
546  libmesh_assert_equal_to(subexpression[assignment_i+1], '=');
547  for (unsigned int i = varname_i+1; i != assignment_i; ++i)
548  libmesh_assert_equal_to(subexpression[i], ' ');
549 
550  std::size_t end_assignment_i =
551  subexpression.find(";", assignment_i+1);
552 
553  libmesh_assert_not_equal_to(end_assignment_i, std::string::npos);
554 
555  std::ostringstream new_subexpression;
556  new_subexpression << subexpression.substr(0, assignment_i+2)
557  << std::setprecision(std::numeric_limits<Output>::digits10+2)
558 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
559  << '(' << newval.real() << '+'
560  << newval.imag() << 'i' << ')'
561 #else
562  << newval
563 #endif
564  << subexpression.substr(end_assignment_i,
565  std::string::npos);
566  _subexpressions[s] = new_subexpression.str();
567  }
568 
569  libmesh_assert(found_var_name);
570 
571  std::string new_expression;
572 
573  for (std::size_t s=0; s != _subexpressions.size(); ++s)
574  {
575  new_expression += '{';
576  new_expression += _subexpressions[s];
577  new_expression += '}';
578  }
579 
580  this->partial_reparse(new_expression);
581 }
void partial_reparse(const std::string &expression)
std::vector< std::string > _subexpressions
libmesh_assert(j)
std::size_t find_name(const std::string &varname, const std::string &expr) const

Member Data Documentation

template<typename Output = Number>
std::vector<std::string> libMesh::ParsedFEMFunction< Output >::_additional_vars
private
template<typename Output = Number>
std::string libMesh::ParsedFEMFunction< Output >::_expression
private
template<typename Output = Number>
std::vector<Output> libMesh::ParsedFEMFunction< Output >::_initial_vals
private
template<typename Output = Number>
unsigned int libMesh::ParsedFEMFunction< Output >::_n_requested_grad_components
private
template<typename Output = Number>
unsigned int libMesh::ParsedFEMFunction< Output >::_n_requested_hess_components
private
template<typename Output = Number>
unsigned int libMesh::ParsedFEMFunction< Output >::_n_requested_vars
private
template<typename Output = Number>
std::vector<bool> libMesh::ParsedFEMFunction< Output >::_need_var
private
template<typename Output = Number>
std::vector<bool> libMesh::ParsedFEMFunction< Output >::_need_var_grad
private
template<typename Output = Number>
std::vector<bool> libMesh::ParsedFEMFunction< Output >::_need_var_hess
private
template<typename Output = Number>
bool libMesh::ParsedFEMFunction< Output >::_requested_normals
private
template<typename Output = Number>
std::vector<Output> libMesh::ParsedFEMFunction< Output >::_spacetime
private
template<typename Output = Number>
std::vector<std::string> libMesh::ParsedFEMFunction< Output >::_subexpressions
private
template<typename Output = Number>
const System& libMesh::ParsedFEMFunction< Output >::_sys
private
template<typename Output = Number>
std::vector<FunctionParserBase<Output> > libMesh::ParsedFEMFunction< Output >::parsers
private
template<typename Output = Number>
std::vector<char> libMesh::ParsedFEMFunction< Output >::parsers
private

Definition at line 175 of file parsed_fem_function.h.


The documentation for this class was generated from the following file: