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=nullptr, const std::vector< Output > *initial_vals=nullptr)
 
ParsedFEMFunctionoperator= (const ParsedFEMFunction &)=delete
 
ParsedFEMFunctionoperator= (ParsedFEMFunction &&)=delete
 
 ParsedFEMFunction (const ParsedFEMFunction &)=default
 
 ParsedFEMFunction (ParsedFEMFunction &&)=default
 
virtual ~ParsedFEMFunction ()=default
 
void reparse (const std::string &expression)
 
virtual void init_context (const FEMContext &c) override
 
virtual std::unique_ptr< FEMFunctionBase< Output > > clone () const override
 
virtual Output operator() (const FEMContext &c, const Point &p, const Real time=0.) override
 
void operator() (const FEMContext &c, const Point &p, const Real time, DenseVector< Output > &output) override
 
virtual Output component (const FEMContext &c, unsigned int i, const Point &p, Real time=0.) 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

◆ ParsedFEMFunction() [1/3]

template<typename Output>
libMesh::ParsedFEMFunction< Output >::ParsedFEMFunction ( const System sys,
const std::string &  expression,
const std::vector< std::string > *  additional_vars = nullptr,
const std::vector< Output > *  initial_vals = nullptr 
)
inlineexplicit

Constructor.

Definition at line 200 of file parsed_fem_function.h.

203  :
204  _sys(sys),
205  _expression (), // overridden by parse()
206  _n_vars(sys.n_vars()),
210  _requested_normals(false),
211  _need_var(_n_vars, false),
212  _need_var_grad(_n_vars*LIBMESH_DIM, false),
213 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
214  _need_var_hess(_n_vars*LIBMESH_DIM*LIBMESH_DIM, false),
215 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
216  _additional_vars (additional_vars ? *additional_vars : std::vector<std::string>()),
217  _initial_vals (initial_vals ? *initial_vals : std::vector<Output>())
218 {
219  this->reparse(expression);
220 }
std::vector< bool > _need_var_grad
std::vector< Output > _initial_vals
std::vector< std::string > _additional_vars
void reparse(const std::string &expression)
std::vector< bool > _need_var_hess
const std::string & expression()

◆ ParsedFEMFunction() [2/3]

template<typename Output = Number>
libMesh::ParsedFEMFunction< Output >::ParsedFEMFunction ( const ParsedFEMFunction< Output > &  )
default

The remaining 5 special functions can be safely defaulted.

Note
The underlying FunctionParserBase class has a copy constructor, so this class should be default-constructible. And, although FunctionParserBase's move constructor is deleted, this class should still be move-constructible because FunctionParserBase only appears in a vector.

◆ ParsedFEMFunction() [3/3]

template<typename Output = Number>
libMesh::ParsedFEMFunction< Output >::ParsedFEMFunction ( ParsedFEMFunction< Output > &&  )
default

◆ ~ParsedFEMFunction()

template<typename Output = Number>
virtual libMesh::ParsedFEMFunction< Output >::~ParsedFEMFunction ( )
virtualdefault

Member Function Documentation

◆ clone()

template<typename Output >
std::unique_ptr< FEMFunctionBase< Output > > libMesh::ParsedFEMFunction< Output >::clone ( ) const
inlineoverridevirtual
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 387 of file parsed_fem_function.h.

388 {
389  return std::unique_ptr<FEMFunctionBase<Output>>
391 }
std::vector< Output > _initial_vals
std::vector< std::string > _additional_vars
ParsedFEMFunction(const System &sys, const std::string &expression, const std::vector< std::string > *additional_vars=nullptr, const std::vector< Output > *initial_vals=nullptr)

◆ component()

template<typename Output >
Output libMesh::ParsedFEMFunction< Output >::component ( const FEMContext context,
unsigned int  i,
const Point p,
Real  time = 0. 
)
inlineoverridevirtual
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 429 of file parsed_fem_function.h.

433 {
434  eval_args(c, p, time);
435 
436  libmesh_assert_less (i, parsers.size());
437  return eval(parsers[i], "f", i);
438 }
Output eval(FunctionParserBase< Output > &parser, const std::string &libmesh_dbg_var(function_name), unsigned int libmesh_dbg_var(component_idx)) const
std::vector< FunctionParserBase< Output > > parsers
void eval_args(const FEMContext &c, const Point &p, const Real time)

