]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add python interface for FunctionManifold 9179/head
authorAlexander Grayver <agrayver@gmail.com>
Tue, 17 Dec 2019 15:27:37 +0000 (16:27 +0100)
committerAlexander Grayver <agrayver@erdw.ethz.ch>
Thu, 9 Jan 2020 08:40:01 +0000 (09:40 +0100)
contrib/python-bindings/include/function_wrapper.h [new file with mode: 0644]
contrib/python-bindings/include/manifold_wrapper.h
contrib/python-bindings/source/export_manifold.cc
contrib/python-bindings/source/manifold_wrapper.cc
contrib/python-bindings/tests/manifold_wrapper.py
contrib/python-bindings/tests/manifold_wrapper.vtk [new file with mode: 0644]
doc/news/changes/minor/20191217Grayver-2 [new file with mode: 0644]
include/deal.II/grid/manifold_lib.h
source/grid/manifold_lib.cc

diff --git a/contrib/python-bindings/include/function_wrapper.h b/contrib/python-bindings/include/function_wrapper.h
new file mode 100644 (file)
index 0000000..5fb1a1d
--- /dev/null
@@ -0,0 +1,70 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 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.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_function_wrapper_h
+#define dealii_function_wrapper_h
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/base/function.h>
+
+#include <boost/python.hpp>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace python
+{
+  template <int dim>
+  class FunctionWrapper : public Function<dim, double>
+  {
+  public:
+    FunctionWrapper(boost::python::object &python_function,
+                    unsigned               n_components)
+      : Function<dim, double>(n_components)
+      , python_function(python_function)
+    {}
+
+    FunctionWrapper(const FunctionWrapper &other)
+    {
+      python_function = other.python_function;
+    }
+
+    virtual void
+    vector_value(const Point<dim> &p, Vector<double> &values) const override
+    {
+      boost::python::list p_list_in;
+
+      for (int d = 0; d < dim; ++d)
+        p_list_in.append(p[d]);
+
+      boost::python::list p_list_out =
+        boost::python::extract<boost::python::list>(python_function(p_list_in));
+
+      for (size_t i = 0; i < this->n_components; ++i)
+        values[i] = boost::python::extract<double>(p_list_out[i]);
+    }
+
+  private:
+    /**
+     * A callback to a python function
+     */
+    boost::python::object python_function;
+  };
+
+} // namespace python
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index c1e6cf21136e15622edf22e6ad795748a13ddf14..2e9f0c8eb12dd95e181eaff56aeaf04f88fc4430 100644 (file)
 
 #include <deal.II/base/config.h>
 
+#include <deal.II/base/function.h>
+
+#include <deal.II/lac/vector.h>
+
 #include <boost/python.hpp>
 
+#include <function_wrapper.h>
 #include <point_wrapper.h>
 
 DEAL_II_NAMESPACE_OPEN
@@ -45,23 +50,39 @@ namespace python
     ~ManifoldWrapper();
 
     /**
-     * Create SphericalManifold
+     * Create SphericalManifold.
      */
     void
     create_spherical(const PointWrapper center);
 
     /**
-     * Create PolarManifold
+     * Create PolarManifold.
      */
     void
     create_polar(const PointWrapper center);
 
     /**
-     * Create CylindricalManifold
+     * Create CylindricalManifold.
      */
     void
     create_cylindrical(const int axis = 0, const double tolerance = 1e-10);
 
+    /**
+     * Create FunctionManifold with string expressions for the push
+     * forward and pull back functions.
+     */
+    void
+    create_function_string(const std::string &push_forward,
+                           const std::string &pull_back);
+
+    /**
+     * Create FunctionManifold with python the push forward and
+     * pull back functions.
+     */
+    void
+    create_function(boost::python::object &push_forward,
+                    boost::python::object &pull_back);
+
     /**
      * Return pointer to an underlying manifold object
      */
index d9225359b32ec1f7befee9a4ed4a13cb20c77957..87a3fe533dcce9dd707c06f1e172d8fdfd56d440 100644 (file)
@@ -28,16 +28,26 @@ namespace python
 
 
   const char create_spherical_docstring[] =
