]> https://gitweb.dealii.org/ - dealii.git/commitdiff
FunctionFromFunctionObjects. 7918/head
authorLuca Heltai <luca.heltai@sissa.it>
Sun, 14 Apr 2019 12:00:03 +0000 (14:00 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Tue, 23 Apr 2019 13:33:30 +0000 (15:33 +0200)
doc/news/changes/minor/20190412LucaHeltai [new file with mode: 0644]
include/deal.II/base/function.h
include/deal.II/base/function.templates.h
source/base/function.inst.in
tests/base/function_from_function_objects_01.cc [new file with mode: 0644]
tests/base/function_from_function_objects_01.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20190412LucaHeltai b/doc/news/changes/minor/20190412LucaHeltai
new file mode 100644 (file)
index 0000000..b2e28e4
--- /dev/null
@@ -0,0 +1,4 @@
+New: The new FunctionFromFunctionObjects class allows one to wrap a vector of
+std::function objects into a Function object, to allow fast prototyping of user codes. 
+<br>
+(Luca Heltai, 2019/04/12)
index 49c42d623b8ab30f8f690e187718f42178f2bbf8..606d7fa1aa3fd5ed6b36a58c1a930b20f866f8c1 100644 (file)
@@ -762,7 +762,8 @@ private:
  *   return 1.0;
  * }
  *
- * VectorFunctionFromScalarFunctionObject<2> component_mask(&one, 1, 3);
+ * VectorFunctionFromScalarFunctionObject<2, RangeNumberType>
+ *   component_mask(&one, 1, 3);
  * @endcode
  * Here, <code>component_mask</code> then represents a Function object that
  * for every point returns the vector $(0, 1, 0)^T$, i.e. a mask function that
@@ -827,6 +828,143 @@ private:
 };
 
 
