]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Function Parser Wrapper.
authorLuca Heltai <luca.heltai@sissa.it>
Mon, 7 Mar 2005 15:51:37 +0000 (15:51 +0000)
committerLuca Heltai <luca.heltai@sissa.it>
Mon, 7 Mar 2005 15:51:37 +0000 (15:51 +0000)
git-svn-id: https://svn.dealii.org/trunk@10030 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/function_parser.h [new file with mode: 0644]

diff --git a/deal.II/base/include/base/function_parser.h b/deal.II/base/include/base/function_parser.h
new file mode 100644 (file)
index 0000000..170a79f
--- /dev/null
@@ -0,0 +1,350 @@
+//----------------------- function_parser.h  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------- function_parser.h  ---------------------------
+
+
+#ifndef __deal2__function_parser_h
+#define __deal2__function_parser_h
+
+#include <base/config.h>
+#include <base/exceptions.h>
+#include <base/function.h>
+#include <vector>
+#include <map>
+#include <functionparser/fparser.h>
+
+template <int dim> class Point;
+template <int rank_, int dim> class Tensor;
+template <int dim> class Tensor<1,dim>;
+template <typename number> class Vector;
+
+/**
+ * Function Parser. Wrapper class for the fparser library (see
+ * http://www.students.tut.fi/~warp/FunctionParser/). This class lets
+ * you evaluate strings such as "sqrt(1-x^2+y^2)" with given values of
+ * 'x' and 'y'. Some of the informations contained here are copied
+ * verbatim from the fparser.txt file that comes with the fparser
+ * library. Please refer also to that file both for clarifications on
+ * how this wrapper works, as well as for any issue regarding the
+ * licence that applies to this class.
+ *
+ * By using this class you indicate that you accept the terms of the
+ * licence that comes with the fparser library.
+ *
+ * The following example shows how to use this class:
+ @verbatim
+
+ // Define some constants that will be used by the function parser
+ std::map<std::string> constants;
+ constants["pi"] = M_PI;
+ // Define the variables that will be used inside the expressions
+ std::string variables = "x,y,z";
+ // Define the expressions of the vector_valued function.
+ std::vector<std::string> expressions(2);
+ expressions[0] = "sin(2*pi*x)+sinh(pi*z)";
+ expressions[1] = "sin(2*pi*y)*exp(x^2)";
+ // Generate an empty function with two components
+ ParsedFunction<3> vector_function(2);
+ // And populate it with the newly created objects.
+ vector_function.initialize(variables,
+                           expressions,
+                           constants);
+ @endverbatim
+ * 
+ * See http://www.students.tut.fi/~warp/FunctionParser/ for an
+ * explanation on how the underlining library works. 
+ *
+ From the fparser.txt file:
+ @verbatim
+
+ The library is intended to be very fast. It byte-compiles the
+ function string at parse time and interpretes this byte-code at
+ evaluation time.  The evaluation is straightforward and no recursions
+ are done (uses stack arithmetic).  Empirical tests show that it indeed
+ is very fast (specially compared to libraries which evaluate functions
+ by just interpreting the raw function string).
+
+ The library is made in ISO C++ and requires a standard-conforming C++
+ compiler.
+
+ @endverbatim
+ * This object overloads for you the virtual methods value() and
+ * vector_value() of the Function base class with the byte compiled
+ * versions of the expressions given to the initialize() methods.
+ *
+ * @author Luca Heltai, 2005
+ */
+template <int dim>
+class FunctionParser : public Function<dim>
+{
+  public:
+  /** Constructor for Parsed functions. Its arguments are the same of
+      the base class Function. The only difference is that this object
+      needs to be initialized with initialize() method before you can
+      use it. If an attempt to use this function is made before the
+      initialize() method has been called, then an exception is
+      thrown. */
+  FunctionParser (const unsigned int n_components = 1,
+                 const double       initial_time = 0.0); 
+  
+  /** Destructor. Explicitly delete the FunctionParser objects (there
+      is one for each component of the function). */
+  ~FunctionParser();
+  
+  /** Type for the constant map. Used by the initialize() method. */
+  typedef std::map<std::string, double> ConstMap;
+  
+  /** Iterator for the constants map. Used by the initialize() method. */
+  typedef ConstMap::iterator ConstMapIterator;
+
+  /** Initialize the function. From the fparser.txt file:
+      @verbatim
+      The function string understood by the class is very similar to the C-syntax.
+      Arithmetic float expressions can be created from float literals, variables
+      or functions using the following operators in this order of precedence:
+
+      ()             expressions in parentheses first
+      -A             unary minus
+      A^B            exponentiation (A raised to the power B)
+      A*B  A/B  A%B  multiplication, division and modulo
+      A+B  A-B       addition and subtraction
+      A=B  A<B  A>B  comparison between A and B (result is either 0 or 1)
+      A&B            result is 1 if int(A) and int(B) differ from 0, else 0.
+      A|B            result is 1 if int(A) or int(B) differ from 0, else 0.
+
+      Since the unary minus has higher precedence than any other operator, for
+      example the following expression is valid: x*-y
+      Note that the '=' comparison can be inaccurate due to floating point
+      precision problems (eg. "sqrt(100)=10" probably returns 0, not 1).
+      
+      The class supports these functions:
+      
+      abs(A)    : Absolute value of A. If A is negative, returns -A otherwise
+                 returns A.
+      acos(A)   : Arc-cosine of A. Returns the angle, measured in radians,
+                 whose cosine is A.
+      acosh(A)  : Same as acos() but for hyperbolic cosine.
+      asin(A)   : Arc-sine of A. Returns the angle, measured in radians, whose
+                 sine is A.
+      asinh(A)  : Same as asin() but for hyperbolic sine.
+      atan(A)   : Arc-tangent of (A). Returns the angle, measured in radians,
+                  whose tangent is (A).
+      atan2(A,B): Arc-tangent of A/B. The two main differences to atan() is
+                  that it will return the right angle depending on the signs of
+                  A and B (atan() can only return values betwen -pi/2 and pi/2),
+                  and that the return value of pi/2 and -pi/2 are possible.
+      atanh(A)  : Same as atan() but for hyperbolic tangent.
+      ceil(A)   : Ceiling of A. Returns the smallest integer greater than A.
+                  Rounds up to the next higher integer.
+      cos(A)    : Cosine of A. Returns the cosine of the angle A, where A is
+                  measured in radians.
+      cosh(A)   : Same as cos() but for hyperbolic cosine.
+      cot(A)    : Cotangent of A (equivalent to 1/tan(A)).
+      csc(A)    : Cosecant of A (equivalent to 1/sin(A)).
+      eval(...) : This a recursive call to the function to be evaluated. The
+                  number of parameters must be the same as the number of parameters
+                 taken by the function. Usually called inside if() to avoid
+                 infinite recursion.
+      exp(A)    : Exponential of A. Returns the value of e raised to the power
+                  A where e is the base of the natural logarithm, i.e. the
+                  non-repeating value approximately equal to 2.71828182846.
+      floor(A)  : Floor of A. Returns the largest integer less than A. Rounds
+                  down to the next lower integer.
+      if(A,B,C) : If int(A) differs from 0, the return value of this function is B,
+                  else C. Only the parameter which needs to be evaluated is
+                  evaluated, the other parameter is skipped; this makes it safe to
+                  use eval() in them.
+      int(A)    : Rounds A to the closest integer. 0.5 is rounded to 1.
+      log(A)    : Natural (base e) logarithm of A.
+      log10(A)  : Base 10 logarithm of A.
+      max(A,B)  : If A>B, the result is A, else B.
+      min(A,B)  : If A<B, the result is A, else B.
+      sec(A)    : Secant of A (equivalent to 1/cos(A)).
+      sin(A)    : Sine of A. Returns the sine of the angle A, where A is
+                  measured in radians.
+      sinh(A)   : Same as sin() but for hyperbolic sine.
+      sqrt(A)   : Square root of A. Returns the value whose square is A.
+      tan(A)    : Tangent of A. Returns the tangent of the angle A, where A
+                  is measured in radians.
+      tanh(A)   : Same as tan() but for hyperbolic tangent.
+
+
+    Examples of function string understood by the class:
+
+    "1+2"
+    "x-1"
+    "-sin(sqrt(x^2+y^2))"
+    "sqrt(XCoord*XCoord + YCoord*YCoord)"
+    
+    An example of a recursive function is the factorial function:
+    
+    "if(n>1, n*eval(n-1), 1)"
+
+  Note that a recursive call has some overhead, which makes it a bit slower
+  than any other operation. It may be a good idea to avoid recursive functions
+  in very time-critical applications. Recursion also takes some memory, so
+  extremely deep recursions should be avoided (eg. millions of nested recursive
+  calls).
+
+  Also note that the if() function is the only place where making a recursive
+  call is safe. In any other place it will cause an infinite recursion (which
+  will make the program eventually run out of memory).
+  @endverbatim
+  
+  This methods accepts the following parameters: 
+  
+  - <b>vars</b>. It contains a string with the variables that will be
+    used by the expressions to be evaluated. Note that the variables
+    can have any name (of course different from the function names
+    defined above!), but the order IS important. The first variable
+    will correspond to the first component of the point in which the
+    function is evaluated, the second variable to the second component
+    and so forth. If this function is also time dependent, then it is
+    necessary to specify it by setting the #time_dependent parameter
+    to true.     
+    An Exception is thrown if the number of variables specified here
+    is different from dim (if this function is not time-dependent) or
+    from dim+1 (if it is time-dependent).
+    
+  - <b>expressions</b>. This is a list of strings containing the
+    expressions that will be byte compiled by the internal parser
+    (FunctionParser). Note that the size of this vector must match
+    exactly the number of components of the FunctionParser. If this is
+    not the case, an exception is thrown.
+    
+  - <b>constants</b>. The map of constants is used to pass any
+    necessary constant that we want to specify in our expressions (in
+    the example above the number pi). An expression is valid if and
+    only if it contains only defined variables and defined constants
+    (other than the functions specified above). If a constant is given
+    whose name is not valid (eg: <tt>constants["sin"] = 1.5;</tt>) an
+    exception is thrown.
+    
+  - <b>time_dependent</b>. If this is a time dependent function, then
+    the last variable is assumed to be the time variable, and
+    this->get_time() is passed to the function. Naturally the number
+    of variables parsed by the initialize() method in this case is
+    dim+1. It defaults to false, i.e.. do not consider time.
+
+  - <b>use_degrees</b>. Parameter to decide if the trigonometric
+    functions work in radians or degrees. The default for this
+    parameter is false, i.e. use radians and not degrees.
+    
+    
+  */
+  void initialize(const std::string vars,
+                 const std::vector<std::string> expressions,
+                 const ConstMap constants = 0,
+                 bool time_dependent = false,
+                 bool use_degrees = false);
+  
+  /** Initialize the function. Same as above, but for scalar
+      functions. An exception is thrown if this method is called and
+      the function has more than one component. 
+      
+      An example of time dependent scalar function is the following:
+      @verbatim
+     
+      // Empty constants object
+      std::map<std::string> constants;
+      
+      // Variables that will be used inside the expressions
+      std::string variables = "x,y,t";
+      
+      // Define the expression of the scalar time dependent function.
+      std::string expression = "exp(y*x)*exp(-t)";
+      
+      // Generate an empty scalar function
+      FunctionParser<2> function;
+      
+      // And populate it with the newly created objects.
+      function.initialize(variables,
+                         expression,
+                         constants,
+                         true);        // This tells the parser that 
+                                       // it is a time-dependent function
+                                       // and there is another variable
+                                       // to be taken into account (t).
+     @endverbatim
+      
+  */
+  void initialize(const std::string vars,
+                 const std::string expression,
+                 const ConstMap constants = 0,
+                 bool time_dependent = false,
+                 bool use_degrees = false);
+  
+
+  /**
+   * Return the value of the function at the given point. Unless there
+   * is only one component (i.e. the function is scalar), you should
+   * state the component you want to have evaluated; it defaults to
+   * zero, i.e. the first component. */
+  virtual double value (const Point<dim>   &p,
+                       const unsigned int  component = 0) const;
+
+  /**
+   * Return all components of a vector-valued function at a given
+   * point.
+   *
+   * <tt>values</tt> shall have the
+   * right size beforehand,
+   * i.e. #n_components.
+   */
+  virtual void vector_value (const Point<dim>   &p,
+                            Vector<double>     &values) const;
+
+  /** @addtogroup Exceptions
+   * @{ */    
+  DeclException2 (ExcParseError,
+                 int, char*, 
+                 << "Parsing Error at Column " << arg1
+                 << ". The parser said: " << arg2);
+    
+  DeclException2 (ExcInvalidExpressionSize,
+                 int, int, 
+                 << "The number of components (" << arg1
+                 << ") is not equal to the number of expressions (" 
+                 << arg2 << ").");
+    
+  //@}
+ private: 
+  /** A pointer to the actual function parsers. */
+ fparser::FunctionParser * fp;
+    
+ /** State of usability. This variable is checked every time the
+     function is called for evaluation. It's set to true in the
+     initialize() methods.*/
+ bool initialized;
+    
+ /** Number of variables. If this is also a function of time, then the
+     number of variables is dim+1, otherwhise it is dim. In the case
+     that this is a time dependent function, the time is supposed to
+     be the last variable. If #n_vars is not identical to the number
+     of the variables parsed by the initialize() method, then an
+     exception is thrown.   
+ */
+ unsigned int n_vars;
+};
+
+
+#endif
+
+  

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.