]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add manifolds to the python bindings
authorAlexander Grayver <agrayver@erdw.ethz.ch>
Thu, 28 Nov 2019 16:50:10 +0000 (17:50 +0100)
committerAlexander Grayver <agrayver@erdw.ethz.ch>
Fri, 29 Nov 2019 11:27:10 +0000 (12:27 +0100)
contrib/python-bindings/include/manifold_wrapper.h [new file with mode: 0644]
contrib/python-bindings/include/triangulation_wrapper.h
contrib/python-bindings/source/CMakeLists.txt
contrib/python-bindings/source/export_manifold.cc [new file with mode: 0644]
contrib/python-bindings/source/export_triangulation.cc
contrib/python-bindings/source/manifold_wrapper.cc [new file with mode: 0644]
contrib/python-bindings/source/triangulation_wrapper.cc
contrib/python-bindings/source/wrappers.cc
doc/news/changes/minor/20191128Grayver [new file with mode: 0644]

diff --git a/contrib/python-bindings/include/manifold_wrapper.h b/contrib/python-bindings/include/manifold_wrapper.h
new file mode 100644 (file)
index 0000000..e72018c
--- /dev/null
@@ -0,0 +1,93 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 - 2017 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_manifold_wrapper_h
+#define dealii_manifold_wrapper_h
+
+#include <deal.II/base/config.h>
+
+#include <boost/python.hpp>
+
+#include <point_wrapper.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace python
+{
+  class ManifoldWrapper
+  {
+  public:
+    /**
+     * Copy constructor.
+     */
+    ManifoldWrapper(const ManifoldWrapper &other);
+
+    /**
+     * Constructor.
+     */
+    ManifoldWrapper(const int dim, const int spacedim);
+
+    /**
+     * Destructor.
+     */
+    ~ManifoldWrapper();
+
+    /**
+     * Create SphericalManifold
+     */
+    void
+    create_spherical(const PointWrapper center);
+
+    /**
+     * Create PolarManifold
+     */
+    void
+    create_polar(const PointWrapper center);
+
+    /**
+     * Create CylindricalManifold
+     */
+    void
+    create_cylindrical(const int axis, const double tolerance);
+
+    /**
+     * Return pointer to an underlying manifold object
+     */
+    void *
+    get_manifold() const;
+
+
+  private:
+    /**
+     * Dimension of the underlying Manifold object.
+     */
+    int dim;
+
+    /**
+     * Space dimension of the underlying Manifold object.
+     */
+    int spacedim;
+
+    /**
+     * Pointer that can be casted to the underlying Mapping object.
+     */
+    void *manifold_ptr;
+  };
+
+} // namespace python
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index 5e3f679be95d4bf1cdc10b99d4922384c646e950..bbab008dee95cbdc3c0901440b0a1e4b17879367 100644 (file)
@@ -20,6 +20,7 @@
 
 #include <boost/python.hpp>
 
+#include <manifold_wrapper.h>
 #include <point_wrapper.h>
 
 #include <string>
@@ -301,6 +302,12 @@ namespace python
     void
     flatten_triangulation(TriangulationWrapper &tria_out);
 
+    /**
+     *
+     */
+    void
+    set_manifold(const int number, ManifoldWrapper &manifold);
+
     /**
      * Refine all the cells @p n times.
      */
index 376dbe35b0ab23268e870fb8f76bacdc26a13887..775a2d7e52269ba5ad85ac7c326fcbfc997f30e2 100644 (file)
@@ -95,11 +95,13 @@ SET(_src
   export_point.cc
   export_triangulation.cc
   export_mapping.cc
+  export_manifold.cc
   cell_accessor_wrapper.cc
   tria_accessor_wrapper.cc
   point_wrapper.cc
   triangulation_wrapper.cc
   mapping_wrapper.cc
+  manifold_wrapper.cc
   )
 
 FOREACH(_build ${DEAL_II_BUILD_TYPES})
