use_degrees);
}
+
template <int dim>
void FunctionParser<dim>::initialize (const std::string &vars,
const std::vector<std::string> &expressions,
const std::map<std::string, double> &constants,
const bool time_dependent)
- {
- initialize(vars, expressions, constants, time_dependent, false);
- }
+{
+ initialize(vars, expressions, constants, time_dependent, false);
+}
+
+
namespace internal
{
{
return static_cast<double>(mu_round(value));
}
+
double mu_ceil(double value)
{
return ceil(value);
}
+
double mu_floor(double value)
{
return floor(value);
}
+
double mu_cot(double value)
{
return 1.0/tan(value);
}
+
double mu_csc(double value)
{
return 1.0/sin(value);
}
- double mu_sec(double value)
+
+ double mu_sec(double value)
{
return 1.0/cos(value);
}
- double mu_log(double value)
+
+ double mu_log(double value)
{
return log(value);
}
-
-
-
}
+
template <int dim>
void FunctionParser<dim>:: init_muparser() const
{
vars.get().resize(var_names.size());
for (unsigned int i=0; i<this->n_components; ++i)
{
- std::map< std::string, double >::const_iterator
- constant = constants.begin(),
- endc = constants.end();
-
- for (; constant != endc; ++constant)
+ for (std::map< std::string, double >::const_iterator constant = constants.begin();
+ constant != constants.end(); ++constant)
{
fp.get()[i].DefineConst(constant->first.c_str(), constant->second);
}
try
{
- fp.get()[i].SetExpr(expressions[i]);
+ // muparser expects that user defined functions have no
+ // space between the name of the function and the opening
+ // parenthesis. this is awkward because it is not backward
+ // compatible to the library we used to use before muparser
+ // (the fparser library) but also makes no real sense.
+ // consequently, in the expressions we set, remove any space
+ // we may find after function names
+ std::string transformed_expression = expressions[i];
+ std::string::size_type pos = 0;
+ while (true)
+ {
+ // try to find any occurrences of 'if'
+ pos = transformed_expression.find ("if", pos);
+ if (pos == std::string::npos)
+ break;
+
+ // replace whitespace until there no longer is any
+ while ((pos+2<transformed_expression.size())
+ &&
+ ((transformed_expression[pos+2] == ' ')
+ ||
+ (transformed_expression[pos+2] == '\t')))
+ transformed_expression.erase (pos+2, pos+3);
+
+ // move the current search position by the size of the
+ // actual 'if'
+ pos += 2;
+ }
+
+ // now use the transformed expression
+ fp.get()[i].SetExpr(transformed_expression);
}
catch (mu::ParserError &e)
{
}
}
+
template <int dim>
void FunctionParser<dim>::initialize (const std::string &variables,
const std::vector<std::string> &expressions,
const bool time_dependent,
const bool use_degrees)
{
-
this->fp.clear(); // this will reset all thread-local objects
this->constants = constants;
initialized = true;
}
+
+
template <int dim>
void FunctionParser<dim>::initialize (const std::string &vars,
- const std::string &expression,
- const std::map<std::string, double> &constants,
- const bool time_dependent)
+ const std::string &expression,
+ const std::map<std::string, double> &constants,
+ const bool time_dependent)
{
initialize(vars, expression, constants, time_dependent, false);
}
+
template <int dim>
void FunctionParser<dim>::initialize (const std::string &variables,
const std::string &expression,
+// ---------------------------------------------------------------------
+// $Id$
+//
+// Copyright (C) 2005 - 2014 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// various checks (compatibility functionparser -> muparser)
+
+
#include "../tests.h"
#include <fstream>
#include <iomanip>
#include <deal.II/lac/vector.h>
#include <deal.II/base/function_parser.h>
-// various checks (compatibility functionparser -> muparser)
void eval(const std::string & exp, const Point<2> & p, double expected)
{
--- /dev/null
+// ---------------------------------------------------------------------
+// $Id$
+//
+// Copyright (C) 2005 - 2014 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// like _06, but ensure that we can not only parse "if(cond,yes,no)" but also
+// "if (cond,yes,no)", i.e., including the whitespace between function and
+// argument list
+
+#include "../tests.h"
+#include <fstream>
+#include <iomanip>
+#include <map>
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/point.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/base/function_parser.h>
+
+
+void eval(const std::string & exp, const Point<2> & p, double expected)
+{
+ std::string variables = "x,y";
+ std::map<std::string,double> constants;
+
+ FunctionParser<2> fp(1);
+ fp.initialize(variables,
+ exp,
+ constants);
+
+ double result = fp.value(p);
+ deallog << "'" << exp << "' @ " << p << " is " << result
+ << " ( expected " << expected << " )" << std::endl;
+ if (fabs(result-expected)>1e-10)
+ deallog << "ERROR!" << std::endl;
+}
+
+
+void test()
+{
+ eval("if (x<0.0,0.0,if(x>1.0,1.0,x))",Point<2>(0.5,0.0), 0.5);
+ eval("if (x<0.0,0.0,if(x>1.0,1.0,x))",Point<2>(-2.0,0.0), 0.0);
+ eval("if (x<0.0,0.0,if(x>1.0,1.0,x))",Point<2>(42.0,0.0), 1.0);
+
+ eval("if (x<1.0 | y<1.0,0,y)",Point<2>(0.5,1.5), 0.0);
+ eval("if (x<1.0 | y<1.0,0,y)",Point<2>(1.5,1.5), 1.5);
+ eval("if (x<1.0 | y<1.0,0,y)",Point<2>(1.5,0.5), 0.0);
+ eval("if (x<1.0 | y<1.0,0,y)",Point<2>(0.5,-2), 0.0);
+
+ eval("if (x<1.0 & y<1.0,0,y)",Point<2>(1.5,-2.0), -2.0);
+
+ double x,y;
+ x=1.0;
+ y=-3.1;
+ eval("atan2(x,y)",Point<2>(x,y), atan2(x,y));
+ x=-1.0;
+ y=3.1;
+ eval("atan2(x,y)",Point<2>(x,y), atan2(x,y));
+
+ eval("if (x==1.0,0,y)",Point<2>(1.0,-2.0), 0.0);
+ eval("if (x==1.0,0,y)",Point<2>(1.1,-2.0), -2.0);
+
+ eval("int(2.1)",Point<2>(1.1,-2.0), 2.0);
+ eval("int(-3.8)",Point<2>(1.1,-2.0), -4.0);
+
+ eval("abs(-2.3)",Point<2>(0,0), 2.3);
+ eval("acos(0.5)",Point<2>(0,0), acos(0.5));
+ eval("acosh(0.5)",Point<2>(0,0), acosh(0.5));
+ eval("asin(0.5)",Point<2>(0,0), asin(0.5));
+ eval("asinh(0.5)",Point<2>(0,0), asinh(0.5));
+ eval("atan(0.5)",Point<2>(0,0), atan(0.5));
+ eval("atanh(0.5)",Point<2>(0,0), atanh(0.5));
+ eval("ceil(0.5)",Point<2>(0,0), ceil(0.5));
+ eval("cos(0.5)",Point<2>(0,0), cos(0.5));
+ eval("cosh(0.5)",Point<2>(0,0), cosh(0.5));
+ eval("cot(0.5)",Point<2>(0,0), 1.0/tan(0.5));
+ eval("csc(0.5)",Point<2>(0,0), 1.0/sin(0.5));
+ eval("exp(0.5)",Point<2>(0,0), exp(0.5));
+ eval("floor(0.5)",Point<2>(0,0), floor(0.5));
+ eval("log(0.5)",Point<2>(0,0), log(0.5));
+ eval("log10(0.5)",Point<2>(0,0), log(0.5)/log(10.0));
+ eval("sec(0.5)",Point<2>(0,0), 1.0/cos(0.5));
+ eval("sin(0.5)",Point<2>(0,0), sin(0.5));
+ eval("sinh(0.5)",Point<2>(0,0), sinh(0.5));
+ eval("sqrt(0.5)",Point<2>(0,0), sqrt(0.5));
+ eval("tan(0.5)",Point<2>(0,0), tan(0.5));
+ eval("tanh(0.5)",Point<2>(0,0), tanh(0.5));
+
+
+}
+
+int main()
+{
+ initlog();
+
+ test();
+}
--- /dev/null
+
+DEAL::'if (x<0.0,0.0,if(x>1.0,1.0,x))' @ 0.500000 0.00000 is 0.500000 ( expected 0.500000 )
+DEAL::'if (x<0.0,0.0,if(x>1.0,1.0,x))' @ -2.00000 0.00000 is 0.00000 ( expected 0.00000 )
+DEAL::'if (x<0.0,0.0,if(x>1.0,1.0,x))' @ 42.0000 0.00000 is 1.00000 ( expected 1.00000 )
+DEAL::'if (x<1.0 | y<1.0,0,y)' @ 0.500000 1.50000 is 0.00000 ( expected 0.00000 )
+DEAL::'if (x<1.0 | y<1.0,0,y)' @ 1.50000 1.50000 is 1.50000 ( expected 1.50000 )
+DEAL::'if (x<1.0 | y<1.0,0,y)' @ 1.50000 0.500000 is 0.00000 ( expected 0.00000 )
+DEAL::'if (x<1.0 | y<1.0,0,y)' @ 0.500000 -2.00000 is 0.00000 ( expected 0.00000 )
+DEAL::'if (x<1.0 & y<1.0,0,y)' @ 1.50000 -2.00000 is -2.00000 ( expected -2.00000 )
+DEAL::'atan2(x,y)' @ 1.00000 -3.10000 is 2.82955 ( expected 2.82955 )
+DEAL::'atan2(x,y)' @ -1.00000 3.10000 is -0.312042 ( expected -0.312042 )
+DEAL::'if (x==1.0,0,y)' @ 1.00000 -2.00000 is 0.00000 ( expected 0.00000 )
+DEAL::'if (x==1.0,0,y)' @ 1.10000 -2.00000 is -2.00000 ( expected -2.00000 )
+DEAL::'int(2.1)' @ 1.10000 -2.00000 is 2.00000 ( expected 2.00000 )
+DEAL::'int(-3.8)' @ 1.10000 -2.00000 is -4.00000 ( expected -4.00000 )
+DEAL::'abs(-2.3)' @ 0.00000 0.00000 is 2.30000 ( expected 2.30000 )
+DEAL::'acos(0.5)' @ 0.00000 0.00000 is 1.04720 ( expected 1.04720 )
+DEAL::'acosh(0.5)' @ 0.00000 0.00000 is -nan ( expected -nan )
+DEAL::'asin(0.5)' @ 0.00000 0.00000 is 0.523599 ( expected 0.523599 )
+DEAL::'asinh(0.5)' @ 0.00000 0.00000 is 0.481212 ( expected 0.481212 )
+DEAL::'atan(0.5)' @ 0.00000 0.00000 is 0.463648 ( expected 0.463648 )
+DEAL::'atanh(0.5)' @ 0.00000 0.00000 is 0.549306 ( expected 0.549306 )
+DEAL::'ceil(0.5)' @ 0.00000 0.00000 is 1.00000 ( expected 1.00000 )
+DEAL::'cos(0.5)' @ 0.00000 0.00000 is 0.877583 ( expected 0.877583 )
+DEAL::'cosh(0.5)' @ 0.00000 0.00000 is 1.12763 ( expected 1.12763 )
+DEAL::'cot(0.5)' @ 0.00000 0.00000 is 1.83049 ( expected 1.83049 )
+DEAL::'csc(0.5)' @ 0.00000 0.00000 is 2.08583 ( expected 2.08583 )
+DEAL::'exp(0.5)' @ 0.00000 0.00000 is 1.64872 ( expected 1.64872 )
+DEAL::'floor(0.5)' @ 0.00000 0.00000 is 0.00000 ( expected 0.00000 )
+DEAL::'log(0.5)' @ 0.00000 0.00000 is -0.693147 ( expected -0.693147 )
+DEAL::'log10(0.5)' @ 0.00000 0.00000 is -0.301030 ( expected -0.301030 )
+DEAL::'sec(0.5)' @ 0.00000 0.00000 is 1.13949 ( expected 1.13949 )
+DEAL::'sin(0.5)' @ 0.00000 0.00000 is 0.479426 ( expected 0.479426 )
+DEAL::'sinh(0.5)' @ 0.00000 0.00000 is 0.521095 ( expected 0.521095 )
+DEAL::'sqrt(0.5)' @ 0.00000 0.00000 is 0.707107 ( expected 0.707107 )
+DEAL::'tan(0.5)' @ 0.00000 0.00000 is 0.546302 ( expected 0.546302 )
+DEAL::'tanh(0.5)' @ 0.00000 0.00000 is 0.462117 ( expected 0.462117 )
projects / dealii.git / commitdiff