--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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 <deal.II/grid/tria_accessor.h>
+
+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
#define dealii__point_wrapper_h
#include <boost/python.hpp>
+
namespace PyDealII
{
class PointWrapper
namespace PyDealII
{
+ class CellAccessorWrapper;
+
class TriangulationWrapper
{
public:
+ typedef std::vector<CellAccessorWrapper>::iterator iterator;
+
/**
* Constructor. Takes a string @p dim with one of the following values
* "2D", "2d", "3D", or "3d".
*/
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.
*/
* Pointer that can be casted to the underlying Triangulation object.
*/
void *triangulation;
+
+ /**
+ * Vector of all the cells of the Triangulation.
+ */
+ std::vector<CellAccessorWrapper> cells;
};
SET(_src
wrappers.cc
+ export_cell_accessor.cc
export_point.cc
export_triangulation.cc
+ cell_accessor_wrapper.cc
point_wrapper.cc
triangulation_wrapper.cc
)
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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 <boost/python.hpp>
+#include <deal.II/python/point_wrapper.h>
+#include <deal.II/python/triangulation_wrapper.h>
+#include <deal.II/python/cell_accessor_wrapper.h>
+
+namespace PyDealII
+{
+ namespace internal
+ {
+ template <int dim, int spacedim>
+ void set_refine_flag(const std::string &refinement_case,
+ void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(cell_accessor);
+
+ std::unique_ptr<dealii::RefinementCase<dim>> ref_case;
+ if (refinement_case.compare("isotropic") == 0)
+ ref_case.reset(new dealii::RefinementCase<dim>(
+ dealii::RefinementPossibilities<dim>::Possibilities::isotropic_refinement));
+ else if (refinement_case.compare("no_refinement") == 0)
+ ref_case.reset(new dealii::RefinementCase<dim>(
+ dealii::RefinementPossibilities<dim>::Possibilities::no_refinement));
+ else if (refinement_case.compare("cut_x") == 0)
+ ref_case.reset(new dealii::RefinementCase<dim>(
+ dealii::RefinementPossibilities<dim>::Possibilities::cut_x));
+ else if (refinement_case.compare("cut_y") == 0)
+ ref_case.reset(new dealii::RefinementCase<dim>(
+ dealii::RefinementPossibilities<dim>::Possibilities::cut_y));
+ else if (refinement_case.compare("cut_xy") == 0)
+ ref_case.reset(new dealii::RefinementCase<dim>(
+ dealii::RefinementPossibilities<dim>::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 <int dim, int spacedim>
+ std::string get_refine_flag(void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(cell_accessor);
+
+ std::string refine_flag;
+ dealii::RefinementCase<dim> ref_case = cell->refine_flag_set();
+ switch (static_cast<int>(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 <int dim, int spacedim>
+ void set_coarsen_flag(const bool coarsen_flag, void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(cell_accessor);
+ if (coarsen_flag == true)
+ cell->set_coarsen_flag();
+ else
+ cell->clear_coarsen_flag();
+ }
+
+
+
+ template <int dim, int spacedim>
+ bool get_coarsen_flag(void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(cell_accessor);
+
+ return cell->coarsen_flag_set();
+ }
+
+
+
+ template <int dim, int spacedim>
+ PointWrapper get_barycenter(void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(cell_accessor);
+ dealii::Point<spacedim> barycenter = cell->barycenter();
+ boost::python::list barycenter_list;
+ for (int i=0; i<dim; ++i)
+ barycenter_list.append(barycenter[i]);
+
+ return PointWrapper(barycenter_list);
+ }
+
+
+
+ template <int dim, int spacedim>
+ void set_material_id(const int material_id, void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(cell_accessor);
+ cell->set_material_id(material_id);
+ }
+
+
+
+ template <int dim, int spacedim>
+ int get_material_id(void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(cell_accessor);
+
+ return cell->material_id();
+ }
+
+
+
+ template <int dim, int spacedim>
+ void set_vertex(const int i,
+ PointWrapper &point_wrapper,
+ void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(cell_accessor);
+ dealii::Point<spacedim> *point =
+ static_cast<dealii::Point<spacedim>*>(point_wrapper.get_point());
+
+ cell->vertex(i) = *point;
+ }
+
+
+
+ template <int dim, int spacedim>
+ PointWrapper get_vertex(const int i, void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(cell_accessor);
+ dealii::Point<spacedim> vertex = cell->vertex(i);
+
+ boost::python::list coordinates;
+ for (int i=0; i<spacedim; ++i)
+ coordinates.append(vertex[i]);
+
+ return PointWrapper(coordinates);
+ }
+
+
+
+ template <int dim, int spacedim>
+ void set_manifold_id(const int manifold_id, void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(cell_accessor);
+ cell->set_manifold_id(manifold_id);
+ }
+
+
+
+ template <int dim, int spacedim>
+ int get_manifold_id(void *cell_accessor)
+ {
+ dealii::CellAccessor<dim,spacedim> *cell =
+ static_cast<dealii::CellAccessor<dim,spacedim>*>(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<dealii::CellAccessor<2,2>*>(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<dealii::CellAccessor<2,3>*>(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<dealii::CellAccessor<3,3>*>(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<dealii::Triangulation<2,2>*> (
+ 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<dealii::Triangulation<2,3>*> (
+ 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<dealii::Triangulation<3,3>*> (
+ 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<dealii::CellAccessor<2,2>*>(cell_accessor);
+ delete tmp;
+ }
+ else if ((dim == 2) && (spacedim == 3))
+ {
+ dealii::CellAccessor<3,3> *tmp =
+ static_cast<dealii::CellAccessor<3,3>*>(cell_accessor);
+ delete tmp;
+ }
+ else
+ {
+ dealii::CellAccessor<3,3> *tmp =
+ static_cast<dealii::CellAccessor<3,3>*>(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<std::pow(2,dim),
+ dealii::ExcMessage("The vertex does not exist."));
+ if ((dim == 2) && (spacedim == 2))
+ internal::set_vertex<2,2>(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<std::pow(2,dim),
+ dealii::ExcMessage("The vertex does not exist."));
+ if ((dim == 2) && (spacedim == 2))
+ return internal::get_vertex<2,2>(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);
+ }
+}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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 <boost/python.hpp>
+#include <deal.II/python/cell_accessor_wrapper.h>
+#include <deal.II/python/triangulation_wrapper.h>
+
+namespace PyDealII
+{
+ void export_cell_accessor()
+ {
+ boost::python::class_<CellAccessorWrapper>("CellAccessor",
+ boost::python::init<TriangulationWrapper &, const int, const int> ())
+ .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"));
+ }
+}
+
#include <boost/python.hpp>
#include <deal.II/python/triangulation_wrapper.h>
+#include <deal.II/python/cell_accessor_wrapper.h>
namespace PyDealII
{
void export_triangulation()
{
boost::python::class_<TriangulationWrapper>("Triangulation", boost::python::init<const std::string &>())
- .def(boost::python::init<const std::string&, 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",
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."));
}
}
// ---------------------------------------------------------------------
#include <deal.II/python/triangulation_wrapper.h>
+#include <deal.II/python/cell_accessor_wrapper.h>
#include <deal.II/grid/tria.h>
#include <deal.II/grid/grid_generator.h>
#include <deal.II/grid/grid_tools.h>
// inside merge_triangulations.
triangulation = tria;
}
+
+
+
+ template <int dim, int spacedim>
+ void build_cells_vector(TriangulationWrapper &triangulation_wrapper,
+ std::vector<CellAccessorWrapper> &cells)
+ {
+ dealii::Triangulation<dim,spacedim> *tria =
+ static_cast<dealii::Triangulation<dim,spacedim>*>(
+ 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");
- TriangulationWrapper::TriangulationWrapper(const std::string &dimension,
+ TriangulationWrapper::TriangulationWrapper(const std::string &dimension,
const std::string &spacedimension)
+ :
+ updated(false)
{
setup(dimension, spacedimension);
}
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<dealii::Triangulation<2,2>*>(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<dealii::Triangulation<2,2>*>(triangulation);
+ delete tmp;
+ }
+ else
+ {
+ dealii::Triangulation<2,3> *tmp =
+ static_cast<dealii::Triangulation<2,3>*>(triangulation);
+ delete tmp;
+ }
+ }
else
- {
- dealii::Triangulation<2,3> *tmp =
- static_cast<dealii::Triangulation<2,3>*>(triangulation);
- delete tmp;
- }
+ {
+ dealii::Triangulation<3,3> *tmp = static_cast<dealii::Triangulation<3,3>*>(triangulation);
+ delete tmp;
+ }
+ triangulation = nullptr;
}
- else
- {
- dealii::Triangulation<3,3> *tmp = static_cast<dealii::Triangulation<3,3>*>(triangulation);
- delete tmp;
- }
- triangulation = nullptr;
- }
dim = -1;
}
+ 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);
oa << *tria;
}
else if ((dim == 2) && (spacedim == 3))
- {
{
- dealii::Triangulation<2,3> *tria =
- static_cast<dealii::Triangulation<2,3>*>(triangulation);
+ {
+ dealii::Triangulation<2,3> *tria =
+ static_cast<dealii::Triangulation<2,3>*>(triangulation);
- oa << *tria;
+ oa << *tria;
+ }
}
- }
else
{
dealii::Triangulation<3,3> *tria =
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;
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;
+ }
}
namespace PyDealII
{
+ void export_cell_accessor();
void export_point();
void export_triangulation();
}
doc_options.enable_py_signatures();
doc_options.disable_cpp_signatures();
+ PyDealII::export_cell_accessor();
PyDealII::export_point();
PyDealII::export_triangulation();
}
--- /dev/null
+# ---------------------------------------------------------------------
+#
+# 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()