◆ eval() [1/2]

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

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

◆ eval() [2/2]

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

◆ eval_args()

template<typename Output >
void libMesh::ParsedFEMFunction< Output >::eval_args ( const FEMContext c,
const Point p,
const Real  time 
)
inlineprotected

Definition at line 673 of file parsed_fem_function.h.

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

◆ expression()

template<typename Output = Number>
const std::string& libMesh::ParsedFEMFunction< Output >::expression ( )
inline

Definition at line 110 of file parsed_fem_function.h.

110 { return _expression; }

◆ find_name()

template<typename Output >
std::size_t libMesh::ParsedFEMFunction< Output >::find_name ( const std::string &  varname,
const std::string &  expr 
) const
inlineprotected

Definition at line 647 of file parsed_fem_function.h.

649 {
650  const std::size_t namesize = varname.size();
651  std::size_t varname_i = expr.find(varname);
652 
653  while ((varname_i != std::string::npos) &&
654  (((varname_i > 0) &&
655  (std::isalnum(expr[varname_i-1]) ||
656  (expr[varname_i-1] == '_'))) ||
657  ((varname_i+namesize < expr.size()) &&
658  (std::isalnum(expr[varname_i+namesize]) ||
659  (expr[varname_i+namesize] == '_')))))
660  {
661  varname_i = expr.find(varname, varname_i+1);
662  }
663 
664  return varname_i;
665 }

◆ get_inline_value()

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

Referenced by libMesh::ParsedFEMFunctionParameter< T >::get().

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

◆ init_context()

template<typename Output >
void libMesh::ParsedFEMFunction< Output >::init_context ( const FEMContext )
inlineoverridevirtual

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

356 {
357  for (unsigned int v=0; v != _n_vars; ++v)
358  {
359  FEBase * elem_fe;
360  c.get_element_fe(v, elem_fe);
361  if (_n_requested_vars)
362  elem_fe->get_phi();
364  elem_fe->get_dphi();
365 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
367  elem_fe->get_d2phi();
368 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
369  }
370 
371  if (_requested_normals)
372  {
373  FEBase * side_fe;
374  c.get_side_fe(0, side_fe);
375 
376  side_fe->get_normals();
377 
378  // FIXME: this is a hack to support normals at quadrature
379  // points; we don't support normals elsewhere.
380  side_fe->get_xyz();
381  }
382 }
FEGenericBase< Real > FEBase

◆ operator()() [1/3]

template<typename Output >
Output libMesh::ParsedFEMFunction< Output >::operator() ( const FEMContext ,
const Point p,
const Real  time = 0. 
)
inlineoverridevirtual
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 396 of file parsed_fem_function.h.

399 {
400  eval_args(c, p, time);
401 
402  return eval(parsers[0], "f", 0);
403 }
Output eval(FunctionParserBase< Output > &parser, const std::string &libmesh_dbg_var(function_name), unsigned int libmesh_dbg_var(component_idx)) const
std::vector< FunctionParserBase< Output > > parsers
void eval_args(const FEMContext &c, const Point &p, const Real time)

◆ operator()() [2/3]

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 145 of file fem_function_base.h.

148 {
149  // Call the time-dependent function with t=0.
150  this->operator()(context, p, 0., output);
151 }
virtual Output operator()(const FEMContext &, const Point &p, const Real time=0.)=0

◆ operator()() [3/3]

template<typename Output>
void libMesh::ParsedFEMFunction< Output >::operator() ( const FEMContext ,
const Point p,
const Real  time,
DenseVector< Output > &  output 
)
inlineoverridevirtual

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

414 {
415  eval_args(c, p, time);
416 
417  unsigned int size = output.size();
418 
419  libmesh_assert_equal_to (size, parsers.size());
420 
421  for (unsigned int i=0; i != size; ++i)
422  output(i) = eval(parsers[i], "f", i);
423 }
virtual unsigned int size() const override
Definition: dense_vector.h:92
Output eval(FunctionParserBase< Output > &parser, const std::string &libmesh_dbg_var(function_name), unsigned int libmesh_dbg_var(component_idx)) const
std::vector< FunctionParserBase< Output > > parsers
void eval_args(const FEMContext &c, const Point &p, const Real time)

◆ operator=() [1/2]