-    " Create spherical manifold with a given center point       \n";
+    " Create spherical manifold with a given center point.      \n";
 
 
   const char create_polar_docstring[] =
-    " Create polar manifold with a given center point           \n";
+    " Create polar manifold with a given center point.          \n";
+
+
+  const char create_function_docstring[] =
+    " Create manifold with the given python push forward and    \n"
+    " pull back functions.                                      \n";
+
+
+  const char create_function_string_docstring[] =
+    " Create manifold with given string expression for the push \n"
+    " forward and pull back functions.                          \n";
 
 
   const char create_cylindrical_docstring[] =
     " Create cylindrical manifold along a given axis            \n"
-    " (0 - x, 1 - y, 2 - z)                                     \n";
+    " (0 - x, 1 - y, 2 - z).                                    \n";
 
 
 
@@ -60,7 +70,15 @@ namespace python
            &ManifoldWrapper::create_cylindrical,
            create_cylindrical_overloads(
              boost::python::args("self", "axis", "tolerance"),
-             create_cylindrical_docstring));
+             create_cylindrical_docstring))
+      .def("create_function",
+           &ManifoldWrapper::create_function,
+           create_function_docstring,
+           boost::python::args("self", "push_forward", "pull_back"))
+      .def("create_function_string",
+           &ManifoldWrapper::create_function_string,
+           create_function_string_docstring,
+           boost::python::args("self", "push_forward", "pull_back"));
   }
 } // namespace python
 
index 14d4888967248f07da2ef605750d47da04981ef1..d364dc014a5279696cfa8adcc4e02159aaaea022 100644 (file)
@@ -13,6 +13,8 @@
 //
 // ---------------------------------------------------------------------
 
+#include <deal.II/base/std_cxx14/memory.h>
+
 #include <deal.II/grid/manifold_lib.h>
 
 #include <manifold_wrapper.h>
@@ -55,6 +57,28 @@ namespace python
 
 
 
+    template <int dim, int spacedim>
+    Manifold<dim, spacedim> *
+    create_function_manifold(const std::string &push_forward,
+                             const std::string &pull_back)
+    {
+      return new FunctionManifold<dim, spacedim>(push_forward, pull_back);
+    }
+
+
+
+    template <int dim, int spacedim>
+    Manifold<dim, spacedim> *
+    create_function_manifold(boost::python::object &push_forward,
+                             boost::python::object &pull_back)
+    {
+      return new FunctionManifold<dim, spacedim>(
+        std_cxx14::make_unique<FunctionWrapper<dim>>(push_forward, spacedim),
+        std_cxx14::make_unique<FunctionWrapper<spacedim>>(pull_back, dim));
+    }
+
+
+
     template <int dim, int spacedim>
     Manifold<dim, spacedim> *
     clone(void *manifold_ptr)
@@ -78,6 +102,12 @@ namespace python
         {
           return new CylindricalManifold<dim, spacedim>(*d);
         }
+      else if (const FunctionManifold<dim, spacedim> *d =
+                 dynamic_cast<const FunctionManifold<dim, spacedim> *>(
+                   manifold))
+        {
+          return new FunctionManifold<dim, spacedim>(*d);
+        }
       else
         ExcMessage("Unsupported manifold type in clone.");
 
@@ -120,7 +150,9 @@ namespace python
         manifold_ptr = internal::clone<3, 3>(other.manifold_ptr);
       }
     else
-      AssertThrow(false, ExcMessage("Wrong dim-spacedim combination."));
+      AssertThrow(false,
+                  ExcMessage(
+                    "Given dim-spacedim combination is not implemented."));
   }
 
 
@@ -171,7 +203,9 @@ namespace python
         manifold_ptr = internal::create_spherical_manifold<3, 3>(center);
       }
     else
-      AssertThrow(false, ExcMessage("Wrong dim-spacedim combination."));
+      AssertThrow(false,
+                  ExcMessage(
+                    "Given dim-spacedim combination is not implemented."));
   }
 
 
@@ -192,7 +226,9 @@ namespace python
         manifold_ptr = internal::create_polar_manifold<3, 3>(center);
       }
     else
