]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add a series of new python wrappers
authorAlexander Grayver <agrayver@erdw.ethz.ch>
Fri, 6 Dec 2019 15:04:59 +0000 (16:04 +0100)
committerAlexander Grayver <agrayver@erdw.ethz.ch>
Sun, 15 Dec 2019 12:44:31 +0000 (13:44 +0100)
contrib/python-bindings/include/cell_accessor_wrapper.h
contrib/python-bindings/include/mapping_wrapper.h
contrib/python-bindings/include/triangulation_wrapper.h
contrib/python-bindings/source/cell_accessor_wrapper.cc
contrib/python-bindings/source/export_cell_accessor.cc
contrib/python-bindings/source/export_triangulation.cc
contrib/python-bindings/source/mapping_wrapper.cc
contrib/python-bindings/source/triangulation_wrapper.cc
doc/news/changes/minor/20191206Grayver [new file with mode: 0644]

index fd2d2bb284e1c7a7094920dbbb75da4b7e4ccd1d..0c7f471a0650e2cea1a8541ea502ef2c2797a51d 100644 (file)
@@ -89,6 +89,13 @@ namespace python
     PointWrapper
     get_barycenter() const;
 
+    /**
+     * Get the barycenter of the cell.
+     */
+    PointWrapper
+    get_center(const bool respect_manifold             = false,
+               const bool interpolate_from_surrounding = false) const;
+
     /**
      * Set the material id.
      */
index b80752ce1480fb1a54e16c2918dfad4b48cb35a0..a51d437b64381096b98d5ffd05cd331012cb7ef2 100644 (file)
@@ -61,6 +61,12 @@ namespace python
     PointWrapper
     transform_real_to_unit_cell(CellAccessorWrapper &cell, PointWrapper &point);
 
