From 3b4aa186f5ed4aba334ec54c89ec0296673be48b Mon Sep 17 00:00:00 2001 From: Bruno Turcksin Date: Fri, 24 Jun 2016 16:47:20 -0400 Subject: [PATCH] Add CellAccessor wrapper. --- .../deal.II/python/cell_accessor_wrapper.h | 125 +++++ include/deal.II/python/point_wrapper.h | 1 + .../deal.II/python/triangulation_wrapper.h | 30 ++ source/python/CMakeLists.txt | 2 + source/python/cell_accessor_wrapper.cc | 488 ++++++++++++++++++ source/python/export_cell_accessor.cc | 49 ++ source/python/export_triangulation.cc | 4 +- source/python/point_wrapper.cc | 2 +- source/python/triangulation_wrapper.cc | 144 ++++-- source/python/wrappers.cc | 2 + tests/python/cell_accessor_wrapper.py | 93 ++++ 11 files changed, 897 insertions(+), 43 deletions(-) create mode 100644 include/deal.II/python/cell_accessor_wrapper.h create mode 100644 source/python/cell_accessor_wrapper.cc create mode 100644 source/python/export_cell_accessor.cc create mode 100644 tests/python/cell_accessor_wrapper.py diff --git a/include/deal.II/python/cell_accessor_wrapper.h b/include/deal.II/python/cell_accessor_wrapper.h new file mode 100644 index 0000000000..86262fdafb --- /dev/null +++ b/include/deal.II/python/cell_accessor_wrapper.h @@ -0,0 +1,125 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#ifndef dealii__cell_accessor_wrapper_h +#define dealii__cell_accessor_wrapper_h + +#include + +namespace PyDealII +{ + class PointWrapper; + class TriangulationWrapper; + + class CellAccessorWrapper + { + public: + /** + * Copy constructor. + */ + CellAccessorWrapper(const CellAccessorWrapper &other); + + /** + * Constructor. Takes a TriangulationWrapper, the level, and the index + * associated to the cell. + */ + CellAccessorWrapper(TriangulationWrapper &triangulation_wrapper, + const int level, + const int index); + + /** + * Destructor. + */ + ~CellAccessorWrapper(); + + /** + * Set the refine flag. The possibilities are: + * - in 2D: isotropic, no_refinement, cut_x, cut_y, and cut_xy. + * - in 3D: isotropic, no_refinement, cut_x, cut_y, cut_z, cut_xy, cut_xz, + * cut_yz, and cut_xyz. + */ + void set_refine_flag(const std::string &refinement_case); + + + /** + * Get the refine flag. + */ + std::string get_refine_flag(); + + /** + * Set the coarsen flag to true or false. + */ + void set_coarsen_flag(const bool coarsen_flag); + + /** + * Get the coarsen flag. + */ + bool get_coarsen_flag(); + + /** + * Get the barycenter of the cell. + */ + PointWrapper get_barycenter(); + + /** + * Set the material id. + */ + void set_material_id(const int material_id); + + /** + * Get the material id. + */ + int get_material_id(); + + /** + * Set the ith vertex of the cell to @p point_wrapper. + */ + void set_vertex(const int i, PointWrapper &point_wrapper); + + /** + * Return the ith vertex of the cell. + */ + PointWrapper get_vertex(const int i); + + /** + * Set the manifold id. + */ + void set_manifold_id(const int manifold_id); + + /** + * Get the manifold id. + */ + int get_manifold_id(); + + private: + /** + * Dimension of the underlying CellAccessor object. + */ + int dim; + + /** + * Space dimension of the underlying CellAccessor object. + */ + int spacedim; + + /** + * Pointer that can be casted to the underlying CellAccessor object. + */ + void *cell_accessor; + + }; +} + +#endif diff --git a/include/deal.II/python/point_wrapper.h b/include/deal.II/python/point_wrapper.h index 7fbeafc112..8d53b7564e 100644 --- a/include/deal.II/python/point_wrapper.h +++ b/include/deal.II/python/point_wrapper.h @@ -17,6 +17,7 @@ #define dealii__point_wrapper_h #include + namespace PyDealII { class PointWrapper diff --git a/include/deal.II/python/triangulation_wrapper.h b/include/deal.II/python/triangulation_wrapper.h index df4cb28617..b63e7472e6 100644 --- a/include/deal.II/python/triangulation_wrapper.h +++ b/include/deal.II/python/triangulation_wrapper.h @@ -22,9 +22,13 @@ namespace PyDealII { + class CellAccessorWrapper; + class TriangulationWrapper { public: + typedef std::vector::iterator iterator; + /** * Constructor. Takes a string @p dim with one of the following values * "2D", "2d", "3D", or "3d". @@ -138,12 +142,33 @@ 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. + */ + bool updated; + /** * Dimension of the underlying Triangulation object. */ @@ -158,6 +183,11 @@ namespace PyDealII * Pointer that can be casted to the underlying Triangulation object. */ void *triangulation; + + /** + * Vector of all the cells of the Triangulation. + */ + std::vector cells; }; diff --git a/source/python/CMakeLists.txt b/source/python/CMakeLists.txt index 6ee5c7261b..9955bd01b6 100644 --- a/source/python/CMakeLists.txt +++ b/source/python/CMakeLists.txt @@ -21,8 +21,10 @@ INCLUDE_DIRECTORIES(SYSTEM ${PYTHON_INCLUDE_DIRS}) SET(_src wrappers.cc + export_cell_accessor.cc export_point.cc export_triangulation.cc + cell_accessor_wrapper.cc point_wrapper.cc triangulation_wrapper.cc ) diff --git a/source/python/cell_accessor_wrapper.cc b/source/python/cell_accessor_wrapper.cc new file mode 100644 index 0000000000..f1afc7e775 --- /dev/null +++ b/source/python/cell_accessor_wrapper.cc @@ -0,0 +1,488 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#include +#include +#include +#include + +namespace PyDealII +{ + namespace internal + { + template + void set_refine_flag(const std::string &refinement_case, + void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + + std::unique_ptr> ref_case; + if (refinement_case.compare("isotropic") == 0) + ref_case.reset(new dealii::RefinementCase( + dealii::RefinementPossibilities::Possibilities::isotropic_refinement)); + else if (refinement_case.compare("no_refinement") == 0) + ref_case.reset(new dealii::RefinementCase( + dealii::RefinementPossibilities::Possibilities::no_refinement)); + else if (refinement_case.compare("cut_x") == 0) + ref_case.reset(new dealii::RefinementCase( + dealii::RefinementPossibilities::Possibilities::cut_x)); + else if (refinement_case.compare("cut_y") == 0) + ref_case.reset(new dealii::RefinementCase( + dealii::RefinementPossibilities::Possibilities::cut_y)); + else if (refinement_case.compare("cut_xy") == 0) + ref_case.reset(new dealii::RefinementCase( + dealii::RefinementPossibilities::Possibilities::cut_xy)); +#if dim==3 + else if (refinement_case.compare("cut_z") == 0) + { + AssertThrow(dim == 3, + dealii::ExcMessage("cut_z can only be use for dim equals to 3.")); + ref_case.reset(new dealii::RefinementCase<3>( + dealii::RefinementPossibilities<3>::Possibilities::cut_z)); + } + else if (refinement_case.compare("cut_xz") == 0) + { + AssertThrow(dim == 3, + dealii::ExcMessage("cut_xz can only be use for dim equals to 3.")); + ref_case.reset(new dealii::RefinementCase<3>( + dealii::RefinementPossibilities<3>::Possibilities::cut_xz)); + } + else if (refinement_case.compare("cut_yz") == 0) + { + AssertThrow(dim == 3, + dealii::ExcMessage("cut_yz can only be use for dim equals to 3.")); + ref_case.reset(new dealii::RefinementCase<3>( + dealii::RefinementPossibilities<3>::Possibilities::cut_yz)); + } + else if (refinement_case.compare("cut_xyz") == 0) + { + AssertThrow(dim == 3, + dealii::ExcMessage("cut_xyz can only be use for dim equals to 3.")); + ref_case.reset(new dealii::RefinementCase<3>( + dealii::RefinementPossibilities<3>::Possibilities::cut_xyz)); + } +#endif + else + AssertThrow(false, + dealii::ExcMessage("Unknown refinement possibility.")); + + cell->set_refine_flag(*ref_case); + } + + + + template + std::string get_refine_flag(void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + + std::string refine_flag; + dealii::RefinementCase ref_case = cell->refine_flag_set(); + switch (static_cast(ref_case)) + { + case (0) : + { + refine_flag = "no_refinement"; + break; + } + case (1) : + { + refine_flag = "cut_x"; + break; + } + case (2) : + { + refine_flag = "cut_y"; + break; + } + case (3) : + { + refine_flag = "cut_xy"; + break; + } + case (4) : + { + refine_flag = "cut_z"; + break; + } + case (5) : + { + refine_flag = "cut_xz"; + break; + } + case (6) : + { + refine_flag = "cut_yz"; + break; + } + case (7) : + { + refine_flag = "cut_xyz"; + break; + } + default : + { + AssertThrow(false, dealii::ExcMessage("Internal error.")); + } + } + + return refine_flag; + } + + + + template + void set_coarsen_flag(const bool coarsen_flag, void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + if (coarsen_flag == true) + cell->set_coarsen_flag(); + else + cell->clear_coarsen_flag(); + } + + + + template + bool get_coarsen_flag(void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + + return cell->coarsen_flag_set(); + } + + + + template + PointWrapper get_barycenter(void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + dealii::Point barycenter = cell->barycenter(); + boost::python::list barycenter_list; + for (int i=0; i + void set_material_id(const int material_id, void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + cell->set_material_id(material_id); + } + + + + template + int get_material_id(void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + + return cell->material_id(); + } + + + + template + void set_vertex(const int i, + PointWrapper &point_wrapper, + void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + dealii::Point *point = + static_cast*>(point_wrapper.get_point()); + + cell->vertex(i) = *point; + } + + + + template + PointWrapper get_vertex(const int i, void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + dealii::Point vertex = cell->vertex(i); + + boost::python::list coordinates; + for (int i=0; i + void set_manifold_id(const int manifold_id, void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + cell->set_manifold_id(manifold_id); + } + + + + template + int get_manifold_id(void *cell_accessor) + { + dealii::CellAccessor *cell = + static_cast*>(cell_accessor); + cell->manifold_id(); + } + } + + + + CellAccessorWrapper::CellAccessorWrapper(const CellAccessorWrapper &other) + : + dim(other.dim), + spacedim(other.dim) + { + if ((dim == 2) && (spacedim == 2)) + { + dealii::CellAccessor<2,2> *other_cell = + static_cast*>(other.cell_accessor); + cell_accessor = new dealii::CellAccessor<2,2>(*other_cell); + } + else if ((dim == 2) && (spacedim == 3)) + { + dealii::CellAccessor<2,3> *other_cell = + static_cast*>(other.cell_accessor); + cell_accessor = new dealii::CellAccessor<2,3>(*other_cell); + } + else if ((dim == 3) && (spacedim == 3)) + { + dealii::CellAccessor<3,3> *other_cell = + static_cast*>(other.cell_accessor); + cell_accessor = new dealii::CellAccessor<3,3>(*other_cell); + } + else + AssertThrow(false, dealii::ExcMessage("Wrong dim-spacedim combination.")); + } + + + + + CellAccessorWrapper::CellAccessorWrapper(TriangulationWrapper &triangulation_wrapper, + const int level, + const int index) + { + dim = triangulation_wrapper.get_dim(); + spacedim = triangulation_wrapper.get_spacedim(); + if ((dim == 2) && (spacedim == 2)) + { + dealii::Triangulation<2,2> *tmp = static_cast*> ( + triangulation_wrapper.get_triangulation()); + cell_accessor = new dealii::CellAccessor<2,2>(tmp ,level, index); + } + else if ((dim == 2) && (spacedim == 3)) + { + dealii::Triangulation<2,3> *tmp = static_cast*> ( + triangulation_wrapper.get_triangulation()); + cell_accessor = new dealii::CellAccessor<2,3>(tmp, level, index); + } + else if ((dim == 3) && (spacedim == 3)) + { + dealii::Triangulation<3,3> *tmp = static_cast*> ( + triangulation_wrapper.get_triangulation()); + cell_accessor = new dealii::CellAccessor<3,3>(tmp, level, index); + } + else + AssertThrow(false, dealii::ExcMessage("Wrong dim-spacedim combination.")); + } + + + + CellAccessorWrapper::~CellAccessorWrapper() + { + if (dim != -1) + { + if ((dim == 2) && (spacedim == 2)) + { + // We cannot call delete on a void pointer so cast the void pointer back + // first. + dealii::CellAccessor<2,2> *tmp = + static_cast*>(cell_accessor); + delete tmp; + } + else if ((dim == 2) && (spacedim == 3)) + { + dealii::CellAccessor<3,3> *tmp = + static_cast*>(cell_accessor); + delete tmp; + } + else + { + dealii::CellAccessor<3,3> *tmp = + static_cast*>(cell_accessor); + delete tmp; + } + + dim = -1; + spacedim = -1; + cell_accessor = nullptr; + } + } + + + + void CellAccessorWrapper::set_refine_flag(const std::string &refinement_case) + { + if ((dim == 2) && (spacedim == 2)) + internal::set_refine_flag<2,2>(refinement_case, cell_accessor); + else if ((dim == 2) && (spacedim == 3)) + internal::set_refine_flag<2,3>(refinement_case, cell_accessor); + else + internal::set_refine_flag<3,3>(refinement_case, cell_accessor); + } + + + + std::string CellAccessorWrapper::get_refine_flag() + { + if ((dim == 2) && (spacedim == 2)) + return internal::get_refine_flag<2,2>(cell_accessor); + else if ((dim == 2) && (spacedim == 3)) + return internal::get_refine_flag<2,3>(cell_accessor); + else + return internal::get_refine_flag<3,3>(cell_accessor); + } + + + + void CellAccessorWrapper::set_coarsen_flag(const bool coarsen_flag) + { + if ((dim == 2) && (spacedim == 2)) + internal::set_coarsen_flag<2,2>(coarsen_flag, cell_accessor); + else if ((dim == 2) && (spacedim == 3)) + internal::set_coarsen_flag<2,3>(coarsen_flag, cell_accessor); + else + internal::set_coarsen_flag<3,3>(coarsen_flag, cell_accessor); + } + + + + bool CellAccessorWrapper::get_coarsen_flag() + { + if ((dim == 2) && (spacedim == 2)) + return internal::get_coarsen_flag<2,2>(cell_accessor); + else if ((dim == 2) && (spacedim == 3)) + return internal::get_coarsen_flag<2,3>(cell_accessor); + else + return internal::get_coarsen_flag<3,3>(cell_accessor); + } + + + PointWrapper CellAccessorWrapper::get_barycenter() + { + if ((dim == 2) && (spacedim == 2)) + return internal::get_barycenter<2,2>(cell_accessor); + else if ((dim == 2) && (spacedim == 3)) + return internal::get_barycenter<2,3>(cell_accessor); + else + return internal::get_barycenter<3,3>(cell_accessor); + } + + + + void CellAccessorWrapper::set_material_id(const int material_id) + { + AssertThrow(material_id < dealii::numbers::invalid_material_id, + dealii::ExcMessage("material_id is too large.")); + if ((dim == 2) && (spacedim == 2)) + return internal::set_material_id<2,2>(material_id, cell_accessor); + else if ((dim == 2) && (spacedim == 3)) + return internal::set_material_id<2,3>(material_id, cell_accessor); + else + return internal::set_material_id<3,3>(material_id, cell_accessor); + } + + + + int CellAccessorWrapper::get_material_id() + { + if ((dim == 2) && (spacedim == 2)) + return internal::get_material_id<2,2>(cell_accessor); + else if ((dim == 2) && (spacedim == 3)) + return internal::get_material_id<2,3>(cell_accessor); + else + return internal::get_material_id<3,3>(cell_accessor); + } + + + + void CellAccessorWrapper::set_vertex(const int i, + PointWrapper &point_wrapper) + { + AssertThrow(i(i, point_wrapper, cell_accessor); + else if ((dim == 2) && (spacedim == 3)) + internal::set_vertex<2,3>(i, point_wrapper, cell_accessor); + else + internal::set_vertex<3,3>(i, point_wrapper, cell_accessor); + } + + + + PointWrapper CellAccessorWrapper::get_vertex(const int i) + { + AssertThrow(i(i, cell_accessor); + else if ((dim == 2) && (spacedim == 3)) + return internal::get_vertex<2,3>(i, cell_accessor); + else + return internal::get_vertex<3,3>(i, cell_accessor); + } + + + + void CellAccessorWrapper::set_manifold_id(const int manifold_id) + { + if ((dim == 2) && (spacedim == 2)) + internal::set_manifold_id<2,2>(manifold_id, cell_accessor); + else if ((dim == 2) && (spacedim == 3)) + internal::set_manifold_id<2,3>(manifold_id, cell_accessor); + else + internal::set_manifold_id<3,3>(manifold_id, cell_accessor); + } + + + + int CellAccessorWrapper::get_manifold_id() + { + if ((dim == 2) && (spacedim == 2)) + internal::get_manifold_id<2,2>(cell_accessor); + else if ((dim== 2) && (spacedim == 3)) + internal::get_manifold_id<2,3>(cell_accessor); + else + internal::get_manifold_id<3,3>(cell_accessor); + } +} diff --git a/source/python/export_cell_accessor.cc b/source/python/export_cell_accessor.cc new file mode 100644 index 0000000000..05c82d9ba4 --- /dev/null +++ b/source/python/export_cell_accessor.cc @@ -0,0 +1,49 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +#include +#include +#include + +namespace PyDealII +{ + void export_cell_accessor() + { + boost::python::class_("CellAccessor", + boost::python::init ()) + .add_property("refine_flag", &CellAccessorWrapper::get_refine_flag, + &CellAccessorWrapper::set_refine_flag, + "Get the refine_flag of the cell. The possibilities are in 2D: isotropic, no_refinement, cut_x, cut_y, and cut_xy and in 3D: isotropic, no_refinement, cut_x, cut_y, cut_z, cut_xy, cut_xz, cut_yz, and cut_xyz.") + .add_property("coarsen_flag", &CellAccessorWrapper::get_coarsen_flag, + &CellAccessorWrapper::set_coarsen_flag, + "Get the coarsen_flag of the cell.") + .add_property("material_id", &CellAccessorWrapper::get_material_id, + &CellAccessorWrapper::set_material_id, + "Get the material_id of the cell.") + .add_property("manifold_id", &CellAccessorWrapper::get_manifold_id, + &CellAccessorWrapper::set_manifold_id, + "Get the manifold_id of the cell.") + .def("barycenter", &CellAccessorWrapper::get_barycenter, + "Return the barycenter of the current cell.", + boost::python::args("self")) + .def("set_vertex", &CellAccessorWrapper::set_vertex, + "Set the ith vertex of the cell to point_wrapper", + boost::python::args("self", "i", "point_wrapper")) + .def("get_vertex", &CellAccessorWrapper::get_vertex, + "Get the ith vertex of the cell", + boost::python::args("self", "i")); + } +} + diff --git a/source/python/export_triangulation.cc b/source/python/export_triangulation.cc index 27628033d8..4b7917c2eb 100644 --- a/source/python/export_triangulation.cc +++ b/source/python/export_triangulation.cc @@ -15,6 +15,7 @@ #include #include +#include namespace PyDealII { @@ -33,7 +34,8 @@ namespace PyDealII void export_triangulation() { boost::python::class_("Triangulation", boost::python::init()) - .def(boost::python::init()) + .def(boost::python::init()) + .def("__iter__", boost::python::iterator()) .def("n_active_cells", &TriangulationWrapper::n_active_cells, "Return the number of active cells", diff --git a/source/python/point_wrapper.cc b/source/python/point_wrapper.cc index cf17461706..bdda73fae6 100644 --- a/source/python/point_wrapper.cc +++ b/source/python/point_wrapper.cc @@ -170,7 +170,7 @@ namespace PyDealII point = new dealii::Point<3> ((*other_point)[0], (*other_point)[1], (*other_point)[2]); } else - AssertThrow(false, + AssertThrow(false, dealii::ExcMessage("The dimension of the point should be 2 or 3.")); } } diff --git a/source/python/triangulation_wrapper.cc b/source/python/triangulation_wrapper.cc index 6f16f30067..d78cdea259 100644 --- a/source/python/triangulation_wrapper.cc +++ b/source/python/triangulation_wrapper.cc @@ -14,6 +14,7 @@ // --------------------------------------------------------------------- #include +#include #include #include #include @@ -166,11 +167,27 @@ namespace PyDealII // inside merge_triangulations. triangulation = tria; } + + + + template + void build_cells_vector(TriangulationWrapper &triangulation_wrapper, + std::vector &cells) + { + dealii::Triangulation *tria = + static_cast*>( + triangulation_wrapper.get_triangulation()); + for (auto cell : tria->active_cell_iterators()) + cells.push_back(CellAccessorWrapper(triangulation_wrapper, cell->level(), + cell->index())); + } } TriangulationWrapper::TriangulationWrapper(const std::string &dimension) + : + updated(false) { if ((dimension.compare("2D")==0) || (dimension.compare("2d")==0)) setup("2D", "2D"); @@ -182,8 +199,10 @@ namespace PyDealII - TriangulationWrapper::TriangulationWrapper(const std::string &dimension, + TriangulationWrapper::TriangulationWrapper(const std::string &dimension, const std::string &spacedimension) + : + updated(false) { setup(dimension, spacedimension); } @@ -193,31 +212,31 @@ namespace PyDealII TriangulationWrapper::~TriangulationWrapper() { if (triangulation != nullptr) - { - if (dim == 2) { - if (spacedim == 2) - { - // We cannot call delete on a void pointer so cast the void pointer back - // first. - dealii::Triangulation<2,2> *tmp = - static_cast*>(triangulation); - delete tmp; - } + if (dim == 2) + { + if (spacedim == 2) + { + // We cannot call delete on a void pointer so cast the void pointer back + // first. + dealii::Triangulation<2,2> *tmp = + static_cast*>(triangulation); + delete tmp; + } + else + { + dealii::Triangulation<2,3> *tmp = + static_cast*>(triangulation); + delete tmp; + } + } else - { - dealii::Triangulation<2,3> *tmp = - static_cast*>(triangulation); - delete tmp; - } + { + dealii::Triangulation<3,3> *tmp = static_cast*>(triangulation); + delete tmp; + } + triangulation = nullptr; } - else - { - dealii::Triangulation<3,3> *tmp = static_cast*>(triangulation); - delete tmp; - } - triangulation = nullptr; - } dim = -1; } @@ -399,6 +418,30 @@ namespace PyDealII + TriangulationWrapper::iterator TriangulationWrapper::begin() + { + // 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(); + } + + + void TriangulationWrapper::save(const std::string &filename) const { std::ofstream ofs(filename); @@ -412,14 +455,14 @@ namespace PyDealII oa << *tria; } else if ((dim == 2) && (spacedim == 3)) - { { - dealii::Triangulation<2,3> *tria = - static_cast*>(triangulation); + { + dealii::Triangulation<2,3> *tria = + static_cast*>(triangulation); - oa << *tria; + oa << *tria; + } } - } else { dealii::Triangulation<3,3> *tria = @@ -467,26 +510,26 @@ namespace PyDealII if ((dimension.compare("2D") == 0) || (dimension.compare("2d") == 0)) { dim = 2; - - if ((spacedimension.compare("2D") == 0) || + + if ((spacedimension.compare("2D") == 0) || (spacedimension.compare("2d") == 0)) - { - spacedim = 2; - triangulation = new dealii::Triangulation<2,2>(); - } - else if ((spacedimension.compare("3D") == 0) || + { + spacedim = 2; + triangulation = new dealii::Triangulation<2,2>(); + } + else if ((spacedimension.compare("3D") == 0) || (spacedimension.compare("3d") == 0)) - { - spacedim = 3; - triangulation = new dealii::Triangulation<2,3>(); - } + { + spacedim = 3; + triangulation = new dealii::Triangulation<2,3>(); + } else AssertThrow(false, dealii::ExcMessage("Spacedimension needs to be 2D or 3D.")); } else if ((dimension.compare("3D") == 0) || (dimension.compare("3d") == 0)) { - if ((spacedimension.compare("3D") != 0) && - (spacedimension.compare("3d") != 0)) + if ((spacedimension.compare("3D") != 0) && + (spacedimension.compare("3d") != 0)) AssertThrow(false, dealii::ExcMessage("Spacedimension needs to be 3D.")); dim = 3; spacedim = 3; @@ -495,4 +538,23 @@ 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; + } } diff --git a/source/python/wrappers.cc b/source/python/wrappers.cc index 3c0b6c3d39..bd3ee4144e 100644 --- a/source/python/wrappers.cc +++ b/source/python/wrappers.cc @@ -17,6 +17,7 @@ namespace PyDealII { + void export_cell_accessor(); void export_point(); void export_triangulation(); } @@ -37,6 +38,7 @@ BOOST_PYTHON_MODULE(PyDealII) doc_options.enable_py_signatures(); doc_options.disable_cpp_signatures(); + PyDealII::export_cell_accessor(); PyDealII::export_point(); PyDealII::export_triangulation(); } diff --git a/tests/python/cell_accessor_wrapper.py b/tests/python/cell_accessor_wrapper.py new file mode 100644 index 0000000000..07792aa311 --- /dev/null +++ b/tests/python/cell_accessor_wrapper.py @@ -0,0 +1,93 @@ +# --------------------------------------------------------------------- +# +# Copyright (C) 2016 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 at +# the top level of the deal.II distribution. +# +# --------------------------------------------------------------------- + +import unittest +from pydealii.PyDealII import * + + +class TestCellAccessorWrapper(unittest.TestCase): + + def setUp(self): + self.triangulation = Triangulation('2D') + self.triangulation.generate_hyper_cube() + self.triangulation.refine_global(1) + + def test_material_id(self): + material_id = 0 + for cell in self.triangulation: + cell.material_id = material_id + material_id += 1 + material_id = 0 + for cell in self.triangulation: + self.assertEqual(cell.material_id, material_id) + material_id += 1 + + def test_manifold_id(self): + manifold_id = 0 + for cell in self.triangulation: + cell.manifold_id = manifold_id + manifold_id += 1 + manifold_id = 0 + for cell in self.triangulation: + 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: + cell.refine_flag = refine_flags[index] + index += 1 + index = 0 + for cell in self.triangulation: + self.assertEqual(cell.refine_flag, refine_flags[index]) + index += 1 + + def test_coarsen_flag(self): + coarsen_flag = True + for cell in self.triangulation: + cell.coarsen_flag = coarsen_flag + coarsen_flag = not coarsen_flag + coarsen_flag = True + for cell in self.triangulation: + 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: + barycenter = cell.barycenter() + self.assertEqual(barycenter.x, centers[index][0]) + self.assertEqual(barycenter.y, centers[index][1]) + index += 1 + + def test_move_vertex(self): + point = Point([0.6, 0.6]) + for cell in self.triangulation: + cell.set_vertex(3, point) + vertex = cell.get_vertex(3) + break + vertices = [3, 2, 1, 0] + index = 0 + for cell in self.triangulation: + vertex = cell.get_vertex(vertices[index]) + self.assertEqual(vertex.x, point.x) + self.assertEqual(vertex.y, point.y) + index += 1 + + +if __name__ == '__main__': + unittest.main() -- 2.39.5