-      AssertThrow(false, ExcMessage("Wrong dim-spacedim combination."));
+      AssertThrow(false,
+                  ExcMessage(
+                    "Given dim-spacedim combination is not implemented."));
   }
 
 
@@ -211,7 +247,63 @@ namespace python
           internal::create_cylindrical_manifold<3, 3>(axis, tolerance);
       }
     else
-      AssertThrow(false, ExcMessage("Wrong dim-spacedim combination."));
+      AssertThrow(false,
+                  ExcMessage(
+                    "Given dim-spacedim combination is not implemented."));
+  }
+
+
+
+  void
+  ManifoldWrapper::create_function_string(const std::string &push_forward,
+                                          const std::string &pull_back)
+  {
+    if ((dim == 2) && (spacedim == 2))
+      {
+        manifold_ptr =
+          internal::create_function_manifold<2, 2>(push_forward, pull_back);
+      }
+    else if ((dim == 2) && (spacedim == 3))
+      {
+        manifold_ptr =
+          internal::create_function_manifold<2, 3>(push_forward, pull_back);
+      }
+    else if ((dim == 3) && (spacedim == 3))
+      {
+        manifold_ptr =
+          internal::create_function_manifold<3, 3>(push_forward, pull_back);
+      }
+    else
+      AssertThrow(false,
+                  ExcMessage(
+                    "Given dim-spacedim combination is not implemented."));
+  }
+
+
+
+  void
+  ManifoldWrapper::create_function(boost::python::object &push_forward,
+                                   boost::python::object &pull_back)
+  {
+    if ((dim == 2) && (spacedim == 2))
+      {
+        manifold_ptr =
+          internal::create_function_manifold<2, 2>(push_forward, pull_back);
+      }
+    else if ((dim == 2) && (spacedim == 3))
+      {
+        manifold_ptr =
+          internal::create_function_manifold<2, 3>(push_forward, pull_back);
+      }
+    else if ((dim == 3) && (spacedim == 3))
+      {
+        manifold_ptr =
+          internal::create_function_manifold<3, 3>(push_forward, pull_back);
+      }
+    else
+      AssertThrow(false,
+                  ExcMessage(
+                    "Given dim-spacedim combination is not implemented."));
   }
 
 
index d3270ccca25c9b27668679dfbe255345280f06fc..f4330bd2adb7316c7e54e70d274ffa8fcc682d33 100644 (file)
@@ -67,5 +67,46 @@ class TestManifoldWrapperBall(unittest.TestCase):
 
         self.assertTrue(abs(volume - 4./3. * math.pi) / (4./3.*math.pi) < 2e-2)
 
+class TestManifoldWrapperFunction(unittest.TestCase):
+
+    def setUp(self):
+        self.manifold_1 = Manifold(dim = 2, spacedim = 2)
+        self.manifold_1.create_function("x^2;y^2", "sqrt(x);sqrt(y)")
+
+        self.manifold_2 = Manifold(dim = 2, spacedim = 2)
+        self.manifold_2.create_function(lambda p: [p[0]**2., p[1]**2.],\
+                                        lambda p: [math.sqrt(p[0]), math.sqrt(p[1])] )
+
+        self.tria_reference = Triangulation('2D')
+        self.tria_reference.read('manifold_wrapper.vtk', 'vtk')
+
+    def test_manifold_str(self):
+        self.triangulation = Triangulation('2D')
+        p_center = Point([0., 0., 0.])
+        self.triangulation.generate_hyper_cube()
+        self.triangulation.reset_manifold(number = 0)
+        self.triangulation.set_manifold(number = 0, manifold = self.manifold_2)
+        for cell in self.triangulation.active_cells():
+            cell.set_all_manifold_ids(0)
+
+        self.triangulation.refine_global(2)
+
+        for cell_ref, cell in zip(self.tria_reference.active_cells(), self.triangulation.active_cells()):
+            self.assertTrue(abs(cell_ref.measure() - cell.measure()) < 1e-8)
+
+    def test_manifold_lambda(self):
+        self.triangulation = Triangulation('2D')
+        p_center = Point([0., 0., 0.])
+        self.triangulation.generate_hyper_cube()
+        self.triangulation.reset_manifold(number = 0)
+        self.triangulation.set_manifold(number = 0, manifold = self.manifold_2)
+        for cell in self.triangulation.active_cells():
+            cell.set_all_manifold_ids(0)
+
+        self.triangulation.refine_global(2)
+
+        for cell_ref, cell in zip(self.tria_reference.active_cells(), self.triangulation.active_cells()):
+            self.assertTrue(abs(cell_ref.measure() - cell.measure()) < 1e-8)
+
 if __name__ == '__main__':
     unittest.main()
