]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Replace __iter__ by active_cells.
authorBruno Turcksin <bruno.turcksin@gmail.com>
Tue, 28 Jun 2016 21:19:06 +0000 (17:19 -0400)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Wed, 3 Aug 2016 20:34:02 +0000 (16:34 -0400)
include/deal.II/python/triangulation_wrapper.h
source/python/export_triangulation.cc
source/python/triangulation_wrapper.cc
tests/python/cell_accessor_wrapper.py
tests/python/triangulation_wrapper.py

index fceaef00803ca91c76aa18685b87940300f98a00..c29a242eea658e9f014cc86c2518ef3dbb006a66 100644 (file)
@@ -19,6 +19,7 @@
 #include <boost/python.hpp>
 #include <deal.II/python/point_wrapper.h>
 #include <string>
+#include <vector>
 
 namespace PyDealII
 {
@@ -121,6 +122,12 @@ namespace PyDealII
      */
     void execute_coarsening_and_refinement();
 
+    /**
+     * Return the list of active cell accessors associated to the underlying
+     * Triangulation.
+     */
+    boost::python::list active_cells();
+
     /**
      * Write mesh to the output file @filename according to the given data format.
      */
@@ -152,28 +159,12 @@ namespace PyDealII
      */
     void *get_triangulation();
 
-    /**
-     * Iterator to the first active cell of the underlying Triangulation.
-     */
-    iterator begin();
-
-    /**
-     * Iterator past the end.
-     */
-    iterator end();
-
   private:
     /**
      * Helper function for the contructors.
      */
     void setup(const std::string &dimension, const std::string &spacedimension);
 
-    /**
-     * Build the vector of CellAccessorWrapper. This allows access to the
-     * CellAccessor associated to the underlying Triangulation.
-     */
-    void build_cells_vector();
-
     /**
      * Flag is true if the vector cells has been updated.
      */
@@ -193,11 +184,6 @@ namespace PyDealII
      * Pointer that can be casted to the underlying Triangulation object.
      */
     void *triangulation;
-
-    /**
-     * Vector of all the cells of the Triangulation.
-     */
-    std::vector<CellAccessorWrapper> cells;
   };
 
 
