]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add a function that converts a lambda function into a vector Function object.
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 19 Jun 2025 00:06:18 +0000 (18:06 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 19 Jun 2025 00:32:35 +0000 (18:32 -0600)
include/deal.II/base/function.h
include/deal.II/base/function.templates.h
source/base/function.inst.in

index 2f555f2ff3a80cdf41710c1bbc89614f30a88d2d..acb9cc0072b496d357ea74f9ed841218e74fcb24 100644 (file)
@@ -1126,12 +1126,12 @@ public:
    * By default, create a Vector object of the same size as
    * <tt>tensor_function</tt> returns, i.e., with <tt>dim</tt> components.
    *
-   * @param tensor_function The TensorFunction that will form one component of
+   * @param tensor_function The TensorFunction that will form `dim` components of
    * the resulting Vector Function object.
    * @param n_components The total number of vector components of the
-   * resulting TensorFunction object.
+   * resulting TensorFunction object. This clearly has to be at least `dim`.
    * @param selected_component The first component that should be filled by
-   * the first argument.  This should be such that the entire tensor_function
+   * the first argument.  This should be such that the entire `tensor_function`
    * fits inside the <tt>n_component</tt> length return vector.
    */
   explicit VectorFunctionFromTensorFunction(
@@ -1243,6 +1243,118 @@ private:
 };
 
 
+/**
+ * This class is built as a means of translating the <code>Tensor<1,dim,
+ * RangeNumberType> </code> values produced by function objects that
+ * for a given point return a tensor,
+ * and returning them as a multiple component version of the same thing as a
+ * Vector for use in, for example, the VectorTools::interpolate or the many
+ * other functions taking Function objects. It allows the user to place the
+ * desired components into an <tt>n_components</tt> long vector starting at
+ * the <tt>selected_component</tt> location in that vector and have all other
+ * components be 0.
+ *
+ * For example: Say you created a function object that returns the gravity
+ * (here, a radially inward pointing vector of magnitude 9.81):
+ * @code
+ *   const auto gravity
+ *     = [](const Point<dim> &p) -> Tensor<1,dim> { return -9.81 * (p /
+ * p.norm()); }
+ * @endcode
+ * and you want to interpolate this onto your mesh using the
+ * VectorTools::interpolate function, with a finite element for the
+ * DoFHandler object has 3 copies of a finite element with <tt>dim</tt>
+ * components, for a total of 3*dim components. To interpolate onto that
+ * DoFHandler, you need an object of type Function that has 3*dim vector
+ * components. Creating such an object from the existing <code>gravity</code>
+ * object is done using this piece of code:
+ * @code
+ * VectorFunctionFromTensorFunctionObject<dim, RangeNumberType>
+ *   gravity_vector_function(gravity, 0, 3*dim);
+ * @endcode
+ * where the <code>dim</code> components of the `gravity` function are placed
+ * into the first <code>dim</code> components of the function object.
+ *
+ * @ingroup functions
+ */
+template <int dim, typename RangeNumberType = double>
+class VectorFunctionFromTensorFunctionObject
+  : public Function<dim, RangeNumberType>
+{
+public:
+  /**
+   * Given a function object that takes a <tt>Point</tt> and returns a
+   * <tt>Tensor<1,dim, RangeNumberType></tt> value, convert this into an object
+   * that matches the Function@<dim@> interface.
+   *
+   * By default, create a Vector object of the same size as
+   * <tt>tensor_function</tt> returns, i.e., with <tt>dim</tt> components.
+   *
+   * @param tensor_function_object The TensorFunction that will form `dim` components of
+   * the resulting Vector Function object.
+   * @param n_components The total number of vector components of the
+   * resulting TensorFunction object. This clearly has to be at least `dim`.
+   * @param selected_component The first component that should be filled by
+   * the first argument.  This should be such that the entire tensor_function
+   * fits inside the <tt>n_component</tt> length return vector.
+   */
+  explicit VectorFunctionFromTensorFunctionObject(
+    const std::function<Tensor<1, dim, RangeNumberType>(const Point<dim> &)>
+                      &tensor_function_object,
+    const unsigned int selected_component = 0,
+    const unsigned int n_components       = dim);
+
+  /**
+   * This destructor is defined as virtual so as to coincide with all other
+   * aspects of class.
+   */
+  virtual ~VectorFunctionFromTensorFunctionObject() override = default;
+
+  /**
+   * Return a single component of a vector-valued function at a given point.
+   */
+  virtual RangeNumberType
+  value(const Point<dim> &p, const unsigned int component = 0) const override;
+
+  /**
+   * Return all components of a vector-valued function at a given point.
+   *
+   * <tt>values</tt> shall have the right size beforehand, i.e. #n_components.
+   */
+  virtual void
+  vector_value(const Point<dim>        &p,
+               Vector<RangeNumberType> &values) const override;
+
+  /**
+   * Return all components of a vector-valued function at a list of points.
+   *
+   * <tt>value_list</tt> shall be the same size as <tt>points</tt> and each
+   * element of the vector will be passed to vector_value() to evaluate the
+   * function
+   */
+  virtual void
+  vector_value_list(
+    const std::vector<Point<dim>>        &points,
+    std::vector<Vector<RangeNumberType>> &value_list) const override;
+
+private:
+  /**
+   * The TensorFunction object which we call when this class's vector_value()
+   * or vector_value_list() functions are called.
+   */
+  const std::function<Tensor<1, dim, RangeNumberType>(const Point<dim> &)>
+    tensor_function_object;
+
+  /**
+   * The first vector component whose value is to be filled by the given
+   * TensorFunction.  The values will be placed in components
+   * selected_component to selected_component+dim-1 for a
+   * <tt>TensorFunction<1,dim, RangeNumberType></tt> object.
+   */
+  const unsigned int selected_component;
+};
+
+
 #ifndef DOXYGEN
 // icc 2018 complains about an undefined reference
 // if we put this in the templates.h file
index d263832b0b35c13731116e8b617f5197171cbf00..371946e42fd0dc28b1e1f2daacd738223f92e904 100644 (file)
@@ -867,6 +867,8 @@ VectorFunctionFromTensorFunction<dim, RangeNumberType>::vector_value_list(
       points[p], value_list[p]);
 }
 