diff --git a/contrib/python-bindings/tests/manifold_wrapper.vtk b/contrib/python-bindings/tests/manifold_wrapper.vtk
new file mode 100644 (file)
index 0000000..7fd1552
--- /dev/null
@@ -0,0 +1,106 @@
+# vtk DataFile Version 3.0
+Triangulation generated with deal.II
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 25 double
+0 0 0
+1 0 0
+0 1 0
+1 1 0
+0.25 0 0
+0 0.25 0
+1 0.25 0
+0.25 1 0
+0.25 0.25 0
+0.0625 0 0
+0.5625 0 0
+0 0.0625 0
+0 0.5625 0
+1 0.0625 0
+1 0.5625 0
+0.0625 1 0
+0.5625 1 0
+0.25 0.0625 0
+0.25 0.5625 0
+0.0625 0.25 0
+0.5625 0.25 0
+0.0625 0.0625 0
+0.5625 0.0625 0
+0.0625 0.5625 0
+0.5625 0.5625 0
+
+CELLS 56 200
+4 0 9 21 11
+4 9 4 17 21
+4 11 21 19 5
+4 21 17 8 19
+4 4 10 22 17
+4 10 1 13 22
+4 17 22 20 8
+4 22 13 6 20
+4 5 19 23 12
+4 19 8 18 23
+4 12 23 15 2
+4 23 18 7 15
+4 8 20 24 18
+4 20 6 14 24
+4 18 24 16 7
+4 24 14 3 16
+2 0 9
+2 9 4
+2 4 10
+2 10 1
+2 0 11
+2 11 5
+2 5 12
+2 12 2
+2 1 13
+2 13 6
+2 6 14
+2 14 3
+2 2 15
+2 15 7
+2 7 16
+2 16 3
+2 4 17
+2 17 8
+2 8 18
+2 18 7
+2 5 19
+2 19 8
+2 8 20
+2 20 6
+2 9 21
+2 21 19
+2 11 21
+2 21 17
+2 10 22
+2 22 20
+2 17 22
+2 22 13
+2 19 23
+2 23 15
+2 12 23
+2 23 18
+2 20 24
+2 24 16
+2 18 24
+2 24 14
+
+CELL_TYPES 56
+9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 
+3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 
+
+
+CELL_DATA 56
+SCALARS MaterialID int 1
+LOOKUP_TABLE default
+0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 
+0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 
+
+
+SCALARS ManifoldID int 1
+LOOKUP_TABLE default
+0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 
+0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 
+
diff --git a/doc/news/changes/minor/20191217Grayver-2 b/doc/news/changes/minor/20191217Grayver-2
new file mode 100644 (file)
index 0000000..411bfd8
--- /dev/null
@@ -0,0 +1,3 @@
+Added python wrappers for the FunctionManifold class
+<br>
+(Alexander Grayver, 2019/12/17)
index 7094b3a113f2347ebc1029515045b434a8849cc3..0a4a3ec3077a970275e7dcea02752a5d2a0029b2 100644 (file)
@@ -600,6 +600,23 @@ public:
     const Tensor<1, chartdim> &periodicity = Tensor<1, chartdim>(),
     const double               tolerance   = 1e-10);
 
