From 62870e0933b44f8adc9aa58ac5a5f1241cbc0679 Mon Sep 17 00:00:00 2001
From: Peter Munch <peterrmuench@gmail.com>
Date: Tue, 15 Oct 2024 10:22:34 +0200
Subject: [PATCH] FunctionFromFunctionObjects: accpet single function

---
 doc/news/changes/minor/20241015Munch      |  5 +++
 include/deal.II/base/function.h           | 20 +++++++++--
 include/deal.II/base/function.templates.h | 41 ++++++++++++++++-------
 3 files changed, 51 insertions(+), 15 deletions(-)
 create mode 100644 doc/news/changes/minor/20241015Munch

diff --git a/doc/news/changes/minor/20241015Munch b/doc/news/changes/minor/20241015Munch
new file mode 100644
index 0000000000..faf1a27e47
--- /dev/null
+++ b/doc/news/changes/minor/20241015Munch
@@ -0,0 +1,5 @@
+Improved: FunctionFromFunctionObjects can now accept
+a single function within which individual components need
+to be handled.
+<br>
+(Peter Munch, 2024/10/15)
diff --git a/include/deal.II/base/function.h b/include/deal.II/base/function.h
index 52c8fcca57..135793555d 100644
--- a/include/deal.II/base/function.h
+++ b/include/deal.II/base/function.h
@@ -991,6 +991,20 @@ public:
                 &values,
     const double initial_time = 0.0);
 
+  /**
+   * Constructor for functions of which you only know the values.
+   *
+   * The resulting function will have a number of components equal @p n_components.
+   * A call to the FunctionFromFunctionObject::gradient()
+   * method will trigger an exception, unless you first call the
+   * set_function_gradients() method.
+   */
+  explicit FunctionFromFunctionObjects(
+    const std::function<RangeNumberType(const Point<dim> &, const unsigned int)>
+                      &values,
+    const unsigned int n_components,
+    const double       initial_time = 0.0);
+
   /**
    * Constructor for functions of which you know both the values and the
    * gradients.
@@ -1052,14 +1066,14 @@ private:
   /**
    * The actual function values.
    */
-  std::vector<std::function<RangeNumberType(const Point<dim> &)>>
+  std::function<RangeNumberType(const Point<dim> &, const unsigned int)>
     function_values;
 
   /**
    * The actual function gradients.
    */
-  std::vector<
-    std::function<Tensor<1, dim, RangeNumberType>(const Point<dim> &)>>
+  std::function<Tensor<1, dim, RangeNumberType>(const Point<dim> &,
+                                                const unsigned int)>
     function_gradients;
 };
 
diff --git a/include/deal.II/base/function.templates.h b/include/deal.II/base/function.templates.h
index 1103012fca..d263832b0b 100644
--- a/include/deal.II/base/function.templates.h
+++ b/include/deal.II/base/function.templates.h
@@ -957,8 +957,6 @@ 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)
 {}
 
 
@@ -968,8 +966,20 @@ 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)
+{
+  this->set_function_values(values);
+}
+
+
+
+template <int dim, typename RangeNumberType>
+FunctionFromFunctionObjects<dim, RangeNumberType>::FunctionFromFunctionObjects(
+  const std::function<RangeNumberType(const Point<dim> &, const unsigned int)>
+                    &values,
+  const unsigned int n_components,
+  const double       initial_time)
+  : Function<dim, RangeNumberType>(n_components, initial_time)
   , function_values(values)
-  , function_gradients(values.size())
 {}
 
 
@@ -982,9 +992,10 @@ FunctionFromFunctionObjects<dim, RangeNumberType>::FunctionFromFunctionObjects(
               &gradients,
   const double initial_time)
   : Function<dim, RangeNumberType>(values.size(), initial_time)
-  , function_values(values)
-  , function_gradients(gradients)
-{}
+{
+  this->set_function_values(values);
+  this->set_function_gradients(gradients);
+}
 
 
 
@@ -995,10 +1006,10 @@ FunctionFromFunctionObjects<dim, RangeNumberType>::value(
   const unsigned int component) const
 {
   AssertIndexRange(component, this->n_components);
-  Assert(function_values[component],
+  Assert(function_values,
          ExcMessage("Accessing value() in FunctionFromFunctionObjects requires "
                     "setting the std::function objects for the value"));
-  return function_values[component](p);
+  return function_values(p, component);
 }
 
 
@@ -1010,11 +1021,11 @@ FunctionFromFunctionObjects<dim, RangeNumberType>::gradient(
   const unsigned int component) const
 {
   AssertIndexRange(component, this->n_components);
-  Assert(function_gradients[component],
+  Assert(function_gradients,
          ExcMessage(
            "Accessing gradient() in FunctionFromFunctionObjects "
            "requires setting the std::function objects for the gradient"));
-  return function_gradients[component](p);
+  return function_gradients(p, component);
 }
 
 
@@ -1025,7 +1036,10 @@ 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;
+  function_values = [values](const auto &p, const auto c) {
+    AssertIndexRange(c, values.size());
+    return values[c](p);
+  };
 }
 
 
@@ -1038,7 +1052,10 @@ FunctionFromFunctionObjects<dim, RangeNumberType>::set_function_gradients(
     &gradients)
 {
   AssertDimension(this->n_components, gradients.size());
-  function_gradients = gradients;
+  function_gradients = [gradients](const auto &p, const auto c) {
+    AssertIndexRange(c, gradients.size());
+    return gradients[c](p);
+  };
 }
 
 DEAL_II_NAMESPACE_CLOSE
-- 
2.39.5