--- /dev/null
+Improvement: Added gradient() implementation to VectorFunctionFromTensorFunction in
+function.h file.
+<br>
+(Abbas Ballout, 2023/07/30)
const std::vector<Point<dim>> & points,
std::vector<Vector<RangeNumberType>> &value_list) const override;
+ /**
+ * Return the gradient of the specified component of the function at the given
+ * point.
+ */
+ virtual Tensor<1, dim, RangeNumberType>
+ gradient(const Point<dim> & p,
+ const unsigned int component = 0) const override;
+
+ /**
+ * Return the gradient of all components of the function at the given point.
+ */
+ virtual void
+ vector_gradient(
+ const Point<dim> & p,
+ std::vector<Tensor<1, dim, RangeNumberType>> &gradients) const override;
+
+ /**
+ * Set <tt>gradients</tt> to the gradients of the specified component of the
+ * function at the <tt>points</tt>. It is assumed that <tt>gradients</tt>
+ * already has the right size, i.e. the same size as the <tt>points</tt>
+ * array.
+ */
+ virtual void
+ gradient_list(const std::vector<Point<dim>> & points,
+ std::vector<Tensor<1, dim, RangeNumberType>> &gradients,
+ const unsigned int component = 0) const override;
+
+ /**
+ * For each component of the function, fill a vector of gradient values, one
+ * for each point.
+ *
+ * The default implementation of this function in Function calls
+ * value_list() for each component. In order to improve performance, this
+ * can be reimplemented in derived classes to speed up performance.
+ */
+ virtual void
+ vector_gradients(const std::vector<Point<dim>> &points,
+ std::vector<std::vector<Tensor<1, dim, RangeNumberType>>>
+ &gradients) const override;
+
+ /**
+ * Set <tt>gradients</tt> to the gradients of the function at the
+ * <tt>points</tt>, for all components. It is assumed that
+ * <tt>gradients</tt> already has the right size, i.e. the same size as the
+ * <tt>points</tt> array.
+ *
+ * The outer loop over <tt>gradients</tt> is over the points in the list,
+ * the inner loop over the different components of the function.
+ */
+ virtual void
+ vector_gradient_list(const std::vector<Point<dim>> &points,
+ std::vector<std::vector<Tensor<1, dim, RangeNumberType>>>
+ &gradients) const override;
+
private:
/**
* The TensorFunction object which we call when this class's vector_value()
points[p], value_list[p]);
}
+template <int dim, typename RangeNumberType>
+inline Tensor<1, dim, RangeNumberType>
+VectorFunctionFromTensorFunction<dim, RangeNumberType>::gradient(
+ 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 Tensor<1, dim>();
+
+ // // 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<2, dim, RangeNumberType> tensor_gradient =
+ tensor_function.gradient(p);
+
+ return tensor_gradient[component - selected_component];
+}
+
+template <int dim, typename RangeNumberType>
+void
+VectorFunctionFromTensorFunction<dim, RangeNumberType>::vector_gradient(
+ const Point<dim> & p,
+ std::vector<Tensor<1, dim, RangeNumberType>> &gradients) const
+{
+ AssertDimension(gradients.size(), this->n_components);
+
+
+ for (unsigned int i = 0; i < this->n_components; ++i)
+ gradients[i] = gradient(p, i);
+}
+
+template <int dim, typename RangeNumberType>
+void
+VectorFunctionFromTensorFunction<dim, RangeNumberType>::gradient_list(
+ const std::vector<Point<dim>> & points,
+ std::vector<Tensor<1, dim, RangeNumberType>> &gradients,
+ const unsigned int component) const
+{
+ Assert(gradients.size() == points.size(),
+ ExcDimensionMismatch(gradients.size(), points.size()));
+
+
+ for (unsigned int i = 0; i < points.size(); ++i)
+ gradients[i] = gradient(points[i], component);
+}
+
+
+
+template <int dim, typename RangeNumberType>
+void
+VectorFunctionFromTensorFunction<dim, RangeNumberType>::vector_gradient_list(
+ const std::vector<Point<dim>> & points,
+ std::vector<std::vector<Tensor<1, dim, RangeNumberType>>> &gradients) const
+{
+ Assert(gradients.size() == points.size(),
+ ExcDimensionMismatch(gradients.size(), points.size()));
+
+
+ for (unsigned int i = 0; i < points.size(); ++i)
+ {
+ Assert(gradients[i].size() == this->n_components,
+ ExcDimensionMismatch(gradients[i].size(), this->n_components));
+ vector_gradient(points[i], gradients[i]);
+ }
+}
+
+template <int dim, typename RangeNumberType>
+void
+VectorFunctionFromTensorFunction<dim, RangeNumberType>::vector_gradients(
+ const std::vector<Point<dim>> & points,
+ std::vector<std::vector<Tensor<1, dim, RangeNumberType>>> &gradients) const
+{
+ const unsigned int n = this->n_components;
+ AssertDimension(gradients.size(), n);
+ for (unsigned int i = 0; i < n; ++i)
+ gradient_list(points, gradients[i], i);
+}
template <int dim, typename RangeNumberType>
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2007 - 2018 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Test VectorFunctionTensorFunction gradient
+
+#include <deal.II/base/function.h>
+#include <deal.II/base/tensor_function.h>
+
+#include <deal.II/lac/vector.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+class X : public TensorFunction<1, dim>
+{
+public:
+ virtual Tensor<2, dim>
+ gradient(const Point<dim> &p) const
+ {
+ Tensor<2, dim> T;
+ for (unsigned int d = 0; d < dim; ++d)
+ T[d] = d * p;
+ return T;
+ }
+};
+
+
+template <int dim>
+void
+check1()
+{
+ X<dim> x;
+ VectorFunctionFromTensorFunction<dim> object(x, 1, dim + 2);
+
+ AssertThrow(object.n_components == dim + 2, ExcInternalError());
+
+ Tensor<1, dim> T;
+
+ for (unsigned int i = 0; i < 10; ++i)
+ {
+ Point<dim> p;
+ for (unsigned int d = 0; d < dim; ++d)
+ p[d] = i + d;
+
+ for (unsigned int c = 0; c < dim + 2; ++c)
+ if (c == 0 || c == dim + 1)
+ {
+ AssertThrow(object.gradient(p, c) == T, ExcInternalError());
+ }
+ else
+ {
+ AssertThrow(object.gradient(p, c) == (c - 1) * p,
+ ExcInternalError());
+ }
+
+ std::vector<Tensor<1, dim>> v(dim + 2);
+ object.vector_gradient(p, v);
+ for (unsigned int c = 0; c < dim + 2; ++c)
+ if (c == 0 || c == dim + 1)
+ {
+ AssertThrow(v[c] == T, ExcInternalError());
+ }
+ else
+ {
+ AssertThrow(v[c] == (c - 1) * p, ExcInternalError());
+ }
+ }
+
+ deallog << "OK" << std::endl;
+}
+
+
+
+int
+main()
+{
+ initlog();
+
+ check1<1>();
+ check1<2>();
+ check1<3>();
+}
--- /dev/null
+
+DEAL::OK
+DEAL::OK
+DEAL::OK
projects / dealii.git / commitdiff