diff --git a/contrib/python-bindings/source/export_manifold.cc b/contrib/python-bindings/source/export_manifold.cc
new file mode 100644 (file)
index 0000000..e752153
--- /dev/null
@@ -0,0 +1,60 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 - 2017 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.
+//
+// ---------------------------------------------------------------------
+
+#include <boost/python.hpp>
+
+#include <manifold_wrapper.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace python
+{
+  const char create_spherical_docstring[] =
+    " Create spherical manifold with a given center point       \n";
+
+
+  const char create_polar_docstring[] =
+    " Create polar manifold with a given center point           \n";
+
+
+  const char create_cylindrical_docstring[] =
+    " Create cylindrical manifold along a given axis            \n"
+    " (0 - x, 1 - y, 2 - z)                                     \n";
+
+
+
+  void
+  export_manifold()
+  {
+    boost::python::class_<ManifoldWrapper>(
+      "ManifoldWrapper",
+      boost::python::init<const int, const int>(
+        boost::python::args("dim", "spacedim")))
+      .def("create_spherical",
+           &ManifoldWrapper::create_spherical,
+           create_spherical_docstring,
+           boost::python::args("self", "center"))
+      .def("create_polar",
+           &ManifoldWrapper::create_polar,
+           create_polar_docstring,
+           boost::python::args("self", "center"))
+      .def("create_cylindrical",
+           &ManifoldWrapper::create_cylindrical,
+           create_cylindrical_docstring,
+           boost::python::args("self", "axis", "tolerance"));
+  }
+} // namespace python
+
+DEAL_II_NAMESPACE_CLOSE
index ed38b6f58f9917e983c11c3cc22179ca8d1f8c7a..f1bab25860306aa240f8f3a7beb39d88feec8225 100644 (file)
@@ -333,6 +333,13 @@ namespace python
 
 
 
+  const char set_manifold_docstring[] =
+    "Assign a manifold object to a certain part of the triangulation.       \n"
+    "The manifold_object is not copied and MUST persist until the           \n"
+    "triangulation is destroyed.                                            \n";
+
+
+
   void
   export_triangulation()
   {
@@ -470,7 +477,11 @@ namespace python
       .def("load",
            &TriangulationWrapper::load,
            load_docstring,
-           boost::python::args("self", "filename"));
+           boost::python::args("self", "filename"))
+      .def("set_manifold",
+           &TriangulationWrapper::set_manifold,
+           set_manifold_docstring,
+           boost::python::args("self", "number", "manifold"));
   }
 } // namespace python
 