index 8bc6cdca319b9b3f43b0ba1244ad5e420e92d266..535119cdecd586af7cb8a266c45c9fcbb3bc1c9e 100644 (file)
@@ -35,7 +35,6 @@ namespace PyDealII
   {
     boost::python::class_<TriangulationWrapper>("Triangulation", boost::python::init<const std::string &>())
     .def(boost::python::init<const std::string &, const std::string &>())
-    .def("__iter__", boost::python::iterator<TriangulationWrapper>())
     .def("n_active_cells",
          &TriangulationWrapper::n_active_cells,
          "Return the number of active cells",
@@ -86,6 +85,10 @@ namespace PyDealII
          &TriangulationWrapper::execute_coarsening_and_refinement,
          "Execute both refinement and coarsening of the Triangulation.",
          boost::python::args("self"))
+    .def("active_cells",
+         &TriangulationWrapper::active_cells,
+         "Return the list of active cell accessors of the Triangulation.",
+         boost::python::args("self"))
     .def("write",
          &TriangulationWrapper::write,
          "Write the mesh to the output file according to the given data format. The possible formats are: none, dx, gnuplot, eps, ucd, xfig, msh, svg, mathgl, vtk, and vtu.",
index 3a749cd5a9d9c8402fd00bf97f7b86a5a09318d4..61efcbcaee59d60b971ce3ce37c90c259fb982ea 100644 (file)
@@ -172,15 +172,17 @@ namespace PyDealII
 
 
     template <int dim, int spacedim>
-    void build_cells_vector(TriangulationWrapper             &triangulation_wrapper,
-                            std::vector<CellAccessorWrapper> &cells)
+    boost::python::list active_cells(TriangulationWrapper &triangulation_wrapper)
     {
       dealii::Triangulation<dim,spacedim> *tria =
         static_cast<dealii::Triangulation<dim,spacedim>*>(
           triangulation_wrapper.get_triangulation());
+      boost::python::list cells;
       for (auto cell : tria->active_cell_iterators())
-        cells.push_back(CellAccessorWrapper(triangulation_wrapper, cell->level(),
-                                            cell->index()));
+        cells.append(CellAccessorWrapper(triangulation_wrapper, cell->level(),
+                                         cell->index()));
+
+      return cells;
     }
 
 
@@ -558,30 +560,17 @@ namespace PyDealII
 
 
 
-  TriangulationWrapper::iterator TriangulationWrapper::begin()
+  boost::python::list TriangulationWrapper::active_cells()
   {
-    // Build the vector containing all the CellAccessorWrapper when the iterator
-    // is called. This way, we don't need to keep track of the changes on the
-    // Triangulation.
-    build_cells_vector();
-
-    return cells.begin();
-  }
-
-
-
-  TriangulationWrapper::iterator TriangulationWrapper::end()
-  {
-    // Python seems to call end() before begin() so call build_cells_vector in
-    // both begin() and end(). The updated flag makes sure that the vector of
-    // CellAccessorWrapper is only build once.
-    build_cells_vector();
-
-    return cells.end();
+    if ((dim == 2) && (spacedim == 2))
+      return internal::active_cells<2,2>(*this);
+    else if ((dim == 2) && (spacedim == 3))
+      return internal::active_cells<2,3>(*this);
+    else
+      return internal::active_cells<3,3>(*this);
   }
 
 
-
   void TriangulationWrapper::setup(const std::string &dimension,
                                    const std::string &spacedimension)
   {
@@ -616,23 +605,4 @@ namespace PyDealII
     else
       AssertThrow(false, dealii::ExcMessage("Dimension needs to be 2D or 3D."));
   }
-
-
-
-  void TriangulationWrapper::build_cells_vector()
-  {
-    if (updated == false)
-      {
-        cells.clear();
-        if ((dim == 2) && (spacedim == 2))
-          internal::build_cells_vector<2,2>(*this, cells);
-        else if ((dim == 2) && (spacedim == 3))
-          internal::build_cells_vector<2,3>(*this, cells);
-        else
-          internal::build_cells_vector<3,3>(*this, cells);
-        updated = true;
-      }
-    else
-      updated = false;
-  }
 }
index 07792aa311d23f2647034f66695dd0bc19391c98..8f8308b885efecaff35fc1fe5b92d4921238b9d0 100644 (file)
@@ -26,49 +26,49 @@ class TestCellAccessorWrapper(unittest.TestCase):
 
     def test_material_id(self):
         material_id = 0
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             cell.material_id = material_id
             material_id += 1
         material_id = 0
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             self.assertEqual(cell.material_id, material_id)
             material_id += 1
 
     def test_manifold_id(self):
         manifold_id = 0
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             cell.manifold_id = manifold_id
             manifold_id += 1
         manifold_id = 0
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             self.assertEqual(cell.manifold_id, manifold_id)
             manifold_id += 1
 
     def test_refine_flag(self):
         index = 0
         refine_flags = ['no_refinement', 'cut_x', 'cut_y', 'cut_xy']
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             cell.refine_flag = refine_flags[index]
             index += 1
         index = 0
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             self.assertEqual(cell.refine_flag, refine_flags[index])
             index += 1
 
     def test_coarsen_flag(self):
         coarsen_flag = True
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             cell.coarsen_flag = coarsen_flag
             coarsen_flag = not coarsen_flag
         coarsen_flag = True
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             self.assertEqual(cell.coarsen_flag, coarsen_flag)
             coarsen_flag = not coarsen_flag
 
     def test_barycenter(self):
         centers = [[0.25, 0.25], [0.75, 0.25], [0.25, 0.75], [0.75, 0.75]]
         index = 0
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             barycenter = cell.barycenter()
             self.assertEqual(barycenter.x, centers[index][0])
             self.assertEqual(barycenter.y, centers[index][1])
@@ -76,13 +76,13 @@ class TestCellAccessorWrapper(unittest.TestCase):
 
     def test_move_vertex(self):
         point = Point([0.6, 0.6])
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             cell.set_vertex(3, point)
             vertex = cell.get_vertex(3)
             break
         vertices = [3, 2, 1, 0]
         index = 0
-        for cell in self.triangulation:
+        for cell in self.triangulation.active_cells():
             vertex = cell.get_vertex(vertices[index])
             self.assertEqual(vertex.x, point.x)
             self.assertEqual(vertex.y, point.y)
index 964c8ddb1aaf1a6e9af89f318f7e7ddb3c42fc8a..52ac71f7d7ed9471ee8c47200bc43c91aeefa1fa 100644 (file)
@@ -127,7 +127,7 @@ class TestTriangulationWrapper(unittest.TestCase):
         for dim in self.dim:
             triangulation = self.build_hyper_cube_triangulation(dim)
             triangulation.refine_global(1)
-            for cell in triangulation:
+            for cell in triangulation.active_cells():
                 cell.refine_flag = 'isotropic'
                 break
             triangulation.execute_coarsening_and_refinement()

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.