ParsedFunction(const unsigned int n_components = 1, const double h = 1e-8);
/**
- * Declare parameters needed by this class. The additional parameter @p
+ * Declare parameters needed by this class. The parameter @p
* n_components is used to generate the right code according to the number
* of components of the function that will parse this ParameterHandler. If
* the number of components which is parsed does not match the number of
* components of this object, an assertion is thrown and the program is
- * aborted. The default behavior for this class is to declare the
- * following entries:
+ * aborted. The additional parameter @p expr is used to declare the default
+ * expression used by the function. The default behavior for this class is
+ * to declare the following entries:
*
* @code
*
*/
static void
declare_parameters(ParameterHandler &prm,
- const unsigned int n_components = 1);
+ const unsigned int n_components = 1,
+ const std::string &input_expr = "");
/**
* Parse parameters needed by this class. If the number of components
template <int dim>
void
ParsedFunction<dim>::declare_parameters(ParameterHandler &prm,
- const unsigned int n_components)
+ const unsigned int n_components,
+ const std::string &input_expr)
{
Assert(n_components > 0, ExcZero());
"variable names in your function expression.");
// The expression of the function
- std::string expr = "0";
- for (unsigned int i = 1; i < n_components; ++i)
- expr += "; 0";
+ // If the string is an empty string, 0 is set for each components.
+ std::string expr = input_expr;
+ if (expr == "")
+ {
+ expr = "0";
+ for (unsigned int i = 1; i < n_components; ++i)
+ expr += "; 0";
+ }
+ else
+ {
+ // If the user specified an input expr, the number of component
+ // specified need to match n_components.
+ AssertDimension((std::count(expr.begin(), expr.end(), ';') + 1) , n_components);
+ }
+
prm.declare_entry(
"Function expression",
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2010 - 2023 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// Part of the source code is dual licensed under Apache-2.0 WITH
+// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
+// governing the source code and code contributions can be found in
+// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
+//
+// ------------------------------------------------------------------------
+
+
+// This program tests the functionality of imposing a default expression
+// when declaring a ParsedFunction.
+
+#include <deal.II/base/parsed_function.h>
+#include <deal.II/base/parameter_handler.h>
+
+#include "../tests.h"
+#include "deal.II/matrix_free/matrix_free.h"
+
+
+void
+test() {
+ // set up problem:
+ ParameterHandler prm;
+ Functions::ParsedFunction<3> fct_1(1);
+ Functions::ParsedFunction<3> fct_2(2);
+ Functions::ParsedFunction<3> fct_3(3);
+
+ std::string expression_1 = "1000";
+ std::string expression_3 = "10;10;10";
+
+ // Declare
+ prm.enter_subsection("Function 1");
+ fct_1.declare_parameters(prm, 1, expression_1);
+ prm.leave_subsection();
+
+ prm.enter_subsection("Function 2");
+ fct_2.declare_parameters(prm, 2);
+ prm.leave_subsection();
+
+ prm.enter_subsection("Function 3");
+ fct_3.declare_parameters(prm, 3, expression_3);
+ prm.leave_subsection();
+
+
+ // Parse
+ prm.enter_subsection("Function 1");
+ fct_1.parse_parameters(prm);
+ prm.leave_subsection();
+
+ prm.enter_subsection("Function 2");
+ fct_2.parse_parameters(prm);
+ prm.leave_subsection();
+
+ prm.enter_subsection("Function 3");
+ fct_3.parse_parameters(prm);
+ prm.leave_subsection();
+
+
+ // Dummy point
+ const Point<3> p;
+
+ const double result_1 = fct_1.value(p);
+ Vector<double> result_2(2);
+ Vector<double> result_3(3);
+
+ fct_2.vector_value(p,result_2);
+ fct_3.vector_value(p,result_3);
+
+ deallog << "Function 1 '" << expression_1 << "' is "
+ << result_1 << std::endl;
+
+ deallog << "Function 2 '' is "
+ << result_2 << std::endl;
+
+ deallog << "Function 3 '" << expression_3 << "' is "
+ << result_3 << std::endl;
+}
+
+int
+main() {
+ initlog();
+
+ test();
+}