+
+
 template <int dim, typename RangeNumberType>
 inline Tensor<1, dim, RangeNumberType>
 VectorFunctionFromTensorFunction<dim, RangeNumberType>::gradient(
@@ -891,6 +893,8 @@ VectorFunctionFromTensorFunction<dim, RangeNumberType>::gradient(
   return tensor_gradient[component - selected_component];
 }
 
+
+
 template <int dim, typename RangeNumberType>
 void
 VectorFunctionFromTensorFunction<dim, RangeNumberType>::vector_gradient(
@@ -904,6 +908,8 @@ VectorFunctionFromTensorFunction<dim, RangeNumberType>::vector_gradient(
     gradients[i] = gradient(p, i);
 }
 
+
+
 template <int dim, typename RangeNumberType>
 void
 VectorFunctionFromTensorFunction<dim, RangeNumberType>::gradient_list(
@@ -939,6 +945,8 @@ VectorFunctionFromTensorFunction<dim, RangeNumberType>::vector_gradient_list(
     }
 }
 
+
+
 template <int dim, typename RangeNumberType>
 void
 VectorFunctionFromTensorFunction<dim, RangeNumberType>::vector_gradients(
@@ -952,6 +960,7 @@ VectorFunctionFromTensorFunction<dim, RangeNumberType>::vector_gradients(
 }
 
 
+
 template <int dim, typename RangeNumberType>
 FunctionFromFunctionObjects<dim, RangeNumberType>::FunctionFromFunctionObjects(
   const unsigned int n_components,
@@ -1058,6 +1067,101 @@ FunctionFromFunctionObjects<dim, RangeNumberType>::set_function_gradients(
   };
 }
 
+
+
+template <int dim, typename RangeNumberType>
+VectorFunctionFromTensorFunctionObject<dim, RangeNumberType>::
+  VectorFunctionFromTensorFunctionObject(
+    const std::function<Tensor<1, dim, RangeNumberType>(const Point<dim> &)>
+                      &tensor_function_object,
+    const unsigned int selected_component,
+    const unsigned int n_components)
+  : Function<dim, RangeNumberType>(n_components)
+  , tensor_function_object(tensor_function_object)
+  , selected_component(selected_component)
+{
+  // Verify that the Tensor<1,dim,RangeNumberType> will fit in the given length
+  // selected_components and not hang over the end of the vector.
+  AssertIndexRange(selected_component + dim - 1, this->n_components);
+}
+
+
+
+template <int dim, typename RangeNumberType>
+inline RangeNumberType
+VectorFunctionFromTensorFunctionObject<dim, RangeNumberType>::value(
+  const Point<dim>  &p,
+  const unsigned int component) const
+{
+  AssertIndexRange(component, this->n_components);
+
+  // if the requested component is out of the range selected, then we can
+  // return early
+  if ((component < selected_component) ||
+      (component >= selected_component + dim))
+    return 0;
+
+  // otherwise retrieve the values from the <tt>tensor_function</tt> to be
+  // placed at the <tt>selected_component</tt> to
+  // <tt>selected_component + dim - 1</tt> elements of the <tt>Vector</tt>
+  // values and pick the correct one
+  const Tensor<1, dim, RangeNumberType> tensor_value =
+    tensor_function_object(p);
+
+  return tensor_value[component - selected_component];
+}
+
+
+template <int dim, typename RangeNumberType>
+inline void
+VectorFunctionFromTensorFunctionObject<dim, RangeNumberType>::vector_value(
+  const Point<dim>        &p,
+  Vector<RangeNumberType> &values) const
+{
+  Assert(values.size() == this->n_components,
+         ExcDimensionMismatch(values.size(), this->n_components));
+
+  // Retrieve the values from the <tt>tensor_function</tt> to be placed at
+  // the <tt>selected_component</tt> to
+  // <tt>selected_component + dim - 1</tt> elements of the <tt>Vector</tt>
+  // values.
+  const Tensor<1, dim, RangeNumberType> tensor_value =
+    tensor_function_object(p);
+
+  // First we make all elements of values = 0
+  values = 0;
+
+  // Second we adjust the desired components to take on the values in
+  // <tt>tensor_value</tt>.
+  for (unsigned int i = 0; i < dim; ++i)
+    values(i + selected_component) = tensor_value[i];
+}
+
+
+/**
+ * Member function <tt>vector_value_list </tt> is the interface for giving a
+ * list of points (<code>vector<Point<dim> ></code>) of which to evaluate
+ * using the <tt>vector_value</tt> member function.  Again, this function is
+ * written so as to not replicate the function definition but passes each
+ * point on to <tt>vector_value</tt> to be evaluated.
+ */
+template <int dim, typename RangeNumberType>
+void
+VectorFunctionFromTensorFunctionObject<dim, RangeNumberType>::vector_value_list(
+  const std::vector<Point<dim>>        &points,
+  std::vector<Vector<RangeNumberType>> &value_list) const
+{
+  Assert(value_list.size() == points.size(),
+         ExcDimensionMismatch(value_list.size(), points.size()));
+
+  const unsigned int n_points = points.size();
+
+  for (unsigned int p = 0; p < n_points; ++p)
+    VectorFunctionFromTensorFunctionObject<dim, RangeNumberType>::vector_value(
+      points[p], value_list[p]);
+}
+
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif /* dealii_function_templates_h */
index 483d1d3ff504148047ddba26feddd1b42ae80c13..442476435ffc466115fc55a12f8cea6b3bbd5003 100644 (file)
@@ -24,6 +24,7 @@ for (S : REAL_SCALARS; dim : SPACE_DIMENSIONS)
     template class VectorFunctionFromScalarFunctionObject<dim, S>;
     template class FunctionFromFunctionObjects<dim, S>;
     template class VectorFunctionFromTensorFunction<dim, S>;
+    template class VectorFunctionFromTensorFunctionObject<dim, S>;
   }
 
 for (S : COMPLEX_SCALARS; dim : SPACE_DIMENSIONS)
@@ -37,4 +38,5 @@ for (S : COMPLEX_SCALARS; dim : SPACE_DIMENSIONS)
     template class VectorFunctionFromScalarFunctionObject<dim, S>;
     template class FunctionFromFunctionObjects<dim, S>;
     template class VectorFunctionFromTensorFunction<dim, S>;
+    template class VectorFunctionFromTensorFunctionObject<dim, S>;
   }

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.