From e667d8548bbd8bacedeb1fb428722748b2e2132e Mon Sep 17 00:00:00 2001 From: David Wells Date: Sat, 21 May 2022 15:04:20 -0400 Subject: [PATCH] Move most of the muParser implementation to the new base class. --- include/deal.II/base/function_parser.h | 20 -- include/deal.II/base/mu_parser_internal.h | 39 ++- include/deal.II/base/tensor_function_parser.h | 10 - source/base/CMakeLists.txt | 1 + source/base/function_parser.cc | 212 +--------------- source/base/mu_parser_internal.cc | 232 +++++++++++++++++- source/base/mu_parser_internal.inst.in | 25 ++ source/base/tensor_function_parser.cc | 184 +------------- 8 files changed, 286 insertions(+), 437 deletions(-) create mode 100644 source/base/mu_parser_internal.inst.in diff --git a/include/deal.II/base/function_parser.h b/include/deal.II/base/function_parser.h index 396f6ae008..4564072f3e 100644 --- a/include/deal.II/base/function_parser.h +++ b/include/deal.II/base/function_parser.h @@ -344,16 +344,6 @@ public: virtual double value(const Point &p, const unsigned int component = 0) const override; - /** - * Return all components of a vector-valued function at the given point @p - * p. - * - * values shall have the right size beforehand, i.e. - * #n_components. - */ - virtual void - vector_value(const Point &p, Vector &values) const override; - /** * Return an array of function expressions (one per component), used to * initialize this function. @@ -379,16 +369,6 @@ public: << ")."); //@} - -private: - /** - * Initialize fp 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 'fp.get().size()==0' (not - * initialized) or >0 (already initialized). - */ - void - init_muparser() const; }; diff --git a/include/deal.II/base/mu_parser_internal.h b/include/deal.II/base/mu_parser_internal.h index a4164ed2ac..562cc6025e 100644 --- a/include/deal.II/base/mu_parser_internal.h +++ b/include/deal.II/base/mu_parser_internal.h @@ -21,19 +21,16 @@ #include +#include +#include #include #include #include #include - DEAL_II_NAMESPACE_OPEN - - -#ifdef DEAL_II_WITH_MUPARSER - namespace internal { namespace FunctionParser @@ -94,6 +91,19 @@ namespace internal */ std::vector get_function_names(); + + /** + * @addtogroup Exceptions + * @{ + */ + DeclException2(ExcParseError, + int, + std::string, + << "Parsing Error at Column " << arg1 + << ". The parser said: " << arg2); + + //@} + /** * deal.II uses muParser as a purely internal dependency. To this end, we do * not include any muParser headers in our own headers (and the bundled @@ -160,6 +170,19 @@ namespace internal mutable Threads::ThreadLocalStorage parser_data; + void + init_muparser() const; + + Number + do_value(const Point &p, + const double time, + unsigned int component) const; + + void + do_all_values(const Point & p, + const double time, + ArrayView &values) const; + /** * An array to keep track of all the constants, required to initialize fp * in each thread. @@ -193,14 +216,8 @@ namespace internal */ unsigned int n_vars; }; - - } // namespace FunctionParser - } // namespace internal -#endif - - DEAL_II_NAMESPACE_CLOSE diff --git a/include/deal.II/base/tensor_function_parser.h b/include/deal.II/base/tensor_function_parser.h index 2a757524c9..1045c106ff 100644 --- a/include/deal.II/base/tensor_function_parser.h +++ b/include/deal.II/base/tensor_function_parser.h @@ -261,17 +261,7 @@ public: << ")."); //@} - private: - /** - * 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; - /** * Number of components is equal dimrank. */ diff --git a/source/base/CMakeLists.txt b/source/base/CMakeLists.txt index 78190e094a..0e58d1df6a 100644 --- a/source/base/CMakeLists.txt +++ b/source/base/CMakeLists.txt @@ -148,6 +148,7 @@ SET(_inst symbolic_function.inst.in tensor_function.inst.in tensor_function_parser.inst.in + mu_parser_internal.inst.in time_stepping.inst.in ) diff --git a/source/base/function_parser.cc b/source/base/function_parser.cc index b64138b88d..30736b08a8 100644 --- a/source/base/function_parser.cc +++ b/source/base/function_parser.cc @@ -25,10 +25,6 @@ #include #include -#ifdef DEAL_II_WITH_MUPARSER -# include -#endif - DEAL_II_NAMESPACE_OPEN @@ -111,145 +107,12 @@ FunctionParser::initialize(const std::string & variables, // 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(); + this->init_muparser(); // finally set the initialization bit this->initialized = true; } -namespace internal -{ - namespace FunctionParserImplementation - { - namespace - { - /** - * PIMPL for mu::Parser. - */ - class Parser : public internal::FunctionParser::muParserBase - { - public: - operator mu::Parser &() - { - return parser; - } - - operator const mu::Parser &() const - { - return parser; - } - - protected: - mu::Parser parser; - }; - } // namespace - } // namespace FunctionParserImplementation -} // namespace internal - - -template -void -FunctionParser::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 - internal::FunctionParser::ParserData &data = this->parser_data.get(); - Assert(data.parsers.size() == 0 && data.vars.size() == 0, ExcInternalError()); - - // initialize the objects for the current thread - data.parsers.reserve(this->n_components); - data.vars.resize(this->var_names.size()); - for (unsigned int component = 0; component < this->n_components; ++component) - { - data.parsers.emplace_back( - new internal::FunctionParserImplementation::Parser()); - mu::Parser &parser = - dynamic_cast( - *data.parsers.back()); - - for (const auto &constant : this->constants) - parser.DefineConst(constant.first, constant.second); - - for (unsigned int iv = 0; iv < this->var_names.size(); ++iv) - parser.DefineVar(this->var_names[iv], &data.vars[iv]); - - // define some compatibility functions: - parser.DefineFun("if", internal::FunctionParser::mu_if, true); - parser.DefineOprt("|", internal::FunctionParser::mu_or, 1); - parser.DefineOprt("&", internal::FunctionParser::mu_and, 2); - parser.DefineFun("int", internal::FunctionParser::mu_int, true); - parser.DefineFun("ceil", internal::FunctionParser::mu_ceil, true); - parser.DefineFun("cot", internal::FunctionParser::mu_cot, true); - parser.DefineFun("csc", internal::FunctionParser::mu_csc, true); - parser.DefineFun("floor", internal::FunctionParser::mu_floor, true); - parser.DefineFun("sec", internal::FunctionParser::mu_sec, true); - parser.DefineFun("log", internal::FunctionParser::mu_log, true); - parser.DefineFun("pow", internal::FunctionParser::mu_pow, true); - parser.DefineFun("erf", internal::FunctionParser::mu_erf, true); - parser.DefineFun("erfc", internal::FunctionParser::mu_erfc, true); - parser.DefineFun("rand_seed", - internal::FunctionParser::mu_rand_seed, - true); - parser.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 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 = this->expressions[component]; - - for (const auto ¤t_function_name : - internal::FunctionParser::get_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 - parser.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 @@ -272,78 +135,7 @@ double FunctionParser::value(const Point & p, const unsigned int component) const { - Assert(this->initialized == true, ExcNotInitialized()); - AssertIndexRange(component, this->n_components); - - // initialize the parser if that hasn't happened yet on the current thread - internal::FunctionParser::ParserData &data = this->parser_data.get(); - if (data.vars.size() == 0) - init_muparser(); - - for (unsigned int i = 0; i < dim; ++i) - data.vars[i] = p(i); - if (dim != this->n_vars) - data.vars[dim] = this->get_time(); - - try - { - Assert(dynamic_cast( - data.parsers[component].get()), - ExcInternalError()); - // using dynamic_cast in the next line is about 6% slower than - // static_cast, so use the assertion above for debugging and disable - // clang-tidy: - mu::Parser &parser = - static_cast( // NOLINT - *data.parsers[component]); - return parser.Eval(); - } - 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 0.0; - } -} - - - -template -void -FunctionParser::vector_value(const Point &p, - Vector & values) const -{ - Assert(this->initialized == true, ExcNotInitialized()); - Assert(values.size() == this->n_components, - ExcDimensionMismatch(values.size(), this->n_components)); - - - // initialize the parser if that hasn't happened yet on the current thread - internal::FunctionParser::ParserData &data = this->parser_data.get(); - if (data.vars.size() == 0) - init_muparser(); - - for (unsigned int i = 0; i < dim; ++i) - data.vars[i] = p(i); - if (dim != this->n_vars) - data.vars[dim] = this->get_time(); - - for (unsigned int component = 0; component < this->n_components; ++component) - { - // Same comment in value() applies here too: - Assert(dynamic_cast( - data.parsers[component].get()), - ExcInternalError()); - mu::Parser &parser = - static_cast( // NOLINT - *data.parsers[component]); - - values(component) = parser.Eval(); - } + return this->do_value(p, this->get_time(), component); } #else diff --git a/source/base/mu_parser_internal.cc b/source/base/mu_parser_internal.cc index 364711d1ff..93c6470bec 100644 --- a/source/base/mu_parser_internal.cc +++ b/source/base/mu_parser_internal.cc @@ -23,13 +23,12 @@ #include #include +#ifdef DEAL_II_WITH_MUPARSER +# include +#endif DEAL_II_NAMESPACE_OPEN - - -#ifdef DEAL_II_WITH_MUPARSER - namespace internal { namespace FunctionParser @@ -200,11 +199,232 @@ namespace internal "rand_seed"}; } - } // namespace FunctionParser +#ifdef DEAL_II_WITH_MUPARSER + /** + * PIMPL for mu::Parser. + */ + class Parser : public muParserBase + { + public: + operator mu::Parser &() + { + return parser; + } + + operator const mu::Parser &() const + { + return parser; + } + + protected: + mu::Parser parser; + }; +#endif -} // namespace internal + + + template + void + ParserImplementation::init_muparser() const + { +#ifdef DEAL_II_WITH_MUPARSER + // 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 + ParserData &data = this->parser_data.get(); + Assert(data.parsers.size() == 0 && data.vars.size() == 0, + ExcInternalError()); + const unsigned int n_components = expressions.size(); + + // initialize the objects for the current thread + data.parsers.reserve(n_components); + data.vars.resize(this->var_names.size()); + for (unsigned int component = 0; component < n_components; ++component) + { + data.parsers.emplace_back(std::make_unique()); + mu::Parser &parser = dynamic_cast(*data.parsers.back()); + + for (const auto &constant : this->constants) + parser.DefineConst(constant.first, constant.second); + + for (unsigned int iv = 0; iv < this->var_names.size(); ++iv) + parser.DefineVar(this->var_names[iv], &data.vars[iv]); + + // define some compatibility functions: + parser.DefineFun("if", mu_if, true); + parser.DefineOprt("|", mu_or, 1); + parser.DefineOprt("&", mu_and, 2); + parser.DefineFun("int", mu_int, true); + parser.DefineFun("ceil", mu_ceil, true); + parser.DefineFun("cot", mu_cot, true); + parser.DefineFun("csc", mu_csc, true); + parser.DefineFun("floor", mu_floor, true); + parser.DefineFun("sec", mu_sec, true); + parser.DefineFun("log", mu_log, true); + parser.DefineFun("pow", mu_pow, true); + parser.DefineFun("erfc", mu_erfc, true); + parser.DefineFun("rand_seed", mu_rand_seed, true); + parser.DefineFun("rand", 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 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 = this->expressions[component]; + + for (const auto ¤t_function_name : get_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 + parser.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())); + } + } +#else + AssertThrow(false, ExcNeedsFunctionparser()); +#endif + } + + template + Number + ParserImplementation::do_value(const Point &p, + const double time, + unsigned int component) const + { +#ifdef DEAL_II_WITH_MUPARSER + Assert(this->initialized == true, ExcNotInitialized()); + + // initialize the parser if that hasn't happened yet on the current + // thread + internal::FunctionParser::ParserData &data = this->parser_data.get(); + if (data.vars.size() == 0) + init_muparser(); + + for (unsigned int i = 0; i < dim; ++i) + data.vars[i] = p(i); + if (dim != this->n_vars) + data.vars[dim] = time; + + try + { + Assert(dynamic_cast(data.parsers[component].get()), + ExcInternalError()); + // NOLINTNEXTLINE don't warn about using static_cast once we check + mu::Parser &parser = static_cast(*data.parsers[component]); + return parser.Eval(); + } // 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())); + } // catch + +#else + AssertThrow(false, ExcNeedsFunctionparser()); +#endif + return std::numeric_limits::signaling_NaN(); + } + + template + void + ParserImplementation::do_all_values( + const Point & p, + const double time, + ArrayView &values) const + { +#ifdef DEAL_II_WITH_MUPARSER + Assert(this->initialized == true, ExcNotInitialized()); + + // initialize the parser if that hasn't happened yet on the current + // thread + internal::FunctionParser::ParserData &data = this->parser_data.get(); + if (data.vars.size() == 0) + init_muparser(); + + for (unsigned int i = 0; i < dim; ++i) + data.vars[i] = p(i); + if (dim != this->n_vars) + data.vars[dim] = time; + + AssertDimension(values.size(), data.parsers.size()); + try + { + for (unsigned int component = 0; component < data.parsers.size(); + ++component) + { + Assert(dynamic_cast(data.parsers[component].get()), + ExcInternalError()); + mu::Parser &parser = + // We just checked that the pointer is valid so suppress the + // clang-tidy check + static_cast(*data.parsers[component]); // NOLINT + values[component] = parser.Eval(); + } + } // 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())); + } // catch +#else + AssertThrow(false, ExcNeedsFunctionparser()); #endif + } +// explicit instantiations +#include "mu_parser_internal.inst" + } // namespace FunctionParser +} // namespace internal DEAL_II_NAMESPACE_CLOSE diff --git a/source/base/mu_parser_internal.inst.in b/source/base/mu_parser_internal.inst.in new file mode 100644 index 0000000000..9993b37291 --- /dev/null +++ b/source/base/mu_parser_internal.inst.in @@ -0,0 +1,25 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2022 - 2022 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; dim : SPACE_DIMENSIONS) + { + template class ParserImplementation; + } + +for (S : COMPLEX_SCALARS; dim : SPACE_DIMENSIONS) + { + template class ParserImplementation; + } diff --git a/source/base/tensor_function_parser.cc b/source/base/tensor_function_parser.cc index cb20799ef1..3d039cdbbe 100644 --- a/source/base/tensor_function_parser.cc +++ b/source/base/tensor_function_parser.cc @@ -26,10 +26,6 @@ #include #include -#ifdef DEAL_II_WITH_MUPARSER -# include -#endif - DEAL_II_NAMESPACE_OPEN @@ -123,7 +119,7 @@ TensorFunctionParser::initialize( // 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(); + this->init_muparser(); // finally set the initialization bit this->initialized = true; @@ -131,142 +127,6 @@ TensorFunctionParser::initialize( -namespace internal -{ - namespace TensorFunctionParserImplementation - { - namespace - { - /** - * PIMPL for mu::Parser. - */ - class Parser : public internal::FunctionParser::muParserBase - { - public: - operator mu::Parser &() - { - return parser; - } - - operator const mu::Parser &() const - { - return parser; - } - - protected: - mu::Parser parser; - }; - } // namespace - } // namespace TensorFunctionParserImplementation -} // namespace internal - - - -template -void -TensorFunctionParser::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 - internal::FunctionParser::ParserData &data = this->parser_data.get(); - Assert(data.parsers.size() == 0 && data.vars.size() == 0, ExcInternalError()); - - // initialize the objects for the current thread - data.parsers.reserve(this->n_components); - data.vars.resize(this->var_names.size()); - for (unsigned int component = 0; component < this->n_components; ++component) - { - data.parsers.emplace_back( - std::make_unique< - internal::TensorFunctionParserImplementation::Parser>()); - mu::Parser &parser = - dynamic_cast( - *data.parsers.back()); - - for (const auto &constant : this->constants) - parser.DefineConst(constant.first, constant.second); - - for (unsigned int iv = 0; iv < this->var_names.size(); ++iv) - parser.DefineVar(this->var_names[iv], &data.vars[iv]); - - // define some compatibility functions: - parser.DefineFun("if", internal::FunctionParser::mu_if, true); - parser.DefineOprt("|", internal::FunctionParser::mu_or, 1); - parser.DefineOprt("&", internal::FunctionParser::mu_and, 2); - parser.DefineFun("int", internal::FunctionParser::mu_int, true); - parser.DefineFun("ceil", internal::FunctionParser::mu_ceil, true); - parser.DefineFun("cot", internal::FunctionParser::mu_cot, true); - parser.DefineFun("csc", internal::FunctionParser::mu_csc, true); - parser.DefineFun("floor", internal::FunctionParser::mu_floor, true); - parser.DefineFun("sec", internal::FunctionParser::mu_sec, true); - parser.DefineFun("log", internal::FunctionParser::mu_log, true); - parser.DefineFun("pow", internal::FunctionParser::mu_pow, true); - parser.DefineFun("erfc", internal::FunctionParser::mu_erfc, true); - parser.DefineFun("rand_seed", - internal::FunctionParser::mu_rand_seed, - true); - parser.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 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 = this->expressions[component]; - - for (const auto ¤t_function_name : - internal::FunctionParser::get_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 - parser.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 void TensorFunctionParser::initialize( @@ -287,46 +147,10 @@ template Tensor TensorFunctionParser::value(const Point &p) const { - Assert(this->initialized == true, ExcNotInitialized()); - - // initialize the parser if that hasn't happened yet on the current thread - internal::FunctionParser::ParserData &data = this->parser_data.get(); - if (data.vars.size() == 0) - init_muparser(); - - for (unsigned int i = 0; i < dim; ++i) - data.vars[i] = p(i); - if (dim != this->n_vars) - data.vars[dim] = this->get_time(); - std::array::n_independent_components> - values; - - try - { - for (unsigned int component = 0; component < values.size(); ++component) - { - Assert(dynamic_cast< - internal::TensorFunctionParserImplementation::Parser *>( - data.parsers[component].get()), - ExcInternalError()); - mu::Parser &parser = static_cast< // NOLINT - internal::TensorFunctionParserImplementation::Parser &>( - *data.parsers[component]); - values[component] = parser.Eval(); - } // 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(); - } // catch + values; + auto values_view = make_array_view(values.begin(), values.end()); + this->do_all_values(p, this->get_time(), values_view); return Tensor( make_array_view(values.begin(), values.end())); -- 2.39.5