diff --git a/contrib/python-bindings/source/manifold_wrapper.cc b/contrib/python-bindings/source/manifold_wrapper.cc
new file mode 100644 (file)
index 0000000..9c78c60
--- /dev/null
@@ -0,0 +1,220 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 - 2017 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.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/grid/manifold_lib.h>
+
+#include <manifold_wrapper.h>
+#include <point_wrapper.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace python
+{
+  namespace internal
+  {
+    template <int dim, int spacedim>
+    Manifold<dim, spacedim> *
+    create_spherical_manifold(PointWrapper center)
+    {
+      const Point<spacedim> *point =
+        static_cast<const Point<spacedim> *>(center.get_point());
+      return new SphericalManifold<dim, spacedim>(*point);
+    }
+
+
+
+    template <int dim, int spacedim>
+    Manifold<dim, spacedim> *
+    create_polar_manifold(PointWrapper center)
+    {
+      const Point<spacedim> *point =
+        static_cast<const Point<spacedim> *>(center.get_point());
+      return new PolarManifold<dim, spacedim>(*point);
+    }
+
+
+
+    template <int dim, int spacedim>
+    Manifold<dim, spacedim> *
+    create_cylindrical_manifold(const int axis, const double tolerance)
+    {
+      return new CylindricalManifold<dim, spacedim>(axis, tolerance);
+    }
+
+
+
+    template <int dim, int spacedim>
+    Manifold<dim, spacedim> *
+    clone(void *manifold_ptr)
+    {
+      const Manifold<dim, spacedim> *manifold =
+        static_cast<const Manifold<dim, spacedim> *>(manifold_ptr);
+
+      if (const SphericalManifold<dim, spacedim> *d =
+            dynamic_cast<const SphericalManifold<dim, spacedim> *>(manifold))
+        {
+          return new SphericalManifold<dim, spacedim>(*d);
+        }
+      else if (const PolarManifold<dim, spacedim> *d =
+                 dynamic_cast<const PolarManifold<dim, spacedim> *>(manifold))
+        {
+          return new PolarManifold<dim, spacedim>(*d);
+        }
+      else
+        ExcMessage("Unsupported manifold type in clone.");
+
+      return nullptr;
+    }
+  } // namespace internal
+
+
+
+  ManifoldWrapper::ManifoldWrapper(const int dim, const int spacedim)
+    : dim(dim)
+    , spacedim(spacedim)
+  {
+    AssertThrow(((dim == 2) && (spacedim == 2)) ||
+                  ((dim == 2) && (spacedim == 3)) ||
+                  ((dim == 3) && (spacedim == 3)),
+                ExcMessage("Wrong dim-spacedim combination."));
+  }
+
+
+
+  ManifoldWrapper::ManifoldWrapper(const ManifoldWrapper &other)
+  {
+    dim      = other.dim;
+    spacedim = other.spacedim;
+
+    AssertThrow(other.manifold_ptr != nullptr,
+                ExcMessage("Underlying manifold does not exist."));
+
+    if ((dim == 2) && (spacedim == 2))
+      {
+        manifold_ptr = internal::clone<2, 2>(other.manifold_ptr);
+      }
+    else if ((dim == 2) && (spacedim == 3))
+      {
+        manifold_ptr = internal::clone<2, 3>(other.manifold_ptr);
+      }
+    else if ((dim == 3) && (spacedim == 3))
+      {
+        manifold_ptr = internal::clone<3, 3>(other.manifold_ptr);
+      }
+    else
+      AssertThrow(false, ExcMessage("Wrong dim-spacedim combination."));
+  }
+
+
+
+  ManifoldWrapper::~ManifoldWrapper()
+  {
+    if (dim != -1)
+      {
+        if ((dim == 2) && (spacedim == 2))
+          {
+            // We cannot call delete on a void pointer so cast the void pointer
+            // back first.
+            Manifold<2, 2> *tmp = static_cast<Manifold<2, 2> *>(manifold_ptr);
+            delete tmp;
+          }
+        else if ((dim == 2) && (spacedim == 3))
+          {
+            Manifold<2, 3> *tmp = static_cast<Manifold<2, 3> *>(manifold_ptr);
+            delete tmp;
+          }
+        else
+          {
+            Manifold<3, 3> *tmp = static_cast<Manifold<3, 3> *>(manifold_ptr);
+            delete tmp;
+          }
+
+        dim          = -1;
+        spacedim     = -1;
+        manifold_ptr = nullptr;
+      }
+  }
+
+
+
+  void
+  ManifoldWrapper::create_spherical(const PointWrapper center)
+  {
+    if ((dim == 2) && (spacedim == 2))
+      {
+        manifold_ptr = internal::create_spherical_manifold<2, 2>(center);
+      }
+    else if ((dim == 2) && (spacedim == 3))
+      {
+        manifold_ptr = internal::create_spherical_manifold<2, 3>(center);
+      }
+    else if ((dim == 3) && (spacedim == 3))
+      {
+        manifold_ptr = internal::create_spherical_manifold<3, 3>(center);
+      }
+    else
+      AssertThrow(false, ExcMessage("Wrong dim-spacedim combination."));
+  }
+
+
+
+  void
+  ManifoldWrapper::create_polar(const PointWrapper center)
+  {
+    if ((dim == 2) && (spacedim == 2))
+      {
+        manifold_ptr = internal::create_polar_manifold<2, 2>(center);
+      }
+    else if ((dim == 2) && (spacedim == 3))
+      {
+        manifold_ptr = internal::create_polar_manifold<2, 3>(center);
+      }
+    else if ((dim == 3) && (spacedim == 3))
+      {
+        manifold_ptr = internal::create_polar_manifold<3, 3>(center);
+      }
+    else
+      AssertThrow(false, ExcMessage("Wrong dim-spacedim combination."));
+  }
+
+
+
+  void
+  ManifoldWrapper::create_cylindrical(const int axis, const double tolerance)
+  {
+    if ((dim == 2) && (spacedim == 3))
+      {
+        manifold_ptr =
+          internal::create_cylindrical_manifold<2, 3>(axis, tolerance);
+      }
+    else if ((dim == 3) && (spacedim == 3))
+      {
+        manifold_ptr =
+          internal::create_cylindrical_manifold<3, 3>(axis, tolerance);
+      }
+    else
+      AssertThrow(false, ExcMessage("Wrong dim-spacedim combination."));
+  }
+
+
+  void *
+  ManifoldWrapper::get_manifold() const
+  {
+    return manifold_ptr;
+  }
+
+} // namespace python
+
+DEAL_II_NAMESPACE_CLOSE
index 3e28e491045bb7f10e39e380c27960d1e672ddf9..45b1024ea820df205d15c151724cb58d19a6b24b 100644 (file)
@@ -1266,6 +1266,39 @@ namespace python
   }
 
 
