#include <boost/python.hpp>
#include <deal.II/python/point_wrapper.h>
#include <string>
+#include <vector>
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.
*/
*/
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.
*/
* Pointer that can be casted to the underlying Triangulation object.
*/
void *triangulation;
-
- /**
- * Vector of all the cells of the Triangulation.
- */
- std::vector<CellAccessorWrapper> cells;
};
{
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",
&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.",
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;
}
- 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)
{
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;
- }
}
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])
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)
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()