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. All overridden virtual functions are documented in fem_function_base.h.

Author
Roy Stogner
Date
2014

Definition at line 55 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 186 of file parsed_fem_function.h.

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

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

189  :
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 }
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:2084
const std::string & expression()
template<typename Output = Number>
virtual libMesh::ParsedFEMFunction< Output >::~ParsedFEMFunction ( )
inlinevirtual

Destructor.

Definition at line 71 of file parsed_fem_function.h.

71 {}

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 373 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().

374 {
375  return UniquePtr<FEMFunctionBase<Output> >
377 }
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 context,
unsigned int  i,
const Point p,
Real  time = 0. 
)
inlinevirtual
Returns
The vector component i at coordinate p and time time.
Note
Subclasses aren't required to override this, since the default implementation is based on the full vector evaluation, which is often correct.
Subclasses are recommended to override this, since the default implementation is based on a vector evaluation, which is usually unnecessarily inefficient.

Reimplemented from libMesh::FEMFunctionBase< Output >.

Definition at line 415 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().

419 {
420  eval_args(c, p, time);
421 
422  libmesh_assert_less (i, parsers.size());
423  return eval(parsers[i], "f", i);
424 }
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 788 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()().

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 }
std::vector< Output > _spacetime
OStreamProxy err(std::cerr)
template<typename Output = Number>
Output libMesh::ParsedFEMFunction< Output >::eval ( char &  libmesh_dbg_varparser,
const std::string &  libmesh_dbg_varfunction_name,
unsigned int   libmesh_dbg_varcomponent_idx 
) const
inlineprotected
template<typename Output >
void libMesh::ParsedFEMFunction< Output >::eval_args ( const FEMContext c,
const Point p,
const Real  time 
)
inlineprotected

Definition at line 659 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()().

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 }
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 633 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().

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 }
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.
Note
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 429 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().

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 }
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 )
inlinevirtual

Prepares a context object for use.

Most problems will want to reimplement this for efficiency, in order to call FE::get_*() as their particular function requires.

Reimplemented from libMesh::FEMFunctionBase< Output >.

Definition at line 341 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().

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 }
FEGenericBase< Real > FEBase
template<typename Output >
Output libMesh::ParsedFEMFunction< Output >::operator() ( const FEMContext ,
const Point p,
const Real  time = 0. 
)
inlinevirtual
Returns
The scalar function value at coordinate p and time time, which defaults to zero.

Pure virtual, so you have to override it.

Implements libMesh::FEMFunctionBase< Output >.

Definition at line 382 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().

385 {
386  eval_args(c, p, time);
387 
388  return eval(parsers[0], "f", 0);
389 }
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 >::operator() ( const FEMContext ,
const Point p,
const Real  time,
DenseVector< Output > &  output 
)
inlinevirtual

Evaluation function for time-dependent vector-valued functions. Sets output values in the passed-in output DenseVector.

Pure virtual, so you have to override it.

Implements libMesh::FEMFunctionBase< Output >.

Definition at line 396 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().

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 }
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

Evaluation function for time-independent vector-valued functions. Sets output values in the passed-in output DenseVector.

Definition at line 143 of file fem_function_base.h.

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

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

Definition at line 569 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().

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 }
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

Re-parse with new expression.

Definition at line 212 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().

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 
336 }
std::vector< bool > _need_var_grad
std::vector< Output > _initial_vals
const std::string & variable_name(const unsigned int i) const
Definition: system.h:2132
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.

Note
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 503 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().

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 }
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 157 of file parsed_fem_function.h.


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