template<typename Output = Number>
ParsedFEMFunction& libMesh::ParsedFEMFunction< Output >::operator= ( const ParsedFEMFunction< Output > &  )
delete

This class contains a const reference so it can't be copy or move-assigned.

◆ operator=() [2/2]

template<typename Output = Number>
ParsedFEMFunction& libMesh::ParsedFEMFunction< Output >::operator= ( ParsedFEMFunction< Output > &&  )
delete

◆ partial_reparse()

template<typename Output >
void libMesh::ParsedFEMFunction< Output >::partial_reparse ( const std::string &  expression)
inlineprotected

Definition at line 583 of file parsed_fem_function.h.

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

◆ reparse()

template<typename Output >
void libMesh::ParsedFEMFunction< Output >::reparse ( const std::string &  expression)
inline

Re-parse with new expression.

Definition at line 226 of file parsed_fem_function.h.

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

◆ set_inline_value()

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

Referenced by libMesh::ParsedFEMFunctionParameter< T >::set().

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

Member Data Documentation

◆ _additional_vars

template<typename Output = Number>
std::vector<std::string> libMesh::ParsedFEMFunction< Output >::_additional_vars
private

Definition at line 191 of file parsed_fem_function.h.

◆ _expression

template<typename Output = Number>
std::string libMesh::ParsedFEMFunction< Output >::_expression
private

Definition at line 161 of file parsed_fem_function.h.

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

◆ _initial_vals

template<typename Output = Number>
std::vector<Output> libMesh::ParsedFEMFunction< Output >::_initial_vals
private

Definition at line 192 of file parsed_fem_function.h.

◆ _n_requested_grad_components

template<typename Output = Number>
unsigned int libMesh::ParsedFEMFunction< Output >::_n_requested_grad_components
private

Definition at line 163 of file parsed_fem_function.h.

◆ _n_requested_hess_components

template<typename Output = Number>
unsigned int libMesh::ParsedFEMFunction< Output >::_n_requested_hess_components
private

Definition at line 163 of file parsed_fem_function.h.

◆ _n_requested_vars

template<typename Output = Number>
unsigned int libMesh::ParsedFEMFunction< Output >::_n_requested_vars
private

Definition at line 163 of file parsed_fem_function.h.

◆ _n_vars

template<typename Output = Number>
unsigned int libMesh::ParsedFEMFunction< Output >::_n_vars
private

Definition at line 163 of file parsed_fem_function.h.

◆ _need_var

template<typename Output = Number>
std::vector<bool> libMesh::ParsedFEMFunction< Output >::_need_var
private

Definition at line 178 of file parsed_fem_function.h.

◆ _need_var_grad

template<typename Output = Number>
std::vector<bool> libMesh::ParsedFEMFunction< Output >::_need_var_grad
private

Definition at line 181 of file parsed_fem_function.h.

◆ _need_var_hess

template<typename Output = Number>
std::vector<bool> libMesh::ParsedFEMFunction< Output >::_need_var_hess
private

Definition at line 186 of file parsed_fem_function.h.

◆ _requested_normals

template<typename Output = Number>
bool libMesh::ParsedFEMFunction< Output >::_requested_normals
private

Definition at line 167 of file parsed_fem_function.h.

◆ _spacetime

template<typename Output = Number>
std::vector<Output> libMesh::ParsedFEMFunction< Output >::_spacetime
private

Definition at line 173 of file parsed_fem_function.h.

◆ _subexpressions

template<typename Output = Number>
std::vector<std::string> libMesh::ParsedFEMFunction< Output >::_subexpressions
private

Definition at line 162 of file parsed_fem_function.h.

◆ _sys

template<typename Output = Number>
const System& libMesh::ParsedFEMFunction< Output >::_sys
private

Definition at line 160 of file parsed_fem_function.h.

◆ parsers [1/2]

template<typename Output = Number>
std::vector<FunctionParserBase<Output> > libMesh::ParsedFEMFunction< Output >::parsers
private

Definition at line 169 of file parsed_fem_function.h.

◆ parsers [2/2]

template<typename Output = Number>
std::vector<char> libMesh::ParsedFEMFunction< Output >::parsers
private

Definition at line 171 of file parsed_fem_function.h.

◆ variables

template<typename Output = Number>
std::string libMesh::ParsedFEMFunction< Output >::variables
private

Definition at line 190 of file parsed_fem_function.h.


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