]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added functionality to the initialize method. Fixed a bug that prevented non
authorLuca Heltai <luca.heltai@sissa.it>
Thu, 31 Mar 2005 09:23:15 +0000 (09:23 +0000)
committerLuca Heltai <luca.heltai@sissa.it>
Thu, 31 Mar 2005 09:23:15 +0000 (09:23 +0000)
vector valued functions to work properly. Corrected destructor.

git-svn-id: https://svn.dealii.org/trunk@10328 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/function_parser.h
deal.II/base/source/function_parser.cc

index 58de5a2e6dd8c6558900c35a2934a20556a3e0c1..ebb308b8e09c4917d6ea603d6e59c41a1c76013a 100644 (file)
@@ -254,9 +254,12 @@ class FunctionParser : public Function<dim>
                  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. 
+  /** Initialize the function. Same as above, but accepts a string
+      rather than a vector of strings. If this is a vector valued
+      function, its componenents are expected to be separated by a
+      semicolon. An exception is thrown if this method is called and
+      the number of components successfully parsed does not match the
+      number of components of the base function.
       
       An example of time dependent scalar function is the following:
       @verbatim
@@ -283,6 +286,29 @@ class FunctionParser : public Function<dim>
                                        // to be taken into account (t).
  
      @endverbatim
+     
+     The following is yet another example of how to instantiate a
+     vector valued function by using a single string:
+      @verbatim
+     
+      // Empty constants object
+      std::map<std::string> constants;
+      
+      // Variables that will be used inside the expressions
+      std::string variables = "x,y";
+      
+      // Define the expression of the vector valued  function.
+      std::string expression = "cos(2*pi*x)*y^2; sin(2*pi*x)*exp(y)";
+      
+      // Generate an empty vector valued function
+      FunctionParser<2> function(2);
+      
+      // And populate it with the newly created objects.
+      function.initialize(variables,
+                         expression,
+                         constants);   
+     @endverbatim
       
   */
   void initialize(const std::string vars,
index a695395d47f6dc1d1b9b68ee7a419ffd31610440..53bbefc95109d12def6fca56a03652ac9f33c67c 100644 (file)
 #include <base/point.h>
 #include <lac/vector.h>
 
+// Utility to split a string given a delimiter.
+unsigned int SplitString(const std::string& input, 
+                        const std::string& delimiter, 
+                        std::vector<std::string>& results)
+{
+  int iPos = 0;
+  int newPos = -1;
+  int sizeS2 = delimiter.size();
+  int isize = input.size();
+
+  std::vector<int> positions;
+
+  newPos = input.find (delimiter, 0);
+
+  if( newPos < 0 ) { return 0; }
+
+  unsigned int numFound = 0;
+
+  while( newPos > iPos )
+  {
+    numFound++;
+    positions.push_back(newPos);
+    iPos = newPos;
+    newPos = input.find (delimiter, iPos+sizeS2+1);
+  }
+
+  for(unsigned int i=0; i <= positions.size(); ++i)
+  {
+    std::string s;
+    if( i == 0 ) { s = input.substr( i, positions[i] ); }
+    int offset = positions[i-1] + sizeS2;
+    if( offset < isize )
+    {
+      if( i == positions.size() )
+      {
+        s = input.substr(offset);
+      }
+      else if( i > 0 )
+      {
+        s = input.substr( positions[i-1] + sizeS2, 
+          positions[i] - positions[i-1] - sizeS2 );
+      }
+    }
+    if( s.size() > 0 )
+    {
+      results.push_back(s);
+    }
+  }
+  return numFound;
+}
+
+
 template <int dim>
 FunctionParser<dim>::FunctionParser(const unsigned int n_components,
                                    const double       initial_time)
@@ -25,11 +77,9 @@ FunctionParser<dim>::FunctionParser(const unsigned int n_components,
 
 template <int dim>
 FunctionParser<dim>::~FunctionParser() {
-  delete fp;
+  delete[] fp;
 }
 
-
-
 template <int dim>
 void FunctionParser<dim>::initialize(const std::string variables,
                                     const std::vector<std::string> expressions,
@@ -37,21 +87,12 @@ void FunctionParser<dim>::initialize(const std::string variables,
                                     bool time_dependent,
                                     bool use_degrees)
 {
-  // We distinguish between two cases: Time dependent problems, and
-  // not time dependent problems. In the first case the number of
-  // variables is given by the number of components plus one. In the
-  // other case, the number of variables is equal to the number of
-  // components. If none of this is the case, then an exception is
-  // thrown.
-  if (time_dependent) 
-    n_vars = this->n_components+1;
-  else 
-    n_vars = this->n_components;
-
+  // We check that the number of components of this function matches
+  // the number of components passed in as a vector of strings.
   AssertThrow(this->n_components == expressions.size(),
              ExcInvalidExpressionSize(this->n_components,
                                       expressions.size()) );
-
+  
   for (unsigned int i=0; i<this->n_components; ++i) {
     // Add the variuous constants to the parser.
     std::map< std::string, double >::const_iterator
@@ -59,22 +100,33 @@ void FunctionParser<dim>::initialize(const std::string variables,
       endc  = constants.end();
     for(;constant != endc; ++constant) {
       AssertThrow( fp[i].AddConstant(constant->first, constant->second),
-             ExcMessage("Invalid Constant Name"));
+                  ExcMessage("Invalid Constant Name"));
     }
-
-
+      
+      
     int ret_value = fp[i].Parse(expressions[i], 
                                variables, 
                                use_degrees);
     AssertThrow(ret_value == -1, 
                ExcParseError(ret_value, fp[i].ErrorMsg()));
-    
+      
     // The fact that the parser did not throw an error does not mean
     // that everything went ok... we can still have problems with the
     // number of variables...
-
   } 
-  
+    
+  // Now we define how many variables we expect to read in.  We
+  // distinguish between two cases: Time dependent problems, and not
+  // time dependent problems. In the first case the number of
+  // variables is given by the dimension plus one. In the other case,
+  // the number of variables is equal to the dimension. Once we parsed
+  // the variables string, if none of this is the case, then an
+  // exception is thrown.
+  if (time_dependent) 
+    n_vars = dim+1;
+  else 
+    n_vars = dim;
+    
   // Let's check if the number of variables is correct...
   AssertThrow(n_vars == fp[0].NVars(), 
              ExcDimensionMismatch(n_vars,fp[0].NVars()));
@@ -82,6 +134,8 @@ void FunctionParser<dim>::initialize(const std::string variables,
   // Now set the initialization bit.
   initialized = true;
 }
+
+
 template <int dim>
 void FunctionParser<dim>::initialize(const std::string variables,
                                     const std::string expression,
@@ -89,8 +143,16 @@ void FunctionParser<dim>::initialize(const std::string variables,
                                     bool time_dependent,
                                     bool use_degrees)
 {
-  std::vector<std::string> expressions(1,expression);
-  initialize(variables, expressions, constants, time_dependent, use_degrees); 
+  // Start up with an empty vector of expressions.
+  std::vector<std::string> expressions;
+  
+  // In this case a vector function is built with the number of
+  // components found between the separators ';' 
+  SplitString(expression, ";", expressions);
+  if (expressions.size() == 0) expressions.push_back(expression);
+  // Now initialize with the things we got.
+  initialize
+    (variables, expressions, constants, time_dependent, use_degrees);  
 }
 
 template <int dim>
@@ -117,7 +179,6 @@ double FunctionParser<dim>::value (const Point<dim> &p,
   
   AssertThrow(fp[component].EvalError() == 0,
              ExcMessage(fp[component].ErrorMsg()));
-  
   return my_value;
 }
 
@@ -125,29 +186,29 @@ double FunctionParser<dim>::value (const Point<dim> &p,
 template <int dim>
 inline
 void FunctionParser<dim>::vector_value (const Point<dim> &p,
-                                          Vector<double>   &values) const 
+                                       Vector<double>   &values) const 
 {
-  AssertThrow(initialized==true, ExcNotInitialized());
-  AssertThrow(values.size() == this->n_components,
-             ExcDimensionMismatch (values.size(), this->n_components));
-  // Statically allocates dim+1 double variables.
-  double vars[dim+1];
-  
-  for(unsigned int i=0; i<dim; ++i)
-    vars[i] = p(i);
-  
-  // We need the time variable only if the number of variables is
-  // different from the dimension. In this case it can only be dim+1,
-  // otherwise an exception would have already been thrown
-  if(dim != n_vars)
-    vars[dim] = this->get_time();
-  
-  for(unsigned int component = 0; component < this->n_components;
-      ++component) {
-    values(component) = fp[component].Eval(vars);
-    Assert(fp[component].EvalError() == 0,
-          ExcMessage(fp[component].ErrorMsg()));
-  }
+    AssertThrow(initialized==true, ExcNotInitialized());
+    AssertThrow(values.size() == this->n_components,
+               ExcDimensionMismatch (values.size(), this->n_components));
+    // Statically allocates dim+1 double variables.
+    double vars[dim+1];
+    
+    for(unsigned int i=0; i<dim; ++i)
+      vars[i] = p(i);
+    
+    // We need the time variable only if the number of variables is
+    // different from the dimension. In this case it can only be dim+1,
+    // otherwise an exception would have already been thrown
+    if(dim != n_vars)
+      vars[dim] = this->get_time();
+    
+    for(unsigned int component = 0; component < this->n_components;
+       ++component) {
+      values(component) = fp[component].Eval(vars);
+      Assert(fp[component].EvalError() == 0,
+            ExcMessage(fp[component].ErrorMsg()));
+    } 
 }
 
 // Explicit Instantiations.

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.