+/**
+ * This class is similar to the ScalarFunctionFromFunctionObject and
+ * VectorFunctionFromFunctionObject classes in that it allows for the easy
+ * conversion of a vector of function objects to something that satisfies the
+ * interface of the Function base class.
+ *
+ * The difference is that here the Function object generated may be vector
+ * valued, and you can specify the gradients of the function. The number of
+ * vector components is deduced from the size of the vector in the constructor.
+ *
+ * To be more concrete, let us consider the following example:
+ * @code
+ * RangeNumberType
+ * first_component(const Point<2> &p)
+ * {
+ *   return 1.0;
+ * }
+ *
+ * RangeNumberType
+ * second_component(const Point<2> &p)
+ * {
+ *   return 2.0;
+ * }
+ *
+ * Tensor<1, 2, RangeNumberType>
+ * zero_gradient(const Point<2> &) {
+ *   return Tensor<1, 2, RangeNumberType>();
+ * }
+ *
+ * FunctionFromFunctionObjects<2, RangeNumberType>
+ *     custom_function({&first_component, &second_component},
+ *                     {&zero_gradient, &zero_gradient});
+ * @endcode
+ *
+ * @author Luca Heltai, 2019
+ */
+template <int dim, typename RangeNumberType = double>
+class FunctionFromFunctionObjects : public Function<dim, RangeNumberType>
+{
+public:
+  /**
+   * Default constructor.
+   *
+   * This constructor does not initialize the internal methods. To have a
+   * usable function, you need to call at least the set_function_values()
+   * method. If you need also the gradients of the solution, then you must
+   * also call the set_function_gradients() method.
+   */
+  explicit FunctionFromFunctionObjects(const unsigned int n_components = 1,
+                                       const double       initial_time = 0);
+
+  /**
+   * Constructor for functions of which you only know the values.
+   *
+   * The resulting function will have a number of components equal to the size
+   * of the vector @p values. A call to the FunctionFromFunctionObject::gradient()
+   * method will trigger an exception, unless you first call the
+   * set_function_gradients() method.
+   */
+  FunctionFromFunctionObjects(
+    const std::vector<std::function<RangeNumberType(const Point<dim> &)>>
+      &          values,
+    const double initial_time = 0.0);
+
+  /**
+   * Constructor for functions of which you know both the values and the
+   * gradients.
+   *
+   * The resulting function will have a number of components equal to the size
+   * of the vector @p values. If the size of @p values and @p gradients does not
+   * match, an exception is triggered.
+   */
+  FunctionFromFunctionObjects(
+    const std::vector<std::function<RangeNumberType(const Point<dim> &)>>
+      &values,
+    const std::vector<
+      std::function<Tensor<1, dim, RangeNumberType>(const Point<dim> &)>>
+      &          gradients,
+    const double initial_time = 0.0);
+
+
+  /**
+   * Return the value of the function at the given point. Unless there is only
+   * one component (i.e. the function is scalar), you should state the
+   * component you want to have evaluated; it defaults to zero, i.e. the first
+   * component.
+   */
+  virtual RangeNumberType
+  value(const Point<dim> &p, const unsigned int component = 0) const override;
+
+  /**
+   * Return the gradient of the function at the given point. Unless there is
+   * only one component (i.e. the function is scalar), you should state the
+   * component you want to have evaluated; it defaults to zero, i.e. the first
+   * component.
+   */
+  virtual Tensor<1, dim, RangeNumberType>
+  gradient(const Point<dim> & p,
+           const unsigned int component = 0) const override;
+
+  /**
+   * Reset the function values of this object. An assertion is thrown if the
+   * size of the @p values parameter does not match the number of components of
+   * this object.
+   */
+  void
+  set_function_values(
+    const std::vector<std::function<RangeNumberType(const Point<dim> &)>>
+      &values);
+
+  /**
+   * Reset the function gradients of this object. An assertion is thrown if the
+   * size of the @p gradients parameter does not match the number of components of
+   * this object.
+   */
+  void
+  set_function_gradients(
+    const std::vector<
+      std::function<Tensor<1, dim, RangeNumberType>(const Point<dim> &)>>
+      &gradients);
+
+private:
+  /**
+   * The actual function values.
+   */
+  std::vector<std::function<RangeNumberType(const Point<dim> &)>>
+    function_values;
+
+  /**
+   * The actual function gradients.
+   */
+  std::vector<
+    std::function<Tensor<1, dim, RangeNumberType>(const Point<dim> &)>>
+    function_gradients;
+};
+
+
 /**
  * This class is built as a means of translating the <code>Tensor<1,dim,
  * RangeNumberType> </code> values produced by objects of type TensorFunction
index 3cab781f8188b9b8f45c835a50bc07335aa39a20..8fd0d49d08acbd531909453828b00e9d0ef19998 100644 (file)
@@ -783,6 +783,96 @@ VectorFunctionFromTensorFunction<dim, RangeNumberType>::vector_value_list(
 }
 
 
+
+template <int dim, typename RangeNumberType>
+FunctionFromFunctionObjects<dim, RangeNumberType>::FunctionFromFunctionObjects(
+  const unsigned int n_components,
+  const double       initial_time)
+  : Function<dim, RangeNumberType>(n_components, initial_time)
+  , function_values(n_components)
+  , function_gradients(n_components)
+{}
+
+
+
+template <int dim, typename RangeNumberType>
+FunctionFromFunctionObjects<dim, RangeNumberType>::FunctionFromFunctionObjects(
+  const std::vector<std::function<RangeNumberType(const Point<dim> &)>> &values,
+  const double initial_time)
+  : Function<dim, RangeNumberType>(values.size(), initial_time)
+  , function_values(values)
+  , function_gradients(values.size())
+{}
+
+
+
+template <int dim, typename RangeNumberType>
+FunctionFromFunctionObjects<dim, RangeNumberType>::FunctionFromFunctionObjects(
+  const std::vector<std::function<RangeNumberType(const Point<dim> &)>> &values,
+  const std::vector<
+    std::function<Tensor<1, dim, RangeNumberType>(const Point<dim> &)>>
+    &          gradients,
+  const double initial_time)
+  : Function<dim, RangeNumberType>(values.size(), initial_time)
+  , function_values(values)
+  , function_gradients(gradients)
+{}
+
+
+
+template <int dim, typename RangeNumberType>
+RangeNumberType
+FunctionFromFunctionObjects<dim, RangeNumberType>::value(
+  const Point<dim> & p,
+  const unsigned int component) const
+{
+  AssertIndexRange(component, this->n_components);
+  Assert(function_values[component],
+         ExcMessage("Accessing value() in FunctionFromFunctionObjects requires "
+                    "setting the std::function objects for the value"));
+  return function_values[component](p);
+}
+
+
+
+template <int dim, typename RangeNumberType>
+Tensor<1, dim, RangeNumberType>
+FunctionFromFunctionObjects<dim, RangeNumberType>::gradient(
+  const Point<dim> & p,
+  const unsigned int component) const
+{
+  AssertIndexRange(component, this->n_components);
+  Assert(function_gradients[component],
+         ExcMessage(
+           "Accessing gradient() in FunctionFromFunctionObjects "
+           "requires setting the std::function objects for the gradient"));
+  return function_gradients[component](p);
+}
+
+
+
+template <int dim, typename RangeNumberType>
+void
+FunctionFromFunctionObjects<dim, RangeNumberType>::set_function_values(
+  const std::vector<std::function<RangeNumberType(const Point<dim> &)>> &values)
+{
+  AssertDimension(this->n_components, values.size());
+  function_values = values;
+}
+
+
+
+template <int dim, typename RangeNumberType>
+void
+FunctionFromFunctionObjects<dim, RangeNumberType>::set_function_gradients(
+  const std::vector<
+    std::function<Tensor<1, dim, RangeNumberType>(const Point<dim> &)>>
+    &gradients)
+{
+  AssertDimension(this->n_components, gradients.size());
+  function_gradients = gradients;
+}
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif /* dealii_function_templates_h */
index 3efdc972bf846ca904dda1ea9f635009d4ee5831..b9177e0a4cc92cbd27af3832a46b60bfbad8191b 100644 (file)
@@ -22,6 +22,7 @@ for (S : REAL_SCALARS; dim : SPACE_DIMENSIONS)
     template class ComponentSelectFunction<dim, S>;
     template class ScalarFunctionFromFunctionObject<dim, S>;
     template class VectorFunctionFromScalarFunctionObject<dim, S>;
