]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add 3D for Triangulation and Point interface.
authorBruno Turcksin <bruno.turcksin@gmail.com>
Fri, 10 Jun 2016 12:02:10 +0000 (08:02 -0400)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Wed, 3 Aug 2016 20:34:02 +0000 (16:34 -0400)
include/deal.II/python/point_wrapper.h
include/deal.II/python/triangulation_wrapper.h [new file with mode: 0644]
source/python/CMakeLists.txt
source/python/export_point.cc
source/python/export_triangulation.cc
source/python/point_wrapper.cc
source/python/triangulation_wrapper.cc [new file with mode: 0644]

index b3f081a786765baf4e1214e7b51f1117783dcde5..73c40ac938d619e7409a3fa1bdcadcac05162e14 100644 (file)
@@ -22,18 +22,79 @@ namespace PyDealII
   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
diff --git a/include/deal.II/python/triangulation_wrapper.h b/include/deal.II/python/triangulation_wrapper.h
new file mode 100644 (file)
index 0000000..cfcd16d
--- /dev/null
@@ -0,0 +1,165 @@
+// ---------------------------------------------------------------------
+//
+// 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 &center,
+                             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
index 7607ad8beca1ab854d773cbe7c80f7ac604d71ac..6ee5c7261b1564709393772757ff9822409623e2 100644 (file)
@@ -24,6 +24,7 @@ SET(_src
   export_point.cc
   export_triangulation.cc
   point_wrapper.cc
+  triangulation_wrapper.cc
   )
 
 FILE(GLOB _header
index b2f9d42cde25222a08b1bf7a020584105850307b..f908f1673fb6b5b61b819962879a850993eca29e 100644 (file)
@@ -22,7 +22,7 @@ namespace PyDealII
 {
   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.");
index 9fe9c4accccec4a234b1a136ea4b9ffac48a9cf0..9ee2fd177c77cb90362507944f81c113479ca25d 100644 (file)
 //
 // ---------------------------------------------------------------------
 
-#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>   &center = 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"));
   }
-
 }
index 07529ac7263db9e6a0fccddb51e0b689a8a516d9..8459bd5ddf5fe10ddb073eaf069284d610e5425b 100644 (file)
 
 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."));
@@ -35,12 +42,35 @@ namespace PyDealII
 
 
 
+  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);
   }
 
 
@@ -48,9 +78,9 @@ namespace PyDealII
   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;
   }
 
 
@@ -58,9 +88,9 @@ namespace PyDealII
   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);
   }
 
 
@@ -68,9 +98,9 @@ namespace PyDealII
   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;
   }
 
 
@@ -78,7 +108,7 @@ namespace PyDealII
   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"));
@@ -89,7 +119,7 @@ namespace PyDealII
   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"));
diff --git a/source/python/triangulation_wrapper.cc b/source/python/triangulation_wrapper.cc
new file mode 100644 (file)
index 0000000..9871759
--- /dev/null
@@ -0,0 +1,412 @@
+// ---------------------------------------------------------------------
+//
+// 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 &center,
+                                                 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> ();
+      }
+  }
+}

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.