]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added gradient() evaluation to VectorFunctionFromTensorFunction 15805/head
authorAbbas Ballout <abb.ballout@gmail.com>
Fri, 28 Jul 2023 09:18:04 +0000 (11:18 +0200)
committerAbbas Ballout <abb.ballout@gmail.com>
Thu, 3 Aug 2023 06:32:27 +0000 (08:32 +0200)
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 <arndtd@ornl.gov>
Update doc/news/changes/minor/20230730AbbasBallout

Co-authored-by: Daniel Arndt <arndtd@ornl.gov>
Update 20230730AbbasBallout doc minor changes

Fix test in base/functions_08_gradient.cc

doc/news/changes/minor/20230730AbbasBallout [new file with mode: 0644]
include/deal.II/base/function.h
include/deal.II/base/function.templates.h
tests/base/functions_08_gradient.cc [new file with mode: 0644]
tests/base/functions_08_gradient.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20230730AbbasBallout b/doc/news/changes/minor/20230730AbbasBallout
new file mode 100644 (file)
index 0000000..30e6639
--- /dev/null
@@ -0,0 +1,4 @@
+Improvement: Added gradient() implementation to VectorFunctionFromTensorFunction in 
+function.h file.
+<br>
+(Abbas Ballout, 2023/07/30)
index bab56d6d274093679fdbb369f3fcd9b11ccc9add..1ad351076bd94742e2c0ff4a00e1ee25e4b1eca0 100644 (file)
@@ -1148,6 +1148,60 @@ public:
     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()
index 67be4c413d05e2a91edb26918b90a51a4ae44710..3ff9abe5cdba10a2a50783c44c95489a8f294f23 100644 (file)
@@ -853,6 +853,89 @@ 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(
+  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>
diff --git a/tests/base/functions_08_gradient.cc b/tests/base/functions_08_gradient.cc
new file mode 100644 (file)
index 0000000..58ebfa0
--- /dev/null
@@ -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 <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>();
+}
diff --git a/tests/base/functions_08_gradient.output b/tests/base/functions_08_gradient.output
new file mode 100644 (file)
index 0000000..fb71de2
--- /dev/null
@@ -0,0 +1,4 @@
+
+DEAL::OK
+DEAL::OK
+DEAL::OK

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.