+    template class FunctionFromFunctionObjects<dim, S>;
     template class VectorFunctionFromTensorFunction<dim, S>;
   }
 
@@ -33,5 +34,6 @@ for (S : COMPLEX_SCALARS; dim : SPACE_DIMENSIONS)
     template class ComponentSelectFunction<dim, S>;
     template class ScalarFunctionFromFunctionObject<dim, S>;
     template class VectorFunctionFromScalarFunctionObject<dim, S>;
+    template class FunctionFromFunctionObjects<dim, S>;
     template class VectorFunctionFromTensorFunction<dim, S>;
   }
diff --git a/tests/base/function_from_function_objects_01.cc b/tests/base/function_from_function_objects_01.cc
new file mode 100644 (file)
index 0000000..5fb9674
--- /dev/null
@@ -0,0 +1,66 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Check FunctionFromFunctionObjects implementation
+
+#include <deal.II/base/function.h>
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/lac/vector.h>
+
+#include <boost/core/demangle.hpp>
+
+#include <string>
+#include <typeinfo>
+#include <vector>
+
+#include "../tests.h"
+
+// Test the lambda function object.
+
+int
+main()
+{
+  initlog();
+
+  // Scalar function:
+  auto val1 = [](const Point<2> &p) { return p.square(); };
+  auto val2 = [](const Point<2> &p) { return 2 * p[0]; };
+
+  auto grad1 = [](const Point<2> &p) { return 2 * Tensor<1, 2>(p); };
+  auto grad2 = [](const Point<2> &) { return 2 * Tensor<1, 2>({1, 0}); };
+
+  FunctionFromFunctionObjects<2> fun1({val1});
+  FunctionFromFunctionObjects<2> fun2({val1, val2});
+
+  FunctionFromFunctionObjects<2> fun3({val1}, {grad1});
+  FunctionFromFunctionObjects<2> fun4({val1, val2}, {grad1, grad2});
+
+
+  Point<2> p(1, 2);
+
+  deallog << fun1.value(p) << std::endl;
+
+  Vector<double> v(2);
+  fun2.vector_value(p, v);
+  deallog << v(0) << " " << v(1) << std::endl;
+
+  deallog << fun3.gradient(p) << std::endl;
+
+  std::vector<Tensor<1, 2>> g(2);
+  fun4.vector_gradient(p, g);
+  deallog << g[0] << ", " << g[1] << std::endl;
+}
diff --git a/tests/base/function_from_function_objects_01.output b/tests/base/function_from_function_objects_01.output
new file mode 100644 (file)
index 0000000..219c3b5
--- /dev/null
@@ -0,0 +1,5 @@
+
+DEAL::5.00000
+DEAL::5.00000 2.00000
+DEAL::2.00000 4.00000
+DEAL::2.00000 4.00000, 2.00000 0.00000

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.