class PointWrapper
{
public:
+ /**
+ * This constructor exists only so that a std::vector of PointWrapper can be
+ * created. This constructor should not be used.
+ */
+ PointWrapper();
+
+ /**
+ * Constructor. Use a list of double to initialize the Point.
+ */
PointWrapper(boost::python::list list);
+
+ /**
+ * Destructpr.
+ */
+ ~PointWrapper();
+
+ /**
+ * Return the first component of the Point.
+ */
double get_x();
+
+ /**
+ * Set the first component of the Point.
+ */
void set_x(double x);
+
+ /**
+ * Return the second component of the Point.
+ */
double get_y();
+
+ /**
+ * Set the second component of the Point.
+ */
void set_y(double y);
+
+ /**
+ * Return the third component of the Point.
+ */
double get_z();
+
+ /**
+ * Set the third component of the Point.
+ */
void set_z(double z);
+ /**
+ * Return a pointer that can be casted to the underlying Point.
+ */
+ void *get_point();
+
private:
+ /**
+ * Dimension of the Point.
+ */
int dim;
- std::shared_ptr<void> point;
+
+ /**
+ * Pointer to the underlying Point object.
+ */
+ void *point;
};
+
+
+ //--------------------- Inline functions ----------------------//
+
+
+
+ inline
+ void *PointWrapper::get_point()
+ {
+ return point;
+ }
}
#endif
--- /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__triangulation_wrapper_h
+#define dealii__triangulation_wrapper_h
+
+#include <boost/python.hpp>
+#include <deal.II/python/point_wrapper.h>
+#include <string>
+
+namespace PyDealII
+{
+ class TriangulationWrapper
+ {
+ public:
+ /**
+ * Constructor. Takes a string @p dim with one of the following values
+ * "2D", "2d", "3D", or "3d".
+ */
+ TriangulationWrapper(const std::string &dim);
+
+ /**
+ * Destructor.
+ */
+ ~TriangulationWrapper();
+
+ /**
+ * Return the number of active cells.
+ */
+ unsigned int n_active_cells();
+
+ /**
+ * Generate a hyper cube (square in 2D and cube in 3D) with exactly one
+ * cell.
+ */
+ void generate_hyper_cube(const double left = 0.,
+ const double right = 1.,
+ const bool colorize = false);
+
+ /**
+ * Generate a simplex with (dim+1) vertices and mesh cells.
+ */
+ void generate_simplex(boost::python::list &vertices);
+
+ /**
+ * Same as hyper_cube but not only one cell is created but each coordinate
+ * direction is subdivided in @p repetitions cells.
+ */
+ void generate_subdivided_hyper_cube(const unsigned int repetitions,
+ const double left = 0.,
+ const double right = 1.);
+
+ /**
+ * Generate a coordinate-parallel brick from the two diagonally opposite
+ * corners points @p1 and @p2.
+ */
+ void generate_hyper_rectangle(PointWrapper &p1,
+ PointWrapper &p2,
+ const bool colorize = false);
+
+ /**
+ * Generate a coordinate-parallel brick from the two diagonally opposite
+ * corners points @p1 and @p2. In direction i, repetitions[i] cells are
+ * created.
+ */
+ void generate_subdivided_hyper_rectangle(boost::python::list &repetitions,
+ PointWrapper &p1,
+ PointWrapper &p2,
+ const bool colorize = false);
+
+ /**
+ * Generate a hyperball, i.e. a circle or a ball around @p center with
+ * given @p radius.
+ */
+ void generate_hyper_ball(PointWrapper ¢er,
+ const double radius = 1.);
+
+ /**
+ * Shift each vertex of the Triangulation by the given @p shift_list.
+ */
+ void shift(boost::python::list &shift_list);
+
+ /**
+ * Given two triangulations, create the triangulation that contains the
+ * cells of both triangulations.
+ */
+ void merge_triangulations(TriangulationWrapper &triangulation_1,
+ TriangulationWrapper &triangulation_2);
+
+ /**
+ * Refine all the cells @p times times.
+ */
+ void refine_global(const unsigned int times);
+
+ /**
+ * Write the Triangulation in file.
+ */
+ void save(const std::string &filename);
+
+ /**
+ * Load the Triangulation from a file.
+ */
+ void load(const std::string &filename);
+
+ /**
+ * Return the dimension of the underlying Triangulation object.
+ */
+ int get_dim();
+
+ /**
+ * Return a pointer that can be casted to the underlying Triangulation
+ * object.
+ */
+ void *get_triangulation();
+
+ private:
+ /**
+ * Reset the underlying Triangulation object.
+ */
+ void reset_triangulation();
+
+ /**
+ * Dimension of the underlying Triangulation object.
+ */
+ int dim;
+
+ /**
+ * Pointer that can be casted to the underlying Triangulation object.
+ */
+ void *triangulation;
+ };
+
+
+//-------------------- inline functions -----------------------//
+
+
+
+ inline
+ int TriangulationWrapper::get_dim()
+ {
+ return dim;
+ }
+
+
+
+ inline
+ void *TriangulationWrapper::get_triangulation()
+ {
+ return triangulation;
+ }
+}
+
+#endif
export_point.cc
export_triangulation.cc
point_wrapper.cc
+ triangulation_wrapper.cc
)
FILE(GLOB _header
{
void export_point()
{
- boost::python::class_<PointWrapper> ("Point",boost::python::init<boost::python::list>())
+ boost::python::class_<PointWrapper> ("Point", boost::python::init<boost::python::list>())
.add_property("x", &PointWrapper::get_x, &PointWrapper::set_x, "Get the x component of the point.")
.add_property("y", &PointWrapper::get_y, &PointWrapper::set_y, "Get the y component of the point.")
.add_property("z", &PointWrapper::get_z, &PointWrapper::set_z, "Get the z component of the point.");
//
// ---------------------------------------------------------------------
-#include <boost/archive/text_oarchive.hpp>
-#include <boost/archive/text_iarchive.hpp>
#include <boost/python.hpp>
-#include <deal.II/grid/tria.h>
-#include <deal.II/grid/grid_generator.h>
-#include <fstream>
+#include <deal.II/python/triangulation_wrapper.h>
namespace PyDealII
{
- template <int dim>
- unsigned int n_active_cells(const dealii::Triangulation<dim> &triangulation)
- {
- return triangulation.n_active_cells();
- }
-
-
-
- template <int dim>
- void generate_hyper_cube(dealii::Triangulation<dim> &triangulation,
- const double left = 0.,
- const double right = 0.,
- const bool colorize = false)
- {
- dealii::GridGenerator::hyper_cube(triangulation, left, right, colorize);
- }
-
-
-
- void generate_simplex(dealii::Triangulation<2> &triangulation,
- dealii::Point<2> vertex_1,
- dealii::Point<2> vertex_2,
- dealii::Point<2> vertex_3)
- {
- std::vector<dealii::Point<2>> vertices(3);
- vertices[0] = vertex_1;
- vertices[1] = vertex_2;
- vertices[2] = vertex_3;
-
- dealii::GridGenerator::simplex(triangulation, vertices);
- }
-
-
-
- template <int dim>
- void generate_subdivided_hyper_cube(dealii::Triangulation<dim> &triangulation,
- const unsigned int repetitions,
- const double left = 0.,
- const double right = 1.)
- {
- dealii::GridGenerator::subdivided_hyper_cube(triangulation, repetitions,
- left, right);
- }
-
-
-
- template <int dim>
- void generate_hyper_rectangle(dealii::Triangulation<dim> &triangulation,
- const dealii::Point<dim> &p1,
- const dealii::Point<dim> &p2,
- const bool colorize = false)
- {
- dealii::GridGenerator::hyper_rectangle(triangulation, p1, p2, colorize);
- }
-
-
-
- void generate_subdivided_hyper_rectangle(dealii::Triangulation<2> &triangulation,
- const unsigned int repetition_x,
- const unsigned int repetition_y,
- const dealii::Point<2> &p1,
- const dealii::Point<2> &p2,
- const bool colorize = false)
- {
- std::vector<unsigned int> repetitions(2);
- repetitions[0] = repetition_x;
- repetitions[1] = repetition_y;
- dealii::GridGenerator::subdivided_hyper_rectangle(triangulation, repetitions,
- p1, p2, colorize);
- }
-
-
-
- template <int dim>
- void generate_hyper_ball(dealii::Triangulation<dim> &triangulation,
- const dealii::Point<dim> ¢er = dealii::Point<dim>(),
- const double radius = 1.)
- {
- dealii::GridGenerator::hyper_ball(triangulation, center, radius);
- }
-
-
-
- template <int dim>
- void save(const dealii::Triangulation<dim> &triangulation,
- const std::string filename)
- {
- std::ofstream ofs(filename);
- boost::archive::text_oarchive oa(ofs);
- oa << triangulation;
- }
-
-
-
- template <int dim>
- void load(dealii::Triangulation<dim> &triangulation,
- const std::string filename)
- {
- std::ifstream ifs(filename);
- boost::archive::text_iarchive ia(ifs);
- ia >> triangulation;
- }
-
-
-
// Macro to enable default arguments
- BOOST_PYTHON_FUNCTION_OVERLOADS(generate_hyper_cube_overloads, generate_hyper_cube, 1, 4)
- BOOST_PYTHON_FUNCTION_OVERLOADS(generate_subdivided_hyper_cube_overloads,
- generate_subdivided_hyper_cube, 2, 4)
- BOOST_PYTHON_FUNCTION_OVERLOADS(generate_hyper_rectangle_overloads,
- generate_hyper_rectangle, 3, 4)
- BOOST_PYTHON_FUNCTION_OVERLOADS(generate_subdivided_hyper_rectangle_overloads,
- generate_subdivided_hyper_rectangle, 5, 6)
- BOOST_PYTHON_FUNCTION_OVERLOADS(generate_hyper_ball_overloads, generate_hyper_ball, 1, 3)
+ BOOST_PYTHON_MEMBER_FUNCTION_OVERLOADS(generate_hyper_cube_overloads, generate_hyper_cube, 0, 3)
+ BOOST_PYTHON_MEMBER_FUNCTION_OVERLOADS(generate_subdivided_hyper_cube_overloads,
+ generate_subdivided_hyper_cube, 1, 3)
+ BOOST_PYTHON_MEMBER_FUNCTION_OVERLOADS(generate_hyper_rectangle_overloads,
+ generate_hyper_rectangle, 2, 3)
+ BOOST_PYTHON_MEMBER_FUNCTION_OVERLOADS(generate_subdivided_hyper_rectangle_overloads,
+ generate_subdivided_hyper_rectangle, 3, 4)
+ BOOST_PYTHON_MEMBER_FUNCTION_OVERLOADS(generate_hyper_ball_overloads, generate_hyper_ball, 1, 2)
void export_triangulation()
{
- boost::python::class_<dealii::Triangulation<2>> ("Triangulation")
- .def("n_active_cells", n_active_cells<2>,
- "Return the number of active cells",
- boost::python::args("self"))
- .def("generate_hyper_cube", generate_hyper_cube<2>,
- generate_hyper_cube_overloads(
- boost::python::args("self", "left", "right", "colorize"),
- "Generate a hyper_cube."))
- .def("generate_simplex", generate_simplex,
- "Generate a simplex.",
- boost::python::args("self", "vertex_1", "vertex_2", "vertex_3"))
- .def("generate_subdivided_hyper_cube", generate_subdivided_hyper_cube<2>,
- generate_subdivided_hyper_cube_overloads(
- boost::python::args("self", "repetitions", "left", "right"),
- "Generate a subdivided hyper_cube."))
- .def("generate_hyper_rectangle", generate_hyper_rectangle<2>,
- generate_hyper_rectangle_overloads(
- boost::python::args("self", "p1", "p2", "colorize"),
- "Generate a hyper_rectangle."))
- .def("generate_subdivided_hyper_rectangle",
- generate_subdivided_hyper_rectangle,
- generate_subdivided_hyper_rectangle_overloads(
- boost::python::args("self", "repetition_x",
- "repetition_y","p1", "p2",
- "colorize"),
- "Generate a subdivided hyper_rectangle."))
- .def("generate_hyper_ball",
- generate_hyper_ball<2>,
- generate_hyper_ball_overloads(
- boost::python::args("self", "center", "radius"),
- "Generate a hyper_ball."))
- .def("refine_global", &dealii::Triangulation<2>::refine_global,
- "Refine the mesh uniformly.",
- boost::python::args("self", "times"))
- .def("save", save<2>, "Serialize and save the triangulation.",
- boost::python::args("self", "filename"))
- .def("load", load<2>, "Load and deserialize a triangulation.",
- boost::python::args("self", "filename"));
+ boost::python::class_<TriangulationWrapper>("Triangulation", boost::python::init<const std::string &>())
+ .def("n_active_cells",
+ &TriangulationWrapper::n_active_cells,
+ "Return the number of active cells",
+ boost::python::args("self"))
+ .def("generate_hyper_cube",
+ &TriangulationWrapper::generate_hyper_cube,
+ generate_hyper_cube_overloads(
+ boost::python::args("self", "left", "right", "colorize"),
+ "Generate a hyper_cube (square in 2D and cube in 3D) with exactly one cell."))
+ .def("generate_simplex",
+ &TriangulationWrapper::generate_simplex,
+ "Generate a simplex with (dim+1) vertices and mesh cells.",
+ boost::python::args("self", "vertices"))
+ .def("generate_subdivided_hyper_cube",
+ &TriangulationWrapper::generate_subdivided_hyper_cube,
+ generate_subdivided_hyper_cube_overloads(
+ boost::python::args("self", "repetitions", "left", "right"),
+ "Same as hyper_cube but not only one cells is created but each coordinate direction is subdivdided in repetitions cells."))
+ .def("generate_hyper_rectangle",
+ &TriangulationWrapper::generate_hyper_rectangle,
+ generate_hyper_rectangle_overloads(
+ boost::python::args("self", "p1", "p2", "colorize"),
+ "Generate a coordinate-parallel brick from the two diagonally opposite corners points p1 and p2."))
+ .def("generate_subdivided_hyper_rectangle",
+ &TriangulationWrapper::generate_subdivided_hyper_rectangle,
+ generate_subdivided_hyper_rectangle_overloads(
+ boost::python::args("self", "repetitions",
+ "p1", "p2", "colorize"),
+ "Generate a coordinate-parallel brick from the two diagonally opposite corners point @p1 and @p2. In direction i, repetitions[i] cells are created."))
+ .def("generate_hyper_ball",
+ &TriangulationWrapper::generate_hyper_ball,
+ generate_hyper_ball_overloads(
+ boost::python::args("self", "center", "radius"),
+ "Generate a hyperball, i.e. a circle or a ball around @p center with given @p radius."))
+ .def("shift",
+ &TriangulationWrapper::shift,
+ "Shift every vertex of the Triangulation by the given shift vector.",
+ boost::python::args("self", "shift"))
+ .def("merge_triangulations",
+ &TriangulationWrapper::merge_triangulations,
+ "Given two triangulations, create the triangulation that contains the cells of both triangulations.",
+ boost::python::args("self", "triangulation_1", "triangulation_2"))
+ .def("refine_global",
+ &TriangulationWrapper::refine_global,
+ "Refine all the cells times times.",
+ boost::python::args("self", "times"))
+ .def("save",
+ &TriangulationWrapper::save,
+ "Write the Triangulation to a file.",
+ boost::python::args("self", "filename"))
+ .def("load",
+ &TriangulationWrapper::load,
+ "Load the Triangulation from a file.",
+ boost::python::args("self", "filename"));
}
-
}
namespace PyDealII
{
+ PointWrapper::PointWrapper()
+ :
+ dim(-1),
+ point(nullptr)
+ {}
+
+
PointWrapper::PointWrapper(boost::python::list coord)
{
dim = boost::python::len(coord);
if (dim == 2)
- point.reset(new dealii::Point<2> (boost::python::extract<double>(coord[0]),
- boost::python::extract<double>(coord[1])));
+ point = new dealii::Point<2> (boost::python::extract<double>(coord[0]),
+ boost::python::extract<double>(coord[1]));
else if (dim == 3)
- point.reset(new dealii::Point<3> (boost::python::extract<double>(coord[0]),
- boost::python::extract<double>(coord[1]),
- boost::python::extract<double>(coord[2])));
+ point = new dealii::Point<3> (boost::python::extract<double>(coord[0]),
+ boost::python::extract<double>(coord[1]),
+ boost::python::extract<double>(coord[2]));
else
AssertThrow(false,
dealii::ExcMessage("The list of coordinates must contain two or three elements."));
+ PointWrapper::~PointWrapper()
+ {
+ if (point != nullptr)
+ {
+ if (dim == 2)
+ {
+ // We cannot call delete on a void pointer so cast the void pointer back
+ // first.
+ dealii::Point<2> *tmp = static_cast<dealii::Point<2>*>(point);
+ delete tmp;
+ }
+ else
+ {
+ dealii::Point<3> *tmp = static_cast<dealii::Point<3>*>(point);
+ delete tmp;
+ }
+ point = nullptr;
+ }
+ dim = -1;
+ }
+
+
+
double PointWrapper::get_x()
{
if (dim == 2)
- return (*std::static_pointer_cast<dealii::Point<2>>(point))(0);
+ return (*static_cast<dealii::Point<2>*>(point))(0);
else
- return (*std::static_pointer_cast<dealii::Point<3>>(point))(0);
+ return (*static_cast<dealii::Point<3>*>(point))(0);
}
void PointWrapper::set_x(double x)
{
if (dim == 2)
- (*std::static_pointer_cast<dealii::Point<2>>(point))(0) = x;
+ (*static_cast<dealii::Point<2>*>(point))(0) = x;
else
- (*std::static_pointer_cast<dealii::Point<3>>(point))(0) = x;
+ (*static_cast<dealii::Point<3>*>(point))(0) = x;
}
double PointWrapper::get_y()
{
if (dim == 2)
- return (*std::static_pointer_cast<dealii::Point<2>>(point))(1);
+ return (*static_cast<dealii::Point<2>*>(point))(1);
else
- return (*std::static_pointer_cast<dealii::Point<3>>(point))(1);
+ return (*static_cast<dealii::Point<3>*>(point))(1);
}
void PointWrapper::set_y(double y)
{
if (dim == 2)
- (*std::static_pointer_cast<dealii::Point<2>>(point))(1) = y;
+ (*static_cast<dealii::Point<2>*>(point))(1) = y;
else
- (*std::static_pointer_cast<dealii::Point<3>>(point))(1) = y;
+ (*static_cast<dealii::Point<3>*>(point))(1) = y;
}
double PointWrapper::get_z()
{
if (dim == 3)
- return (*std::static_pointer_cast<dealii::Point<3>>(point))(2);
+ return (*static_cast<dealii::Point<3>*>(point))(2);
else
AssertThrow(false,
dealii::ExcMessage("The z coordinate is only available for three-dimensional points"));
void PointWrapper::set_z(double z)
{
if (dim == 3)
- (*std::static_pointer_cast<dealii::Point<3>>(point))(2) = z;
+ (*static_cast<dealii::Point<3>*>(point))(2) = z;
else
AssertThrow(false,
dealii::ExcMessage("The z coordinate is only available for three-dimensional points"));
--- /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 <deal.II/python/triangulation_wrapper.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <boost/archive/text_oarchive.hpp>
+#include <boost/archive/text_iarchive.hpp>
+#include <fstream>
+
+namespace PyDealII
+{
+ TriangulationWrapper::TriangulationWrapper(const std::string &dimension)
+ :
+ triangulation(nullptr)
+ {
+ if ((dimension.compare("2D")==0) || (dimension.compare("2d")==0))
+ {
+ dim = 2;
+ }
+ else if ((dimension.compare("3D")==0) || (dimension.compare("3d")==0))
+ {
+ dim = 3;
+ }
+ else
+ AssertThrow(false, dealii::ExcMessage("Dimension needs to be 2D or 3D"));
+ }
+
+
+
+ TriangulationWrapper::~TriangulationWrapper()
+ {
+ if (triangulation != nullptr)
+ {
+ if (dim == 2)
+ {
+ // We cannot call delete on a void pointer so cast the void pointer back
+ // first.
+ dealii::Triangulation<2> *tmp = static_cast<dealii::Triangulation<2>*>(triangulation);
+ delete tmp;
+ }
+ else
+ {
+ dealii::Triangulation<3> *tmp = static_cast<dealii::Triangulation<3>*>(triangulation);
+ delete tmp;
+ }
+ triangulation = nullptr;
+ }
+ dim = -1;
+ }
+
+
+
+ unsigned int TriangulationWrapper::n_active_cells()
+ {
+ if (dim == 2)
+ {
+ return (*static_cast<dealii::Triangulation<2>*>(triangulation)).n_active_cells();
+ }
+ else
+ return (*static_cast<dealii::Triangulation<3>*>(triangulation)).n_active_cells();
+ }
+
+
+
+ void TriangulationWrapper::generate_hyper_cube(const double left,
+ const double right,
+ const bool colorize)
+ {
+ // Create a new Triangulation
+ reset_triangulation();
+
+ if (dim == 2)
+ {
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+ dealii::GridGenerator::hyper_cube(*tria, left, right, colorize);
+ }
+ else
+ {
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(triangulation);
+ dealii::GridGenerator::hyper_cube(*tria, left, right, colorize);
+ }
+ }
+
+
+
+ void TriangulationWrapper::generate_simplex(boost::python::list &vertices)
+ {
+ // Create a new Triangulation
+ reset_triangulation();
+
+ // Extract the PointWrapper object from the python list
+ std::vector<PointWrapper> wrapped_points(dim);
+ for (int i=0; i<dim; ++i)
+ wrapped_points[i] = boost::python::extract<PointWrapper>(vertices[i]);
+
+ if (dim == 2)
+ {
+ // Cast the PointWrapper objects to Point<dim>
+ std::vector<dealii::Point<2>> points(dim);
+ for (int i=0; i<dim; ++i)
+ points[i] = *(static_cast<dealii::Point<2>*>((wrapped_points[i]).get_point()));
+
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+ dealii::GridGenerator::simplex(*tria, points);
+ }
+ else
+ {
+ std::vector<dealii::Point<3>> points(dim);
+ for (int i=0; i<dim; ++i)
+ points[i] = *(static_cast<dealii::Point<3>*>((wrapped_points[i]).get_point()));
+
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(triangulation);
+ dealii::GridGenerator::simplex(*tria, points);
+ }
+ }
+
+
+
+ void TriangulationWrapper::generate_subdivided_hyper_cube(const unsigned int repetitions,
+ const double left,
+ const double right)
+ {
+ // Create a new Triangulation
+ reset_triangulation();
+
+ if (dim == 2)
+ {
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+ dealii::GridGenerator::subdivided_hyper_cube(*tria, repetitions, left, right);
+ }
+ else
+ {
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(triangulation);
+ dealii::GridGenerator::subdivided_hyper_cube(*tria, repetitions, left, right);
+ }
+ }
+
+
+
+ void TriangulationWrapper::generate_hyper_rectangle(PointWrapper &p1,
+ PointWrapper &p2,
+ const bool colorize)
+ {
+ // Create a new Triangulation
+ reset_triangulation();
+
+ if (dim == 2)
+ {
+ // Cast the PointWrapper object to Point<dim>
+ dealii::Point<2> point_1 = *(static_cast<dealii::Point<2>*>(p1.get_point()));
+ dealii::Point<2> point_2 = *(static_cast<dealii::Point<2>*>(p2.get_point()));
+
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+ dealii::GridGenerator::hyper_rectangle(*tria, point_1, point_2, colorize);
+ }
+ else
+ {
+ dealii::Point<3> point_1 = *(static_cast<dealii::Point<3>*>(p1.get_point()));
+ dealii::Point<3> point_2 = *(static_cast<dealii::Point<3>*>(p2.get_point()));
+
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(triangulation);
+ dealii::GridGenerator::hyper_rectangle(*tria, point_1, point_2, colorize);
+ }
+
+ }
+
+
+
+ void TriangulationWrapper::generate_subdivided_hyper_rectangle(boost::python::list &repetition_list,
+ PointWrapper &p1,
+ PointWrapper &p2,
+ const bool colorize)
+ {
+ // Create a new Triangulation
+ reset_triangulation();
+
+ // Extract the repetitions from the python list
+ std::vector<unsigned int> repetitions(dim);
+ for (int i=0; i<dim; ++i)
+ repetitions[i] = boost::python::extract<unsigned int>(repetition_list[i]);
+
+ if (dim == 2)
+ {
+ // Cast the PointWrapper object to Point<dim>
+ dealii::Point<2> point_1 = *(static_cast<dealii::Point<2>*>(p1.get_point()));
+ dealii::Point<2> point_2 = *(static_cast<dealii::Point<2>*>(p2.get_point()));
+
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+ dealii::GridGenerator::subdivided_hyper_rectangle(*tria, repetitions, point_1,
+ point_2, colorize);
+ }
+ else
+ {
+ dealii::Point<3> point_1 = *(static_cast<dealii::Point<3>*>(p1.get_point()));
+ dealii::Point<3> point_2 = *(static_cast<dealii::Point<3>*>(p2.get_point()));
+
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(triangulation);
+ dealii::GridGenerator::subdivided_hyper_rectangle(*tria, repetitions, point_1,
+ point_2, colorize);
+ }
+ }
+
+
+
+ void TriangulationWrapper::generate_hyper_ball(PointWrapper ¢er,
+ const double radius)
+ {
+ // Create a new Triangulation
+ reset_triangulation();
+
+ if (dim == 2)
+ {
+ // Cast the PointWrapper object to Point<dim>
+ dealii::Point<2> center_point = *(static_cast<dealii::Point<2>*>(
+ center.get_point()));
+
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+ dealii::GridGenerator::hyper_ball(*tria, center_point, radius);
+ }
+ else
+ {
+ dealii::Point<3> center_point = *(static_cast<dealii::Point<3>*>(
+ center.get_point()));
+
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(triangulation);
+ dealii::GridGenerator::hyper_ball(*tria, center_point, radius);
+ }
+ }
+
+
+
+ void TriangulationWrapper::shift(boost::python::list &shift_list)
+ {
+ if (dim == 2)
+ {
+ // Extract the shift vector from the python list
+ dealii::Tensor<1,2> shift_vector;
+ for (int i=0; i<dim; ++i)
+ shift_vector[i] = boost::python::extract<double>(shift_list[i]);
+
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+ dealii::GridTools::shift(shift_vector, *tria);
+ }
+ else
+ {
+ dealii::Tensor<1,3> shift_vector;
+ for (int i=0; i<dim; ++i)
+ shift_vector[i] = boost::python::extract<double>(shift_list[i]);
+
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(triangulation);
+ dealii::GridTools::shift(shift_vector, *tria);
+ }
+ }
+
+
+
+ void TriangulationWrapper::merge_triangulations(TriangulationWrapper &triangulation_1,
+ TriangulationWrapper &triangulation_2)
+ {
+ AssertThrow(triangulation_1.get_dim() == triangulation_2.get_dim(),
+ dealii::ExcMessage("Triangulation_1 and Triangulation_2 should have the same dimension."));
+
+ // Create a new Triangulation
+ reset_triangulation();
+
+ if (triangulation_1.get_dim() == 2)
+ {
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+ tria = new dealii::Triangulation<2>();
+ dealii::Triangulation<2> *tria_1 =
+ static_cast<dealii::Triangulation<2>*>(triangulation_1.get_triangulation());
+ dealii::Triangulation<2> *tria_2 =
+ static_cast<dealii::Triangulation<2>*>(triangulation_2.get_triangulation());
+ dealii::GridGenerator::merge_triangulations(*tria_1, *tria_2, *tria);
+ // We need to reassign tria to triangulation because tria was cleared
+ // inside merge_triangulations.
+ triangulation = tria;
+ }
+ else
+ {
+ TriangulationWrapper merged_triangulation("3D");
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(merged_triangulation.get_triangulation());
+ tria = new dealii::Triangulation<3>();
+ dealii::Triangulation<3> *tria_1 =
+ static_cast<dealii::Triangulation<3>*>(triangulation_1.get_triangulation());
+ dealii::Triangulation<3> *tria_2 =
+ static_cast<dealii::Triangulation<3>*>(triangulation_2.get_triangulation());
+ dealii::GridGenerator::merge_triangulations(*tria_1, *tria_2, *tria);
+ triangulation = tria;
+ }
+ }
+
+
+
+ void TriangulationWrapper::save(const std::string &filename)
+ {
+ std::ofstream ofs(filename);
+ boost::archive::text_oarchive oa(ofs);
+
+ if (dim == 2)
+ {
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+
+ oa << *tria;
+ }
+ else
+ {
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(triangulation);
+
+ oa << *tria;
+ }
+ }
+
+
+
+ void TriangulationWrapper::load(const std::string &filename)
+ {
+ std::ifstream ifs(filename);
+ boost::archive::text_iarchive ia(ifs);
+
+ if (dim == 2)
+ {
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+
+ ia >> *tria;
+ }
+ else
+ {
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(triangulation);
+
+ ia >> *tria;
+ }
+ }
+
+
+
+ void TriangulationWrapper::refine_global(const unsigned int times)
+ {
+ if (dim == 2)
+ {
+ dealii::Triangulation<2> *tria =
+ static_cast<dealii::Triangulation<2>*>(triangulation);
+ tria->refine_global(times);
+ }
+ else
+ {
+ dealii::Triangulation<3> *tria =
+ static_cast<dealii::Triangulation<3>*>(triangulation);
+ tria->refine_global(times);
+ }
+ }
+
+
+
+ void TriangulationWrapper::reset_triangulation()
+ {
+ if (dim == 2)
+ {
+ if (triangulation != nullptr)
+ {
+ // We cannot call delete on a void pointer so cast the void pointer back
+ // first.
+ dealii::Triangulation<2> *tmp = static_cast<dealii::Triangulation<2>*>(triangulation);
+ delete tmp;
+ }
+ triangulation = new dealii::Triangulation<2>();
+ }
+ else
+ {
+ if (triangulation != nullptr)
+ {
+ dealii::Triangulation<3> *tmp = static_cast<dealii::Triangulation<3>*>(triangulation);
+ delete tmp;
+ }
+ triangulation = new dealii::Triangulation<3> ();
+ }
+ }
+}