From eb3cf6ba7d253662a5bcdc954a1f4481abee47e4 Mon Sep 17 00:00:00 2001 From: Abbas Ballout Date: Fri, 28 Jul 2023 11:18:04 +0200 Subject: [PATCH] Added gradient() evaluation to VectorFunctionFromTensorFunction Added gradient() evaluation to VectorFunctionFromTensorFunction Added gradient() evaluation to VectorFunctionFromTensorFunction Added gradient() evaluation to VectorFunctionFromTensorFunction Added gradient() evaluation to VectorFunctionFromTensorFunction Added gradient() evaluation to VectorFunctionFromTensorFunction Add test functions_08_gradient Removed value() functions_08_gradient test Added doc/news/changes/minor/ Update doc/news/changes/minor/20230730AbbasBallout Co-authored-by: Daniel Arndt Update doc/news/changes/minor/20230730AbbasBallout Co-authored-by: Daniel Arndt Update 20230730AbbasBallout doc minor changes Fix test in base/functions_08_gradient.cc --- doc/news/changes/minor/20230730AbbasBallout | 4 + include/deal.II/base/function.h | 54 ++++++++++++ include/deal.II/base/function.templates.h | 83 ++++++++++++++++++ tests/base/functions_08_gradient.cc | 96 +++++++++++++++++++++ tests/base/functions_08_gradient.output | 4 + 5 files changed, 241 insertions(+) create mode 100644 doc/news/changes/minor/20230730AbbasBallout create mode 100644 tests/base/functions_08_gradient.cc create mode 100644 tests/base/functions_08_gradient.output diff --git a/doc/news/changes/minor/20230730AbbasBallout b/doc/news/changes/minor/20230730AbbasBallout new file mode 100644 index 0000000000..30e6639ec0 --- /dev/null +++ b/doc/news/changes/minor/20230730AbbasBallout @@ -0,0 +1,4 @@ +Improvement: Added gradient() implementation to VectorFunctionFromTensorFunction in +function.h file. +
+(Abbas Ballout, 2023/07/30) diff --git a/include/deal.II/base/function.h b/include/deal.II/base/function.h index bab56d6d27..1ad351076b 100644 --- a/include/deal.II/base/function.h +++ b/include/deal.II/base/function.h @@ -1148,6 +1148,60 @@ public: const std::vector> & points, std::vector> &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 & 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 & p, + std::vector> &gradients) const override; + + /** + * Set gradients to the gradients of the specified component of the + * function at the points. It is assumed that gradients + * already has the right size, i.e. the same size as the points + * array. + */ + virtual void + gradient_list(const std::vector> & points, + std::vector> &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> &points, + std::vector>> + &gradients) const override; + + /** + * Set gradients to the gradients of the function at the + * points, for all components. It is assumed that + * gradients already has the right size, i.e. the same size as the + * points array. + * + * The outer loop over gradients 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> &points, + std::vector>> + &gradients) const override; + private: /** * The TensorFunction object which we call when this class's vector_value() diff --git a/include/deal.II/base/function.templates.h b/include/deal.II/base/function.templates.h index 67be4c413d..3ff9abe5cd 100644 --- a/include/deal.II/base/function.templates.h +++ b/include/deal.II/base/function.templates.h @@ -853,6 +853,89 @@ VectorFunctionFromTensorFunction::vector_value_list( points[p], value_list[p]); } +template +inline Tensor<1, dim, RangeNumberType> +VectorFunctionFromTensorFunction::gradient( + const Point & 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 tensor_function to be + // // placed at the selected_component to + // // selected_component + dim - 1 elements of the Vector + // // values and pick the correct one + const Tensor<2, dim, RangeNumberType> tensor_gradient = + tensor_function.gradient(p); + + return tensor_gradient[component - selected_component]; +} + +template +void +VectorFunctionFromTensorFunction::vector_gradient( + const Point & p, + std::vector> &gradients) const +{ + AssertDimension(gradients.size(), this->n_components); + + + for (unsigned int i = 0; i < this->n_components; ++i) + gradients[i] = gradient(p, i); +} + +template +void +VectorFunctionFromTensorFunction::gradient_list( + const std::vector> & points, + std::vector> &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 +void +VectorFunctionFromTensorFunction::vector_gradient_list( + const std::vector> & points, + std::vector>> &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 +void +VectorFunctionFromTensorFunction::vector_gradients( + const std::vector> & points, + std::vector>> &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 diff --git a/tests/base/functions_08_gradient.cc b/tests/base/functions_08_gradient.cc new file mode 100644 index 0000000000..58ebfa02b3 --- /dev/null +++ b/tests/base/functions_08_gradient.cc @@ -0,0 +1,96 @@ +// --------------------------------------------------------------------- +// +// 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 +#include + +#include + +#include "../tests.h" + + +template +class X : public TensorFunction<1, dim> +{ +public: + virtual Tensor<2, dim> + gradient(const Point &p) const + { + Tensor<2, dim> T; + for (unsigned int d = 0; d < dim; ++d) + T[d] = d * p; + return T; + } +}; + + +template +void +check1() +{ + X x; + VectorFunctionFromTensorFunction 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 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> 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>(); +} diff --git a/tests/base/functions_08_gradient.output b/tests/base/functions_08_gradient.output new file mode 100644 index 0000000000..fb71de2867 --- /dev/null +++ b/tests/base/functions_08_gradient.output @@ -0,0 +1,4 @@ + +DEAL::OK +DEAL::OK +DEAL::OK -- 2.39.5