+    /**
+     * Get the underlying mapping
+     */
+    void *
+    get_mapping();
+
   private:
     /**
      * Dimension of the underlying Mapping object.
index 6910143fc7354f01cfcc8679be0e2956c703ad83..63f28e0d7f32240125c42c556ccfe34da4d7c103 100644 (file)
@@ -21,6 +21,7 @@
 #include <boost/python.hpp>
 
 #include <manifold_wrapper.h>
+#include <mapping_wrapper.h>
 #include <point_wrapper.h>
 
 #include <string>
@@ -316,6 +317,44 @@ namespace python
     void
     flatten_triangulation(TriangulationWrapper &tria_out);
 
+    /**
+     * Take a 2d Triangulation that is being extruded in z direction by
+     * the total height of height using n_slices slices (minimum is 2).
+     * The boundary indicators of the faces of input are going to be
+     * assigned to the corresponding side walls in z direction. The
+     * bottom and top get the next two free boundary indicators.
+     */
+    void
+    extrude_triangulation(const unsigned int    n_slices,
+                          const double          height,
+                          TriangulationWrapper &tria_out);
+
+    /**
+     * Distort the given triangulation by randomly moving around all the
+     * vertices of the grid. The direction of movement of each vertex is
+     * random, while the length of the shift vector has a value of factor
+     * times the minimal length of the active edges adjacent to this vertex.
+     * Note that factor should obviously be well below 0.5.
+     */
+    void
+    distort_random(const double factor, const bool keep_boundary = true);
+
+    /**
+     * Transform the vertices of the given triangulation by applying the
+     * function object provided as first argument to all its vertices.
+     */
+    void
+    transform(boost::python::object &transformation);
+
+    /**
+     * Find and return an active cell that surrounds a given point p.
+     * The mapping used to determine whether the given point is inside a given
+     * cell.
+     */
+    CellAccessorWrapper
+    find_active_cell_around_point(PointWrapper &          p,
+                                  MappingQGenericWrapper &mapping);
+
     /**
      * Assign a manifold object to a certain part of the triangulation.
      * If an object with manifold number is refined, this object
index 73991421da3531bcd81a7b10ebf987b169e6d9de..d8ede1bafb6573c32219ac92511775bd04a53536 100644 (file)
@@ -173,6 +173,22 @@ namespace python
       return PointWrapper(barycenter_list);
     }
 
+    template <int dim, int spacedim>
+    PointWrapper
+    get_center(const bool  respect_manifold,
+               const bool  interpolate_from_surrounding,
+               const void *cell_accessor)
+    {
+      const CellAccessor<dim, spacedim> *cell =
+        static_cast<const CellAccessor<dim, spacedim> *>(cell_accessor);
+      Point<spacedim> center =
+        cell->center(respect_manifold, interpolate_from_surrounding);
+      boost::python::list center_list;
+      for (int i = 0; i < spacedim; ++i)
+        center_list.append(center[i]);
+
+      return PointWrapper(center_list);
+    }
 
 
     template <int dim, int spacedim>
@@ -496,6 +512,26 @@ namespace python
 
 
 
+  PointWrapper
+  CellAccessorWrapper::get_center(const bool respect_manifold,
+                                  const bool interpolate_from_surrounding) const
+  {
+    if ((dim == 2) && (spacedim == 2))
+      return internal::get_center<2, 2>(respect_manifold,
+                                        interpolate_from_surrounding,
+                                        cell_accessor);
+    else if ((dim == 2) && (spacedim == 3))
+      return internal::get_center<2, 3>(respect_manifold,
+                                        interpolate_from_surrounding,
+                                        cell_accessor);
+    else
+      return internal::get_center<3, 3>(respect_manifold,
+                                        interpolate_from_surrounding,
+                                        cell_accessor);
+  }
+
+
+
   void
   CellAccessorWrapper::set_material_id(const int material_id)
   {
index 877da21b98fba6ffb4b2caf02298fcaa3944118b..6e088447db81ffe3b91c107df6f2f7fd193047bb 100644 (file)
@@ -22,6 +22,8 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace python
 {
+  BOOST_PYTHON_MEMBER_FUNCTION_OVERLOADS(get_center_overloads, get_center, 0, 2)
+
   const char refine_flag_docstring[] =
     "Get/Set the refine_flag of the cell. In 2D, the possibilities are: \n"
     "  - isotropic                                                      \n"
@@ -68,6 +70,12 @@ namespace python
 
 
 
+  const char center_docstring[] =
+    "Center of the object. This method can take into account manifold   \n"
+    "attached to a cell.                                                \n";
+
+
+
   const char set_vertex_docstring[] =
     " Set the ith vertex of the cell to point_wrapper                   \n";
 
@@ -138,6 +146,13 @@ namespace python
            &CellAccessorWrapper::get_barycenter,
            barycenter_docstring,
            boost::python::args("self"))
+      .def("center",
+           &CellAccessorWrapper::get_center,
+           get_center_overloads(
+             boost::python::args("self",
+                                 "respect_manifold",
+                                 "interpolate_from_surrounding"),
+             center_docstring))
       .def("set_vertex",
            &CellAccessorWrapper::set_vertex,
            set_vertex_docstring,
index 59195628996eecad72c4b7a3aa872269d52ae644..65835a57bb5e5cdcbffe48cc59ae687d52751597 100644 (file)
@@ -97,6 +97,10 @@ namespace python
                                          generate_hyper_shell,
                                          3,
                                          5)
+  BOOST_PYTHON_MEMBER_FUNCTION_OVERLOADS(distort_random_overloads,
+                                         distort_random,
+                                         1,
+                                         2)
 
   const char n_active_cells_docstring[] =
     "Return the number of active cells                                      \n";
@@ -283,6 +287,15 @@ namespace python
 
 
 
+  const char extrude_docstring[] =
+    "Take a 2d Triangulation that is being extruded in z direction by       \n"
+    "the total height of height using n_slices slices (minimum is 2).       \n"
+    "The boundary indicators of the faces of input are going to be          \n"
+    "assigned to the corresponding side walls in z direction. The           \n"
+    "bottom and top get the next two free boundary indicators.              \n";
+
+
+
   const char flatten_triangulation_docstring[] =
     "Create a new flat triangulation out_tria which contains a single       \n"
     "level with all active cells of the input triangulation. If the spacedim\n"
@@ -303,6 +316,15 @@ namespace python
 
 
 
+  const char distort_random_docstring[] =
+    "Distort the given triangulation by randomly moving around all the      \n"
+    "vertices of the grid. The direction of movement of each vertex is      \n"
+    "random, while the length of the shift vector has a value of factor     \n"
+    "times the minimal length of the active edges adjacent to this vertex.  \n"
+    "Note that factor should obviously be well below 0.5.                   \n";
+
+
+
   const char refine_global_docstring[] =
     "Refine all the cells times time                                        \n";
 
@@ -365,6 +387,17 @@ namespace python
 
 
 
+  const char transform_docstring[] =
+    "Transform the vertices of the given triangulation by applying the      \n"
+    "function object provided as first argument to all its vertices.        \n";
+
+
+
+  const char find_active_cell_around_point_docstring[] =
+    "Find and return an active cell that surrounds a given point p.         \n";
+
+
+
   void
   export_triangulation()
   {
@@ -484,10 +517,27 @@ namespace python
            &TriangulationWrapper::merge_triangulations,
            merge_docstring,
            boost::python::args("self", "triangulation_1", "triangulation_2"))
+      .def("extrude_triangulation",
+           &TriangulationWrapper::extrude_triangulation,
+           extrude_docstring,
+           boost::python::args("self", "n_slices", "depth", "tria_out"))
       .def("flatten_triangulation",
            &TriangulationWrapper::flatten_triangulation,
            flatten_triangulation_docstring,
            boost::python::args("self", "tria_out"))
+      .def("distort_random",
+           &TriangulationWrapper::distort_random,
+           distort_random_overloads(
+             boost::python::args("self", "factor", "keep_boundary"),
+             distort_random_docstring))
+      .def("transform",
+           &TriangulationWrapper::transform,
+           transform_docstring,
+           boost::python::args("self", "transformation"))
+      .def("find_active_cell_around_point",
+           &TriangulationWrapper::find_active_cell_around_point,
+           find_active_cell_around_point_docstring,
+           boost::python::args("self", "point", "mapping"))
       .def("refine_global",
            &TriangulationWrapper::refine_global,
            refine_global_docstring,
index 0c7ad75c675cb5e2d4ff35cc6d4d128a13c34bdf..97c317fe243d2f8b4f2791651f2bbbf522860dcd 100644 (file)
@@ -222,6 +222,14 @@ namespace python
                                                          p.point);
   }
 
+
+
+  void *
+  MappingQGenericWrapper::get_mapping()
+  {
+    return mapping_ptr;
+  }
+
 } // namespace python
 
 DEAL_II_NAMESPACE_CLOSE
index 1123f7e42b0bc0e37c3ed3e3e532a369e8c5bea0..a9adb51806eef4830eec8bad1fbbcfc940a62a8e 100644 (file)
@@ -484,6 +484,73 @@ namespace python
 
 
 
+    template <int dim, int spacedim>
+    void
+    distort_random(const double factor,
+                   const bool   keep_boundary,
+                   void *       triangulation)
+    {
+      Triangulation<dim, spacedim> *tria =
+        static_cast<Triangulation<dim, spacedim> *>(triangulation);
+
+      GridTools::distort_random(factor, *tria, keep_boundary);
+    }
+
+
+
+    template <int dim, int spacedim>
+    void
+    transform(boost::python::object &transformation, void *triangulation)
+    {
+      Triangulation<dim, spacedim> *tria =
+        static_cast<Triangulation<dim, spacedim> *>(triangulation);
+
+      auto func = [&transformation](const Point<spacedim> &p_in) {
+        boost::python::list p_list_in;
+
+        for (int d = 0; d < spacedim; ++d)
+          p_list_in.append(p_in[d]);
+
+        boost::python::list p_list_out =
+          boost::python::extract<boost::python::list>(
+            transformation(p_list_in));
+
+        Point<spacedim> p_out;
+
+        for (int d = 0; d < spacedim; ++d)
+          p_out[d] = boost::python::extract<double>(p_list_out[d]);
+
+        return p_out;
+      };
+
+      GridTools::transform(func, *tria);
+    }
+
+
+
+    template <int dim, int spacedim>
+    std::pair<int, int>
+    find_active_cell_around_point(PointWrapper &          p,
+                                  MappingQGenericWrapper &mapping_wrapper,
+                                  void *                  triangulation)
+    {
+      Triangulation<dim, spacedim> *tria =
+        static_cast<Triangulation<dim, spacedim> *>(triangulation);
+
+      Point<spacedim> point = *(static_cast<Point<spacedim> *>(p.get_point()));
+
+      const MappingQGeneric<dim, spacedim> *mapping =
+        static_cast<const MappingQGeneric<dim, spacedim> *>(
+          mapping_wrapper.get_mapping());
+
+      auto cell_pair =
+        GridTools::find_active_cell_around_point(*mapping, *tria, point);
+
+      return std::make_pair(cell_pair.first->level(), cell_pair.first->index());
+    }
+
+
+
     template <int dim, int spacedim>
     boost::python::list
     active_cells(TriangulationWrapper &triangulation_wrapper)
@@ -1199,6 +1266,91 @@ namespace python
 
 
 
+  void
+  TriangulationWrapper::extrude_triangulation(
+    const unsigned int    n_slices,
+    const double          height,
+    TriangulationWrapper &triangulation_out)
+  {
+    AssertThrow((dim == 2) && (spacedim == 2),
+                ExcMessage(
+                  "Extrude can only be applied to the dim and spacedim two."));
+    AssertThrow((triangulation_out.get_dim() == 3) &&
+                  (triangulation_out.get_spacedim() == 3),
+                ExcMessage(
+                  "The output Triangulation must be of dimension three"));
+
+    if ((dim == 2) && (spacedim == 2))
+      {
+        Triangulation<2, 2> *tria =
+          static_cast<Triangulation<2, 2> *>(triangulation);
+        Triangulation<3, 3> *tria_out = static_cast<Triangulation<3, 3> *>(
+          triangulation_out.get_triangulation());
+        GridGenerator::extrude_triangulation(*tria,
+                                             n_slices,
+                                             height,
+                                             *tria_out);
+      }
+  }
+
+
+
+  void
+  TriangulationWrapper::distort_random(const double factor,
+                                       const bool   keep_boundary)
+  {
+    if ((dim == 2) && (spacedim == 2))
+      internal::distort_random<2, 2>(factor, keep_boundary, triangulation);
+    else if ((dim == 2) && (spacedim == 3))
+      internal::distort_random<2, 3>(factor, keep_boundary, triangulation);
+    else
+      internal::distort_random<3, 3>(factor, keep_boundary, triangulation);
+  }
+
+
+
+  void
+  TriangulationWrapper::transform(boost::python::object &transformation)
+  {
+    if ((dim == 2) && (spacedim == 2))
+      internal::transform<2, 2>(transformation, triangulation);
+    else if ((dim == 2) && (spacedim == 3))
+      internal::transform<2, 3>(transformation, triangulation);
+    else
+      internal::transform<3, 3>(transformation, triangulation);
+  }
+
+
+
+  CellAccessorWrapper
+  TriangulationWrapper::find_active_cell_around_point(
+    PointWrapper &          p,
+    MappingQGenericWrapper &mapping)
+  {
+    std::pair<int, int> level_index_pair;
+    if ((dim == 2) && (spacedim == 2))
+      level_index_pair =
+        internal::find_active_cell_around_point<2, 2>(p,
+                                                      mapping,
+                                                      triangulation);
+    else if ((dim == 2) && (spacedim == 3))
+      level_index_pair =
+        internal::find_active_cell_around_point<2, 3>(p,
+                                                      mapping,
+                                                      triangulation);
+    else
+      level_index_pair =
+        internal::find_active_cell_around_point<3, 3>(p,
+                                                      mapping,
+                                                      triangulation);
+
+    return CellAccessorWrapper(*this,
+                               level_index_pair.first,
+                               level_index_pair.second);
+  }
+
+
+
   void
   TriangulationWrapper::refine_global(const unsigned int n)
   {
diff --git a/doc/news/changes/minor/20191206Grayver b/doc/news/changes/minor/20191206Grayver
new file mode 100644 (file)
index 0000000..4b1de45
--- /dev/null
@@ -0,0 +1,5 @@
+Added: Python wrappers for the CellAccessor::center, GridTools::transform,
+GridTools::distort_random, GridGenerator::extrude_triangulation, 
+GridTools::find_active_cell_around_point
+<br>
+(Alexander Grayver, 2019/12/06)

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.