+  /**
+   * Same as previous, except this constructor takes ownership of the Function
+   * objects passed as first and second argument, and is ultimately in charge
+   * of deleting the pointers when the FunctionManifold object is destroyed.
+   *
+   * This constructor is useful because it allows creating function objects at
+   * the place of calling the constructor without having to name and later
+   * delete these objects. This allows the following idiom:
+   * FunctionManifold<dim> manifold(std_cxx14::make_unique<MyPushForward>(...),
+   *                                std_cxx14::make_unique<MyPullBack>(...));
+   */
+  FunctionManifold(
+    std::unique_ptr<Function<chartdim>> push_forward_function,
+    std::unique_ptr<Function<spacedim>> pull_back_function,
+    const Tensor<1, chartdim> &         periodicity = Tensor<1, chartdim>(),
+    const double                        tolerance   = 1e-10);
+
   /**
    * Expressions constructor. Takes the expressions of the push_forward
    * function of spacedim components, and of the pull_back function of @p
@@ -711,7 +728,7 @@ private:
    * destructor will delete the function objects pointed to be the two
    * pointers.
    */
-  const bool owns_pointers;
+  bool owns_pointers;
 
   /**
    * The expression used to construct the push_forward function.
index d9723307e9b5fa8e3d186ae42068505a85baaf6d..135c2ca15df61d234fb9a765504a2c0c9cd1a5ec 100644 (file)
@@ -1321,6 +1321,26 @@ FunctionManifold<dim, spacedim, chartdim>::FunctionManifold(
 
 
 
+template <int dim, int spacedim, int chartdim>
+FunctionManifold<dim, spacedim, chartdim>::FunctionManifold(
+  std::unique_ptr<Function<chartdim>> push_forward_function,
+  std::unique_ptr<Function<spacedim>> pull_back_function,
+  const Tensor<1, chartdim> &         periodicity,
+  const double                        tolerance)
+  : ChartManifold<dim, spacedim, chartdim>(periodicity)
+  , const_map()
+  , push_forward_function(push_forward_function.release())
+  , pull_back_function(pull_back_function.release())
+  , tolerance(tolerance)
+  , owns_pointers(true)
+  , finite_difference_step(0)
+{
+  AssertDimension(push_forward_function->n_components, spacedim);
+  AssertDimension(pull_back_function->n_components, chartdim);
+}
+
+
+
 template <int dim, int spacedim, int chartdim>
 FunctionManifold<dim, spacedim, chartdim>::FunctionManifold(
   const std::string                                 push_forward_expression,
@@ -1376,13 +1396,15 @@ FunctionManifold<dim, spacedim, chartdim>::clone() const
   // push forward and the pull back charts, or by providing two Function
   // objects. In the first case, the push_forward and pull_back functions are
   // created internally in FunctionManifold, and destroyed when this object is
-  // deleted. We need to make sure that our cloned object is constructed in the
+  // deleted. In the second case, the function objects are destroyed if they
+  // are passed as pointers upon construction.
+  // We need to make sure that our cloned object is constructed in the
   // same way this class was constructed, and that its internal Function
   // pointers point either to the same Function objects used to construct this
-  // function (owns_pointers == false) or that the newly generated manifold
-  // creates internally the push_forward and pull_back functions using the same
-  // expressions that were used to construct this class (own_pointers == true).
-  if (owns_pointers == true)
+  // function or that the newly generated manifold creates internally the
+  // push_forward and pull_back functions using the same expressions that were
+  // used to construct this class.
+  if (!(push_forward_expression.empty() && pull_back_expression.empty()))
     {
       return std_cxx14::make_unique<FunctionManifold<dim, spacedim, chartdim>>(
         push_forward_expression,
@@ -1395,11 +1417,18 @@ FunctionManifold<dim, spacedim, chartdim>::clone() const
         finite_difference_step);
     }
   else
-    return std_cxx14::make_unique<FunctionManifold<dim, spacedim, chartdim>>(
-      *push_forward_function,
-      *pull_back_function,
-      this->get_periodicity(),
-      tolerance);
+    {
+      auto manifold_clone =
+        std_cxx14::make_unique<FunctionManifold<dim, spacedim, chartdim>>(
+          *push_forward_function,
+          *pull_back_function,
+          this->get_periodicity(),
+          tolerance);
+
+      Assert(owns_pointers == false, ExcNotImplemented());
+
+      return manifold_clone;
+    }
 }
 
 

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.