]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add tests 9148/head
authorAlexander Grayver <agrayver@erdw.ethz.ch>
Fri, 13 Dec 2019 14:21:01 +0000 (15:21 +0100)
committerAlexander Grayver <agrayver@erdw.ethz.ch>
Sun, 15 Dec 2019 19:37:11 +0000 (20:37 +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_triangulation.cc
contrib/python-bindings/source/mapping_wrapper.cc
contrib/python-bindings/source/triangulation_wrapper.cc
contrib/python-bindings/tests/triangulation_wrapper.py

index 0c7f471a0650e2cea1a8541ea502ef2c2797a51d..ef578822562311ba772122f2c5a29f4677e31ec1 100644 (file)
@@ -90,7 +90,7 @@ namespace python
     get_barycenter() const;
 
     /**
-     * Get the barycenter of the cell.
+     * Get the center of the cell.
      */
     PointWrapper
     get_center(const bool respect_manifold             = false,
index a51d437b64381096b98d5ffd05cd331012cb7ef2..f9d3d8e1c7f4e02bc15dab39c8d9e30e0e6befab 100644 (file)
@@ -37,6 +37,11 @@ namespace python
      */
     MappingQGenericWrapper(const MappingQGenericWrapper &other);
 
+    /**
+     * Default constructor.
+     */
+    MappingQGenericWrapper();
+
     /**
      * Constructor.
      */
@@ -62,7 +67,7 @@ namespace python
     transform_real_to_unit_cell(CellAccessorWrapper &cell, PointWrapper &point);
 
     /**
-     * Get the underlying mapping
+     * Get the underlying mapping.
      */
     void *
     get_mapping();
index 63f28e0d7f32240125c42c556ccfe34da4d7c103..bb2505c62847508645b6b1b39613496fb92529e4 100644 (file)
@@ -348,12 +348,13 @@ namespace python
 
     /**
      * 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.
+     * The mapping is used to determine whether the given point is inside a
+     * given cell.
      */
     CellAccessorWrapper
-    find_active_cell_around_point(PointWrapper &          p,
-                                  MappingQGenericWrapper &mapping);
+    find_active_cell_around_point(
+      PointWrapper &         p,
+      MappingQGenericWrapper mapping = MappingQGenericWrapper());
 
     /**
      * Assign a manifold object to a certain part of the triangulation.
index d8ede1bafb6573c32219ac92511775bd04a53536..c17981ed45eb31daa541b46f7d6a93e6ca6eed50 100644 (file)
@@ -173,6 +173,8 @@ namespace python
       return PointWrapper(barycenter_list);
     }
 
+
+
     template <int dim, int spacedim>
     PointWrapper
     get_center(const bool  respect_manifold,
index 65835a57bb5e5cdcbffe48cc59ae687d52751597..cb1e00097325aaef20d0908f5b467454de591d15 100644 (file)
@@ -101,6 +101,11 @@ namespace python
                                          distort_random,
                                          1,
                                          2)
+  BOOST_PYTHON_MEMBER_FUNCTION_OVERLOADS(
+    find_active_cell_around_point_overloads,
+    find_active_cell_around_point,
+    1,
+    2)
 
   const char n_active_cells_docstring[] =
     "Return the number of active cells                                      \n";
@@ -536,8 +541,9 @@ namespace python
            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"))
+           find_active_cell_around_point_overloads(
+             boost::python::args("self", "point", "mapping"),
+             find_active_cell_around_point_docstring))
       .def("refine_global",
            &TriangulationWrapper::refine_global,
            refine_global_docstring,
index 97c317fe243d2f8b4f2791651f2bbbf522860dcd..edf344b3d6d925d0d8525d41f8ff5b6c4a8e697f 100644 (file)
@@ -87,6 +87,15 @@ namespace python
 
 
 
+  MappingQGenericWrapper::MappingQGenericWrapper()
+    : dim(-1)
+    , spacedim(-1)
+    , degree(-1)
+    , mapping_ptr(nullptr)
+  {}
+
+
+
   MappingQGenericWrapper::MappingQGenericWrapper(const int dim,
                                                  const int spacedim,
                                                  const int degree)
index a9adb51806eef4830eec8bad1fbbcfc940a62a8e..0805fd8e56d0fa190e55e6a2fe5ebf81416c49e5 100644 (file)
@@ -539,14 +539,22 @@ namespace python
 
       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());
+      if (mapping_wrapper.get_mapping() != nullptr)
+        {
+          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());
+        }
+      else
+        {
+          auto cell = GridTools::find_active_cell_around_point(*tria, point);
+          return std::make_pair(cell->level(), cell->index());
+        }
     }
 
 
@@ -1280,17 +1288,12 @@ namespace python
                 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);
-      }
+
+    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);
   }
 
 
@@ -1324,8 +1327,8 @@ namespace python
 
   CellAccessorWrapper
   TriangulationWrapper::find_active_cell_around_point(
-    PointWrapper &          p,
-    MappingQGenericWrapper &mapping)
+    PointWrapper &         p,
+    MappingQGenericWrapper mapping)
   {
     std::pair<int, int> level_index_pair;
     if ((dim == 2) && (spacedim == 2))
index d23d8c9aadee3c3c84554066153f00fb0efff72c..890f9050661d06bea6e689d9c79fc02a697eaf29 100644 (file)
@@ -16,7 +16,6 @@
 import unittest
 from PyDealII.Debug import *
 
-
 class TestTriangulationWrapper(unittest.TestCase):
 
     def setUp(self):
@@ -375,6 +374,35 @@ class TestTriangulationWrapper(unittest.TestCase):
             else:
                 self.assertEqual(n_cells, 64)
 
+
+    def test_transform(self):
+        for dim in self.dim:
+            triangulation_1 = self.build_hyper_cube_triangulation(dim)
+            triangulation_1.refine_global(1)
+            triangulation_2 = self.build_hyper_cube_triangulation(dim)
+            triangulation_2.refine_global(1)
+
+            triangulation_1.transform(lambda p: [v + 1. for v in p])
+
+            if dim[1] == '3D':
+                offset = Point([1., 1., 1.])
+            else:
+                offset = Point([1., 1.])
+
+            for (cell_1, cell_2) in zip(triangulation_1.active_cells(), triangulation_2.active_cells()):
+                self.assertTrue(cell_1.center().distance(cell_2.center() + offset) < 1e-8)
+
+
+    def test_find_active_cell_around_point(self):
+        for dim in self.dim:
+            triangulation = self.build_hyper_cube_triangulation(dim)
+            triangulation.refine_global(2)
+
+            for cell in triangulation.active_cells():
+                cell_ret = triangulation.find_active_cell_around_point(cell.center())
+                self.assertTrue(cell.center().distance(cell_ret.center()) < 1e-8)
+
+
     def test_save_load(self):
         for dim in self.dim:
             triangulation_1 = self.build_hyper_cube_triangulation(dim)

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.