and has an interface similar to FunctionParser. Some stuff to iterface muParser was
moved to different files mu_parser_internal.h and *.cc. The documentation provided
is self contained apart from some places where we refer to the documentation of
FuntionParser. Changelog file provided. (All commits are squashed into this one.)
--- /dev/null
+New: The TensorFunctionParser class allows to read
+in a function similarly to the FunctionParser class using the MuParser.
+<br>
+(Konrad Simon, 2019/12/17)
\ No newline at end of file
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_mu_parser_internal_h
+#define dealii_mu_parser_internal_h
+
+// This file contains functions used internally by the FunctionParser
+// and the TensorFunctionParser class.
+
+#include <deal.II/base/config.h>
+
+#include <string>
+#include <vector>
+
+
+DEAL_II_NAMESPACE_OPEN
+
+
+
+#ifdef DEAL_II_WITH_MUPARSER
+
+namespace internal
+{
+ namespace FunctionParser
+ {
+ int
+ mu_round(double val);
+
+ double
+ mu_if(double condition, double thenvalue, double elsevalue);
+
+ double
+ mu_or(double left, double right);
+
+ double
+ mu_and(double left, double right);
+
+ double
+ mu_int(double value);
+
+ double
+ mu_ceil(double value);
+
+ double
+ mu_floor(double value);
+
+ double
+ mu_cot(double value);
+
+ double
+ mu_csc(double value);
+
+ double
+ mu_sec(double value);
+
+ double
+ mu_log(double value);
+
+ double
+ mu_pow(double a, double b);
+
+ double
+ mu_erfc(double value);
+
+ // returns a random value in the range [0,1] initializing the generator
+ // with the given seed
+ double
+ mu_rand_seed(double seed);
+
+ // returns a random value in the range [0,1]
+ double
+ mu_rand();
+
+ extern std::vector<std::string> function_names;
+
+ } // namespace FunctionParser
+
+} // namespace internal
+#endif
+
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_tensor_function_parser_h
+#define dealii_tensor_function_parser_h
+
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/base/exceptions.h>
+#include <deal.II/base/point.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_function.h>
+#include <deal.II/base/thread_local_storage.h>
+
+#include <map>
+#include <memory>
+#include <vector>
+
+namespace mu
+{
+ class Parser;
+}
+
+DEAL_II_NAMESPACE_OPEN
+
+// Forward declaration
+#ifndef DOXYGEN
+template <typename>
+class Vector;
+#endif
+
+/**
+ * This class implements a tensor function object that gets its value by parsing
+ * a string describing this function. It is a wrapper class for the muparser
+ * library (see http://muparser.beltoforion.de/). This class is essentially an
+ * extension of the FunctionParser class to read in a TensorFunction. The class
+ * reads in an expression of length dim<sup>rank</sup> (separated by a
+ * semicolon) where the components of the tensor function are filled according
+ * to the C++ convention (fastest index is the most right one).
+ *
+ * @note In contrast to the FunctionParser class the TensorFunctionParser class does not support
+ * automatic differentiation.
+ *
+ * A minimal example for the usage of the class would be:
+ * @code
+ * // set up time dependent tensor function:
+ * const std::string variables = "x,y,t";
+ * const std::string expression =
+ * "exp(-t)*cos(x+y);-sin(pi*x*y-t);sin(pi*x*y-t);exp(t)*cos(x+y)";
+ * std::map<std::string,double> constants;
+ * constants["pi"] = numbers::PI;
+ *
+ * // TensorFunctionParser with 2+1 variables (space + time) in 2D of rank 2.
+ * // It is necessary to tell the parser that there is an additional variable
+ * // to be taken into account (t).
+ * TensorFunctionParser<2,2> tfp;
+ * tfp.initialize(variables,
+ * expression,
+ * constants,
+ * true); // flag for time dependence
+ *
+ * // Point at which we want to evaluate the function
+ * Point<2> point(0.0, 1.0);
+ *
+ * // evaluate the expression at 'point':
+ * double result = tfp.value(point);
+ *
+ * deallog << "Function '" << expression << "'"
+ * << " @ " << point
+ * << " is: "
+ * << std::endl
+ * << result[0][0] << " " << result[0][1] << std::endl
+ * << result[1][0] << " " << result[1][1]
+ * << std::endl;
+ * @endcode
+ *
+ * See also the documentation of the FunctionParser class.
+ *
+ * This class overloads the virtual method value() and value_list() of the
+ * TensorFunction base class with the byte compiled versions of the expressions
+ * given to the initialize() methods. Note that the class will not work unless
+ * you first call the initialize() method that accepts the text description of
+ * the function as an argument (among other things).
+ *
+ * The syntax to describe a function follows usual programming practice, and
+ * is explained in detail at the homepage of the underlying muparser library
+ * at http://muparser.beltoforion.de/ .
+ *
+ *
+ * Vector-valued functions can either be declared using strings where the
+ * function components are separated by semicolons, or using a vector of
+ * strings each defining one vector component.
+ *
+ *
+ * @ingroup functions
+ * @author Konrad Simon, 2019
+ */
+template <int rank, int dim, typename Number = double>
+class TensorFunctionParser : public TensorFunction<rank, dim, Number>
+{
+public:
+ /**
+ * Standard constructor. Only set initial time. This object needs to be
+ * initialized with the 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.
+ */
+ TensorFunctionParser(const double initial_time = 0.0);
+
+ /**
+ * Constructor for parsed functions. This object needs to be initialized
+ * with the 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.
+ * Takes a semicolon separated list of expressions (one for each component
+ * of the tensor function), an optional comma-separated list of constants.
+ */
+ TensorFunctionParser(
+ const std::string &expression,
+ const std::string &constants = "",
+ const std::string &variable_names = default_variable_names() + ",t");
+
+ /**
+ * Copy constructor. Objects of this type can not be copied, and
+ * consequently this constructor is deleted.
+ */
+ TensorFunctionParser(const TensorFunctionParser &) = delete;
+
+ /**
+ * Move constructor. Objects of this type can not be moved, and
+ * consequently this constructor is deleted.
+ */
+ TensorFunctionParser(TensorFunctionParser &&) = delete;
+
+ /**
+ * Destructor.
+ */
+ virtual ~TensorFunctionParser() override;
+
+ /**
+ * Copy operator. Objects of this type can not be copied, and
+ * consequently this operator is deleted.
+ */
+ TensorFunctionParser &
+ operator=(const TensorFunctionParser &) = delete;
+
+ /**
+ * Move operator. Objects of this type can not be moved, and
+ * consequently this operator is deleted.
+ */
+ TensorFunctionParser &
+ operator=(TensorFunctionParser &&) = delete;
+
+ /**
+ * Type for the constant map. Used by the initialize() method.
+ */
+ using ConstMap = std::map<std::string, double>;
+
+
+ /**
+ * Initialize the tensor function. This method accepts the following
+ * parameters:
+ *
+ * @param[in] vars 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
+ * <code>time_dependent</code> 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).
+ *
+ * @param[in] expressions A vector of strings containing the expressions that
+ * will be byte compiled by the internal parser (TensorFunctionParser). Note
+ * that the size of this vector must match exactly the number of components of
+ * the TensorFunctionParser, as declared in the constructor. If this is not
+ * the case, an exception is thrown.
+ *
+ * @param[in] constants A map of constants 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:
+ * <code>constants["sin"] = 1.5;</code>) an exception is thrown.
+ *
+ * @param[in] time_dependent If this is a time dependent function, then the
+ * last variable declared in <b>vars</b> is assumed to be the time variable,
+ * and this->get_time() is used to initialize it when evaluating the
+ * function. Naturally the number of variables parsed by the initialize()
+ * method in this case is dim+1. The value of this parameter defaults to
+ * false, i.e. do not consider time.
+ */
+ void
+ initialize(const std::string & vars,
+ const std::vector<std::string> &expressions,
+ const ConstMap & constants,
+ const bool time_dependent = false);
+
+ /**
+ * Initialize the function. Same as above, but accepts a string rather than
+ * a vector of strings. If this is a vector valued function, its components
+ * 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.
+ */
+ void
+ initialize(const std::string &vars,
+ const std::string &expression,
+ const ConstMap & constants,
+ const bool time_dependent = false);
+
+ /**
+ * A function that returns default names for variables, to be used in the
+ * first argument of the initialize() functions: it returns "x" in 1d, "x,y"
+ * in 2d, and "x,y,z" in 3d.
+ */
+ static std::string
+ default_variable_names();
+
+ /**
+ * Return the value of the tensor function at the given point.
+ */
+ virtual Tensor<rank, dim, Number>
+ value(const Point<dim> &p) const override;
+
+ /**
+ * Return the value of the tensor function at the given point.
+ */
+ virtual void
+ value_list(const std::vector<Point<dim>> & p,
+ std::vector<Tensor<rank, dim, Number>> &values) const override;
+
+ /**
+ * Return an array of function expressions (one per component), used to
+ * initialize this function.
+ */
+ const std::vector<std::string> &
+ get_expressions() const;
+
+ /**
+ * @addtogroup Exceptions
+ * @{
+ */
+ DeclException2(ExcParseError,
+ int,
+ std::string,
+ << "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:
+#ifdef DEAL_II_WITH_MUPARSER
+ /**
+ * Place for the variables for each thread
+ */
+ mutable Threads::ThreadLocalStorage<std::vector<double>> vars;
+
+ /**
+ * The muParser objects for each thread (and one for each component). We are
+ * storing a unique_ptr so that we don't need to include the definition of
+ * mu::Parser in this header.
+ */
+ mutable Threads::ThreadLocalStorage<std::vector<std::unique_ptr<mu::Parser>>>
+ tfp;
+
+ /**
+ * An array to keep track of all the constants, required to initialize tfp in
+ * each thread.
+ */
+ std::map<std::string, double> constants;
+
+ /**
+ * An array for the variable names, required to initialize tfp in each
+ * thread.
+ */
+ std::vector<std::string> var_names;
+
+ /**
+ * Initialize tfp and vars on the current thread. This function may only be
+ * called once per thread. A thread can test whether the function has
+ * already been called by testing whether 'tfp.get().size()==0' (not
+ * initialized) or >0 (already initialized).
+ */
+ void
+ init_muparser() const;
+#endif
+
+ /**
+ * An array of function expressions (one per component), required to
+ * initialize tfp in each thread.
+ */
+ std::vector<std::string> expressions;
+
+ /**
+ * 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, otherwise 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;
+
+ /**
+ * Number of components is equal dim<sup>rank<\sup>.
+ */
+ unsigned int n_components;
+};
+
+
+template <int rank, int dim, typename Number>
+std::string
+TensorFunctionParser<rank, dim, Number>::default_variable_names()
+{
+ switch (dim)
+ {
+ case 1:
+ return "x";
+ case 2:
+ return "x,y";
+ case 3:
+ return "x,y,z";
+ default:
+ Assert(false, ExcNotImplemented());
+ }
+ return "";
+}
+
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
logstream.cc
hdf5.cc
mpi.cc
+ mu_parser_internal.cc
multithread_info.cc
named_selection.cc
numbers.cc
subscriptor.cc
table_handler.cc
tensor_function.cc
+ tensor_function_parser.cc
tensor_product_polynomials.cc
tensor_product_polynomials_bubbles.cc
tensor_product_polynomials_const.cc
symmetric_tensor.inst.in
tensor.inst.in
tensor_function.inst.in
+ tensor_function_parser.inst.in
time_stepping.inst.in
)
#include <deal.II/base/function_parser.h>
+#include <deal.II/base/mu_parser_internal.h>
#include <deal.II/base/patterns.h>
#include <deal.II/base/thread_management.h>
#include <deal.II/base/utilities.h>
-namespace internal
-{
- // convert double into int
- int
- mu_round(double val)
- {
- return static_cast<int>(val + ((val >= 0.0) ? 0.5 : -0.5));
- }
-
- double
- mu_if(double condition, double thenvalue, double elsevalue)
- {
- if (mu_round(condition))
- return thenvalue;
- else
- return elsevalue;
- }
-
- double
- mu_or(double left, double right)
- {
- return (mu_round(left)) || (mu_round(right));
- }
-
- double
- mu_and(double left, double right)
- {
- return (mu_round(left)) && (mu_round(right));
- }
-
- double
- mu_int(double value)
- {
- return static_cast<double>(mu_round(value));
- }
-
- double
- mu_ceil(double value)
- {
- return std::ceil(value);
- }
-
- double
- mu_floor(double value)
- {
- return std::floor(value);
- }
-
- double
- mu_cot(double value)
- {
- return 1.0 / std::tan(value);
- }
-
- double
- mu_csc(double value)
- {
- return 1.0 / std::sin(value);
- }
-
- double
- mu_sec(double value)
- {
- return 1.0 / std::cos(value);
- }
-
- double
- mu_log(double value)
- {
- return std::log(value);
- }
-
- double
- mu_pow(double a, double b)
- {
- return std::pow(a, b);
- }
-
- double
- mu_erfc(double value)
- {
- return std::erfc(value);
- }
-
- // returns a random value in the range [0,1] initializing the generator
- // with the given seed
- double
- mu_rand_seed(double seed)
- {
- static Threads::Mutex rand_mutex;
- std::lock_guard<std::mutex> lock(rand_mutex);
-
- static boost::random::uniform_real_distribution<> uniform_distribution(0,
- 1);
-
- // for each seed an unique random number generator is created,
- // which is initialized with the seed itself
- // we could use std::mt19937 but doing so results in compiler-dependent
- // output.
- static std::map<double, boost::random::mt19937> rng_map;
-
- if (rng_map.find(seed) == rng_map.end())
- rng_map[seed] = boost::random::mt19937(static_cast<unsigned int>(seed));
-
- return uniform_distribution(rng_map[seed]);
- }
-
- // returns a random value in the range [0,1]
- double
- mu_rand()
- {
- static Threads::Mutex rand_mutex;
- std::lock_guard<std::mutex> lock(rand_mutex);
- static boost::random::uniform_real_distribution<> uniform_distribution(0,
- 1);
- static boost::random::mt19937 rng(
- static_cast<unsigned long>(std::time(nullptr)));
- return uniform_distribution(rng);
- }
-
-} // namespace internal
-
-
-
template <int dim>
void
FunctionParser<dim>::init_muparser() const
fp.get()[component]->DefineVar(var_names[iv], &vars.get()[iv]);
// define some compatibility functions:
- fp.get()[component]->DefineFun("if", internal::mu_if, true);
- fp.get()[component]->DefineOprt("|", internal::mu_or, 1);
- fp.get()[component]->DefineOprt("&", internal::mu_and, 2);
- fp.get()[component]->DefineFun("int", internal::mu_int, true);
- fp.get()[component]->DefineFun("ceil", internal::mu_ceil, true);
- fp.get()[component]->DefineFun("cot", internal::mu_cot, true);
- fp.get()[component]->DefineFun("csc", internal::mu_csc, true);
- fp.get()[component]->DefineFun("floor", internal::mu_floor, true);
- fp.get()[component]->DefineFun("sec", internal::mu_sec, true);
- fp.get()[component]->DefineFun("log", internal::mu_log, true);
- fp.get()[component]->DefineFun("pow", internal::mu_pow, true);
- fp.get()[component]->DefineFun("erfc", internal::mu_erfc, true);
- fp.get()[component]->DefineFun("rand_seed", internal::mu_rand_seed, true);
- fp.get()[component]->DefineFun("rand", internal::mu_rand, true);
+ fp.get()[component]->DefineFun("if",
+ internal::FunctionParser::mu_if,
+ true);
+ fp.get()[component]->DefineOprt("|", internal::FunctionParser::mu_or, 1);
+ fp.get()[component]->DefineOprt("&", internal::FunctionParser::mu_and, 2);
+ fp.get()[component]->DefineFun("int",
+ internal::FunctionParser::mu_int,
+ true);
+ fp.get()[component]->DefineFun("ceil",
+ internal::FunctionParser::mu_ceil,
+ true);
+ fp.get()[component]->DefineFun("cot",
+ internal::FunctionParser::mu_cot,
+ true);
+ fp.get()[component]->DefineFun("csc",
+ internal::FunctionParser::mu_csc,
+ true);
+ fp.get()[component]->DefineFun("floor",
+ internal::FunctionParser::mu_floor,
+ true);
+ fp.get()[component]->DefineFun("sec",
+ internal::FunctionParser::mu_sec,
+ true);
+ fp.get()[component]->DefineFun("log",
+ internal::FunctionParser::mu_log,
+ true);
+ fp.get()[component]->DefineFun("pow",
+ internal::FunctionParser::mu_pow,
+ true);
+ fp.get()[component]->DefineFun("erfc",
+ internal::FunctionParser::mu_erfc,
+ true);
+ fp.get()[component]->DefineFun("rand_seed",
+ internal::FunctionParser::mu_rand_seed,
+ true);
+ fp.get()[component]->DefineFun("rand",
+ internal::FunctionParser::mu_rand,
+ true);
try
{
// we may find after function names
std::string transformed_expression = expressions[component];
- const char *function_names[] = {// functions predefined by muparser
- "sin",
- "cos",
- "tan",
- "asin",
- "acos",
- "atan",
- "sinh",
- "cosh",
- "tanh",
- "asinh",
- "acosh",
- "atanh",
- "atan2",
- "log2",
- "log10",
- "log",
- "ln",
- "exp",
- "sqrt",
- "sign",
- "rint",
- "abs",
- "min",
- "max",
- "sum",
- "avg",
- // functions we define ourselves above
- "if",
- "int",
- "ceil",
- "cot",
- "csc",
- "floor",
- "sec",
- "pow",
- "erfc",
- "rand",
- "rand_seed"};
- for (const auto &function_name_c_string : function_names)
+ for (const auto ¤t_function_name :
+ internal::FunctionParser::function_names)
{
- const std::string function_name = function_name_c_string;
- const unsigned int function_name_length = function_name.size();
+ const unsigned int function_name_length =
+ current_function_name.size();
std::string::size_type pos = 0;
while (true)
{
// try to find any occurrences of the function name
- pos = transformed_expression.find(function_name, pos);
+ pos = transformed_expression.find(current_function_name, pos);
if (pos == std::string::npos)
break;
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/mu_parser_internal.h>
+#include <deal.II/base/thread_management.h>
+
+#include <boost/random.hpp>
+
+#include <cmath>
+#include <ctime>
+#include <map>
+#include <mutex>
+#include <vector>
+
+
+DEAL_II_NAMESPACE_OPEN
+
+
+
+#ifdef DEAL_II_WITH_MUPARSER
+
+namespace internal
+{
+ namespace FunctionParser
+ {
+ int
+ mu_round(double val)
+ {
+ return static_cast<int>(val + ((val >= 0.0) ? 0.5 : -0.5));
+ }
+
+ double
+ mu_if(double condition, double thenvalue, double elsevalue)
+ {
+ if (mu_round(condition))
+ return thenvalue;
+ else
+ return elsevalue;
+ }
+
+ double
+ mu_or(double left, double right)
+ {
+ return (mu_round(left)) || (mu_round(right));
+ }
+
+ double
+ mu_and(double left, double right)
+ {
+ return (mu_round(left)) && (mu_round(right));
+ }
+
+ double
+ mu_int(double value)
+ {
+ return static_cast<double>(mu_round(value));
+ }
+
+ double
+ mu_ceil(double value)
+ {
+ return std::ceil(value);
+ }
+
+ double
+ mu_floor(double value)
+ {
+ return std::floor(value);
+ }
+
+ double
+ mu_cot(double value)
+ {
+ return 1.0 / std::tan(value);
+ }
+
+ double
+ mu_csc(double value)
+ {
+ return 1.0 / std::sin(value);
+ }
+
+ double
+ mu_sec(double value)
+ {
+ return 1.0 / std::cos(value);
+ }
+
+ double
+ mu_log(double value)
+ {
+ return std::log(value);
+ }
+
+ double
+ mu_pow(double a, double b)
+ {
+ return std::pow(a, b);
+ }
+
+ double
+ mu_erfc(double value)
+ {
+ return std::erfc(value);
+ }
+
+ // returns a random value in the range [0,1] initializing the generator
+ // with the given seed
+ double
+ mu_rand_seed(double seed)
+ {
+ static Threads::Mutex rand_mutex;
+ std::lock_guard<std::mutex> lock(rand_mutex);
+
+ static boost::random::uniform_real_distribution<> uniform_distribution(0,
+ 1);
+
+ // for each seed an unique random number generator is created,
+ // which is initialized with the seed itself
+ // we could use std::mt19937 but doing so results in compiler-dependent
+ // output.
+ static std::map<double, boost::random::mt19937> rng_map;
+
+ if (rng_map.find(seed) == rng_map.end())
+ rng_map[seed] = boost::random::mt19937(static_cast<unsigned int>(seed));
+
+ return uniform_distribution(rng_map[seed]);
+ }
+
+ // returns a random value in the range [0,1]
+ double
+ mu_rand()
+ {
+ static Threads::Mutex rand_mutex;
+ std::lock_guard<std::mutex> lock(rand_mutex);
+ static boost::random::uniform_real_distribution<> uniform_distribution(0,
+ 1);
+ static boost::random::mt19937 rng(
+ static_cast<unsigned long>(std::time(nullptr)));
+ return uniform_distribution(rng);
+ }
+
+ std::vector<std::string> function_names = {
+ // functions predefined by muparser
+ "sin",
+ "cos",
+ "tan",
+ "asin",
+ "acos",
+ "atan",
+ "sinh",
+ "cosh",
+ "tanh",
+ "asinh",
+ "acosh",
+ "atanh",
+ "atan2",
+ "log2",
+ "log10",
+ "log",
+ "ln",
+ "exp",
+ "sqrt",
+ "sign",
+ "rint",
+ "abs",
+ "min",
+ "max",
+ "sum",
+ "avg",
+ // functions we define ourselves above
+ "if",
+ "int",
+ "ceil",
+ "cot",
+ "csc",
+ "floor",
+ "sec",
+ "pow",
+ "erfc",
+ "rand",
+ "rand_seed"};
+
+ } // namespace FunctionParser
+
+} // namespace internal
+#endif
+
+
+
+DEAL_II_NAMESPACE_CLOSE
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/base/mu_parser_internal.h>
+#include <deal.II/base/patterns.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_function.h>
+#include <deal.II/base/tensor_function_parser.h>
+#include <deal.II/base/thread_management.h>
+#include <deal.II/base/utilities.h>
+
+#include <boost/random.hpp>
+
+#include <cmath>
+#include <map>
+
+#ifdef DEAL_II_WITH_MUPARSER
+# include <muParser.h>
+#endif
+
+DEAL_II_NAMESPACE_OPEN
+
+
+template <int rank, int dim, typename Number>
+const std::vector<std::string> &
+TensorFunctionParser<rank, dim, Number>::get_expressions() const
+{
+ return expressions;
+}
+
+
+
+template <int rank, int dim, typename Number>
+TensorFunctionParser<rank, dim, Number>::TensorFunctionParser(
+ const double initial_time)
+ : TensorFunction<rank, dim, Number>(initial_time)
+ , initialized(false)
+ , n_vars(0)
+ , n_components(Utilities::pow(dim, rank))
+{}
+
+
+template <int rank, int dim, typename Number>
+TensorFunctionParser<rank, dim, Number>::TensorFunctionParser(
+ const std::string &expression,
+ const std::string &constants,
+ const std::string &variable_names)
+ : TensorFunction<rank, dim, Number>()
+ , initialized(false)
+ , n_vars(0)
+ , n_components(Utilities::pow(dim, rank))
+{
+ auto constants_map = Patterns::Tools::Convert<ConstMap>::to_value(
+ constants,
+ std_cxx14::make_unique<Patterns::Map>(Patterns::Anything(),
+ Patterns::Double(),
+ 0,
+ Patterns::Map::max_int_value,
+ ",",
+ "="));
+ initialize(variable_names,
+ expression,
+ constants_map,
+ Utilities::split_string_list(variable_names, ",").size() ==
+ dim + 1);
+}
+
+
+// We deliberately delay the definition of the default destructor
+// so that we don't need to include the definition of mu::Parser
+// in the header file.
+template <int rank, int dim, typename Number>
+TensorFunctionParser<rank, dim, Number>::~TensorFunctionParser() = default;
+
+
+#ifdef DEAL_II_WITH_MUPARSER
+
+template <int rank, int dim, typename Number>
+void
+TensorFunctionParser<rank, dim, Number>::initialize(
+ const std::string & variables,
+ const std::vector<std::string> & expressions,
+ const std::map<std::string, double> &constants,
+ const bool time_dependent)
+{
+ this->tfp.clear(); // this will reset all thread-local objects
+
+ this->constants = constants;
+ this->var_names = Utilities::split_string_list(variables, ',');
+ this->expressions = expressions;
+ AssertThrow(((time_dependent) ? dim + 1 : dim) == var_names.size(),
+ ExcMessage("Wrong number of variables"));
+
+ // 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()));
+
+ // 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;
+
+ // create a parser object for the current thread we can then query
+ // in value() and vector_value(). this is not strictly necessary
+ // because a user may never call these functions on the current
+ // thread, but it gets us error messages about wrong formulas right
+ // away
+ init_muparser();
+
+ // finally set the initialization bit
+ initialized = true;
+}
+
+
+
+template <int rank, int dim, typename Number>
+void
+TensorFunctionParser<rank, dim, Number>::init_muparser() const
+{
+ // check that we have not already initialized the parser on the
+ // current thread, i.e., that the current function is only called
+ // once per thread
+ Assert(tfp.get().size() == 0, ExcInternalError());
+
+ // initialize the objects for the current thread (tfp.get() and
+ // vars.get())
+ tfp.get().reserve(this->n_components);
+ vars.get().resize(var_names.size());
+ for (unsigned int component = 0; component < this->n_components; ++component)
+ {
+ tfp.get().emplace_back(new mu::Parser());
+
+ for (const auto &constant : constants)
+ {
+ tfp.get()[component]->DefineConst(constant.first, constant.second);
+ }
+
+ for (unsigned int iv = 0; iv < var_names.size(); ++iv)
+ tfp.get()[component]->DefineVar(var_names[iv], &vars.get()[iv]);
+
+ // define some compatibility functions:
+ tfp.get()[component]->DefineFun("if",
+ internal::FunctionParser::mu_if,
+ true);
+ tfp.get()[component]->DefineOprt("|", internal::FunctionParser::mu_or, 1);
+ tfp.get()[component]->DefineOprt("&",
+ internal::FunctionParser::mu_and,
+ 2);
+ tfp.get()[component]->DefineFun("int",
+ internal::FunctionParser::mu_int,
+ true);
+ tfp.get()[component]->DefineFun("ceil",
+ internal::FunctionParser::mu_ceil,
+ true);
+ tfp.get()[component]->DefineFun("cot",
+ internal::FunctionParser::mu_cot,
+ true);
+ tfp.get()[component]->DefineFun("csc",
+ internal::FunctionParser::mu_csc,
+ true);
+ tfp.get()[component]->DefineFun("floor",
+ internal::FunctionParser::mu_floor,
+ true);
+ tfp.get()[component]->DefineFun("sec",
+ internal::FunctionParser::mu_sec,
+ true);
+ tfp.get()[component]->DefineFun("log",
+ internal::FunctionParser::mu_log,
+ true);
+ tfp.get()[component]->DefineFun("pow",
+ internal::FunctionParser::mu_pow,
+ true);
+ tfp.get()[component]->DefineFun("erfc",
+ internal::FunctionParser::mu_erfc,
+ true);
+ tfp.get()[component]->DefineFun("rand_seed",
+ internal::FunctionParser::mu_rand_seed,
+ true);
+ tfp.get()[component]->DefineFun("rand",
+ internal::FunctionParser::mu_rand,
+ true);
+
+ try
+ {
+ // muparser expects that 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 tfparser 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[component];
+
+ for (const auto ¤t_function_name :
+ internal::FunctionParser::function_names)
+ {
+ const unsigned int function_name_length =
+ current_function_name.size();
+
+ std::string::size_type pos = 0;
+ while (true)
+ {
+ // try to find any occurrences of the function name
+ pos = transformed_expression.find(current_function_name, pos);
+ if (pos == std::string::npos)
+ break;
+
+ // replace whitespace until there no longer is any
+ while ((pos + function_name_length <
+ transformed_expression.size()) &&
+ ((transformed_expression[pos + function_name_length] ==
+ ' ') ||
+ (transformed_expression[pos + function_name_length] ==
+ '\t')))
+ transformed_expression.erase(
+ transformed_expression.begin() + pos +
+ function_name_length);
+
+ // move the current search position by the size of the
+ // actual function name
+ pos += function_name_length;
+ }
+ }
+
+ // now use the transformed expression
+ tfp.get()[component]->SetExpr(transformed_expression);
+ }
+ catch (mu::ParserError &e)
+ {
+ std::cerr << "Message: <" << e.GetMsg() << ">\n";
+ std::cerr << "Formula: <" << e.GetExpr() << ">\n";
+ std::cerr << "Token: <" << e.GetToken() << ">\n";
+ std::cerr << "Position: <" << e.GetPos() << ">\n";
+ std::cerr << "Errc: <" << e.GetCode() << ">" << std::endl;
+ AssertThrow(false, ExcParseError(e.GetCode(), e.GetMsg()));
+ }
+ }
+}
+
+
+
+template <int rank, int dim, typename Number>
+void
+TensorFunctionParser<rank, dim, Number>::initialize(
+ const std::string & vars,
+ const std::string & expression,
+ const std::map<std::string, double> &constants,
+ const bool time_dependent)
+{
+ initialize(vars,
+ Utilities::split_string_list(expression, ';'),
+ constants,
+ time_dependent);
+}
+
+
+
+template <int rank, int dim, typename Number>
+Tensor<rank, dim, Number>
+TensorFunctionParser<rank, dim, Number>::value(const Point<dim> &p) const
+{
+ Assert(initialized == true, ExcNotInitialized());
+
+ // initialize the parser if that hasn't happened yet on the current thread
+ if (tfp.get().size() == 0)
+ init_muparser();
+
+ for (unsigned int i = 0; i < dim; ++i)
+ vars.get()[i] = p(i);
+ if (dim != n_vars)
+ vars.get()[dim] = this->get_time();
+
+ // initialize tensor with zeros
+ Tensor<rank, dim, Number> value;
+
+ try
+ {
+ unsigned int component = 0;
+ for (Number *value_ptr = value.begin_raw(); value_ptr != value.end_raw();
+ ++value_ptr)
+ {
+ *value_ptr = tfp.get()[component]->Eval();
+ ++component;
+ } // for
+ } // try
+ catch (mu::ParserError &e)
+ {
+ std::cerr << "Message: <" << e.GetMsg() << ">\n";
+ std::cerr << "Formula: <" << e.GetExpr() << ">\n";
+ std::cerr << "Token: <" << e.GetToken() << ">\n";
+ std::cerr << "Position: <" << e.GetPos() << ">\n";
+ std::cerr << "Errc: <" << e.GetCode() << ">" << std::endl;
+ AssertThrow(false, ExcParseError(e.GetCode(), e.GetMsg()));
+
+ return Tensor<rank, dim, Number>();
+ } // catch
+
+ return value;
+}
+
+
+
+template <int rank, int dim, typename Number>
+void
+TensorFunctionParser<rank, dim, Number>::value_list(
+ const std::vector<Point<dim>> & p,
+ std::vector<Tensor<rank, dim, Number>> &values) const
+{
+ Assert(p.size() == values.size(),
+ ExcDimensionMismatch(p.size(), values.size()));
+
+ for (unsigned int i = 0; i < p.size(); ++i)
+ {
+ values[i] = value(p[i]);
+ }
+}
+
+#else
+
+
+template <int rank, int dim, typename Number>
+void
+TensorFunctionParser<rank, dim, Number>::initialize(
+ const std::string &,
+ const std::vector<std::string> &,
+ const std::map<std::string, double> &,
+ const bool)
+{
+ Assert(false, ExcNeedsFunctionparser());
+}
+
+template <int rank, int dim, typename Number>
+void
+TensorFunctionParser<rank, dim, Number>::initialize(
+ const std::string &,
+ const std::string &,
+ const std::map<std::string, double> &,
+ const bool)
+{
+ Assert(false, ExcNeedsFunctionparser());
+}
+
+
+
+template <int rank, int dim, typename Number>
+Tensor<rank, dim, Number>
+TensorFunctionParser<rank, dim, Number>::value(const Point<dim> &) const
+{
+ Assert(false, ExcNeedsFunctionparser());
+ return Tensor<rank, dim, Number>();
+}
+
+
+
+template <int rank, int dim, typename Number>
+void
+TensorFunctionParser<rank, dim, Number>::value_list(
+ const std::vector<Point<dim>> & p,
+ std::vector<Tensor<rank, dim, Number>> &values) const
+{
+ Assert(false, ExcNeedsFunctionparser());
+}
+
+
+#endif
+
+// explicit instantiations
+#include "tensor_function_parser.inst"
+
+DEAL_II_NAMESPACE_CLOSE
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 - 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+for (S : REAL_SCALARS; rank : RANKS; dim : SPACE_DIMENSIONS)
+ {
+ template class TensorFunctionParser<rank, dim, S>;
+ }
+
+for (S : COMPLEX_SCALARS; rank : RANKS; dim : SPACE_DIMENSIONS)
+ {
+ template class TensorFunctionParser<rank, dim, S>;
+ }
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// This program tests the functionality of the function parser
+// wrapper.
+
+#include <deal.II/base/point.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_function_parser.h>
+
+#include "../tests.h"
+
+
+// Test initializtaion and evaluation of timedependent function in 2d.
+int
+main()
+{
+ initlog();
+
+ // Define some constants that will be used by the function parser
+ std::map<std::string, double> constants;
+ constants["pi"] = numbers::PI;
+
+ // Define the variables that will be used inside the expressions
+ std::string variables = "x,y,t";
+
+ // Define the expressions of the vector_valued function.
+ std::vector<std::string> expressions;
+ expressions.push_back("cos(2*pi*(x*y+t))");
+ expressions.push_back("sin(2*pi*(x*y+t))");
+ expressions.push_back("-sin(2*pi*(x*y+t))");
+ expressions.push_back("cos(2*pi*(x*y+t))");
+
+ // Test time dependent function
+ try
+ {
+ {
+ TensorFunctionParser<2, 2> tensor_function;
+ tensor_function.initialize(variables,
+ expressions,
+ constants,
+ /* time dependent */ true);
+
+ deallog << "Initialize Succeeded with dim = 2, rank = 2, "
+ << expressions.size() << " expressions, " << variables
+ << " as variables." << std::endl;
+ }
+ }
+ catch (...)
+ {
+ deallog << "Initialization or Evaluation Failed with dim = 2, rank = 2, "
+ << expressions.size() << " expressions, " << variables
+ << " as variables." << std::endl;
+ }
+}
--- /dev/null
+
+DEAL::Initialize Succeeded with dim = 2, rank = 2, 4 expressions, x,y,t as variables.
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// This program tests the functionality of the function parser
+// wrapper.
+
+#include <deal.II/base/point.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_function_parser.h>
+
+#include "../tests.h"
+
+
+// Test initializtaion and evaluation of timedependent function in 2d.
+int
+main()
+{
+ initlog();
+
+ // Define some constants that will be used by the function parser
+ std::map<std::string, double> constants;
+ constants["pi"] = numbers::PI;
+
+ // Define the variables that will be used inside the expressions
+ std::string variables = "x,y";
+
+ // Define the expressions of the vector_valued function.
+ std::vector<std::string> expressions;
+ expressions.push_back("cos(2*pi*(x*y))");
+ expressions.push_back("sin(2*pi*(x*y))");
+ expressions.push_back("-sin(2*pi*(x*y))");
+ expressions.push_back("cos(2*pi*(x*y))");
+
+ // Test time dependent function
+ try
+ {
+ {
+ TensorFunctionParser<2, 2> tensor_function;
+ tensor_function.initialize(variables,
+ expressions,
+ constants,
+ /* time dependent */ false);
+
+ deallog << "Initialize Succeeded with dim = 2, rank = 2, "
+ << expressions.size() << " expressions, " << variables
+ << " as variables." << std::endl;
+ }
+ }
+ catch (...)
+ {
+ deallog << "Initialization or Evaluation Failed with dim = 2, rank = 2, "
+ << expressions.size() << " expressions, " << variables
+ << " as variables." << std::endl;
+ }
+}
--- /dev/null
+
+DEAL::Initialize Succeeded with dim = 2, rank = 2, 4 expressions, x,y as variables.
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// This program tests the functionality of the function parser
+// wrapper.
+
+#include <deal.II/base/point.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_function_parser.h>
+
+#include "../tests.h"
+
+
+// Test initializtaion and evaluation of timedependent function in 2d.
+int
+main()
+{
+ initlog();
+
+ // Define some constants that will be used by the function parser
+ std::map<std::string, double> constants;
+ constants["pi"] = numbers::PI;
+
+ // Define the variables that will be used inside the expressions
+ std::string variables = "x,y,t";
+
+ // Define the expressions of the vector_valued function.
+ std::vector<std::string> expressions;
+ expressions.push_back("cos(2*pi*(x*y+t))");
+ expressions.push_back("sin(2*pi*(x*y+t))");
+ expressions.push_back("-sin(2*pi*(x*y+t))");
+ expressions.push_back("cos(2*pi*(x*y+t))");
+
+ Point<2> p(0.5, 0.5); // evaluation point
+
+ // Test time dependent function
+ try
+ {
+ {
+ TensorFunctionParser<2, 2> tensor_function;
+ tensor_function.initialize(variables,
+ expressions,
+ constants,
+ /* time dependent */ true);
+
+ deallog << "Initialize Succeeded with dim = 2, rank = 2, "
+ << expressions.size() << " expressions, " << variables
+ << " as variables." << std::endl;
+
+ for (unsigned int t = 0; t < 5; ++t)
+ {
+ tensor_function.set_time(t);
+ deallog << "Value at p = (0.5,0.5): " << std::endl
+ << tensor_function.value(p)[0][0] << " "
+ << tensor_function.value(p)[0][1] << std::endl
+ << tensor_function.value(p)[1][0] << " "
+ << tensor_function.value(p)[1][1] << " at time t = " << t
+ << std::endl;
+ }
+ }
+ }
+ catch (...)
+ {
+ deallog << "Initialization or Evaluation Failed with dim = 2, rank = 2, "
+ << expressions.size() << " expressions, " << variables
+ << " as variables." << std::endl;
+ }
+}
--- /dev/null
+
+DEAL::Initialize Succeeded with dim = 2, rank = 2, 4 expressions, x,y,t as variables.
+DEAL::Value at p = (0.5,0.5):
+DEAL::6.12323e-17 1.00000
+DEAL::-1.00000 6.12323e-17 at time t = 0
+DEAL::Value at p = (0.5,0.5):
+DEAL::3.06162e-16 1.00000
+DEAL::-1.00000 3.06162e-16 at time t = 1
+DEAL::Value at p = (0.5,0.5):
+DEAL::5.51091e-16 1.00000
+DEAL::-1.00000 5.51091e-16 at time t = 2
+DEAL::Value at p = (0.5,0.5):
+DEAL::-9.80336e-16 1.00000
+DEAL::-1.00000 -9.80336e-16 at time t = 3
+DEAL::Value at p = (0.5,0.5):
+DEAL::-7.35407e-16 1.00000
+DEAL::-1.00000 -7.35407e-16 at time t = 4
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// This program tests the functionality of the function parser
+// wrapper.
+
+#include <deal.II/base/point.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_function_parser.h>
+
+#include "../tests.h"
+
+
+// Test initializtaion and evaluation of timedependent function in 2d.
+int
+main()
+{
+ initlog();
+
+ // Define some constants that will be used by the function parser
+ std::map<std::string, double> constants;
+ constants["pi"] = numbers::PI;
+
+ // Define the variables that will be used inside the expressions
+ std::string variables = "x,y,t";
+
+ // Define the expressions of the vector_valued function.
+ std::vector<std::string> expressions;
+ expressions.push_back("cos(2*pi*(x*y+t+0.5))");
+ expressions.push_back("sin(2*pi*(x*y-t+0.5))");
+
+ Point<2> p(0.5, 0.5); // evaluation point
+
+ // Test time dependent function
+ try
+ {
+ {
+ TensorFunctionParser<1, 2> tensor_function;
+ tensor_function.initialize(variables,
+ expressions,
+ constants,
+ /* time dependent */ true);
+
+ deallog << "Initialize Succeeded with dim = 2, rank = 1, "
+ << expressions.size() << " expressions, " << variables
+ << " as variables." << std::endl;
+
+ for (unsigned int t = 0; t < 5; ++t)
+ {
+ tensor_function.set_time(t);
+ deallog << "Value at p = (0.5,0.5): " << std::endl
+ << tensor_function.value(p)[0] << " "
+ << tensor_function.value(p)[1] << " at time t = " << t
+ << std::endl;
+ }
+ }
+ }
+ catch (...)
+ {
+ deallog << "Initialization or Evaluation Failed with dim = 2, rank = 2, "
+ << expressions.size() << " expressions, " << variables
+ << " as variables." << std::endl;
+ }
+}
--- /dev/null
+
+DEAL::Initialize Succeeded with dim = 2, rank = 1, 2 expressions, x,y,t as variables.
+DEAL::Value at p = (0.5,0.5):
+DEAL::-1.83697e-16 -1.00000 at time t = 0
+DEAL::Value at p = (0.5,0.5):
+DEAL::-4.28626e-16 -1.00000 at time t = 1
+DEAL::Value at p = (0.5,0.5):
+DEAL::-2.44991e-15 -1.00000 at time t = 2
+DEAL::Value at p = (0.5,0.5):
+DEAL::-2.69484e-15 -1.00000 at time t = 3
+DEAL::Value at p = (0.5,0.5):
+DEAL::-2.93977e-15 -1.00000 at time t = 4
projects / dealii.git / commitdiff