+
+  void
+  TriangulationWrapper::set_manifold(const int        number,
+                                     ManifoldWrapper &manifold)
+  {
+    if ((dim == 2) && (spacedim == 2))
+      {
+        Triangulation<2, 2> *tria =
+          static_cast<Triangulation<2, 2> *>(triangulation);
+        Manifold<2, 2> *m =
+          static_cast<Manifold<2, 2> *>(manifold.get_manifold());
+        tria->set_manifold(number, *m);
+      }
+    else if ((dim == 2) && (spacedim == 3))
+      {
+        Triangulation<2, 3> *tria =
+          static_cast<Triangulation<2, 3> *>(triangulation);
+        Manifold<2, 3> *m =
+          static_cast<Manifold<2, 3> *>(manifold.get_manifold());
+        tria->set_manifold(number, *m);
+      }
+    else
+      {
+        Triangulation<3, 3> *tria =
+          static_cast<Triangulation<3, 3> *>(triangulation);
+        Manifold<3, 3> *m =
+          static_cast<Manifold<3, 3> *>(manifold.get_manifold());
+        tria->set_manifold(number, *m);
+      }
+  }
+
+
+
   void
   TriangulationWrapper::setup(const std::string &dimension,
                               const std::string &spacedimension)
index 5d26719cd0c2ebaa6908133a712f67bbb7789e5f..aa18628b452de32fba0f4bda4cc32c0e06b58c68 100644 (file)
@@ -33,6 +33,8 @@ namespace python
   export_triangulation();
   void
   export_mapping();
+  void
+  export_manifold();
 } // namespace python
 
 DEAL_II_NAMESPACE_CLOSE
@@ -66,6 +68,7 @@ BOOST_PYTHON_MODULE(Debug)
   dealii::python::export_point();
   dealii::python::export_triangulation();
   dealii::python::export_mapping();
+  dealii::python::export_manifold();
 }
 
 #else
@@ -84,6 +87,7 @@ BOOST_PYTHON_MODULE(Release)
   dealii::python::export_point();
   dealii::python::export_triangulation();
   dealii::python::export_mapping();
+  dealii::python::export_manifold();
 }
 
 #endif
diff --git a/doc/news/changes/minor/20191128Grayver b/doc/news/changes/minor/20191128Grayver
new file mode 100644 (file)
index 0000000..dfb0c89
--- /dev/null
@@ -0,0 +1,5 @@
+Added: New python wrappers for the Manifold class
+with a possibility to assign manifold object to a
+python triangulation object.
+<br>
+(Alexander Grayver, 2019/11/28)

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.