]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add more GridGenerator.
authorBruno Turcksin <bruno.turcksin@gmail.com>
Wed, 1 Jun 2016 13:40:50 +0000 (09:40 -0400)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Wed, 3 Aug 2016 20:34:02 +0000 (16:34 -0400)
source/python/CMakeLists.txt
source/python/export_point.cc [new file with mode: 0644]
source/python/export_triangulation.cc [new file with mode: 0644]
source/python/wrappers.cc

index baf781342ee180a9602834842596d9e2e9903831..886dce057e7a3d77df08c6fbf4bfa938dda27409 100644 (file)
@@ -21,6 +21,8 @@ INCLUDE_DIRECTORIES(SYSTEM ${PYTHON_INCLUDE_DIRS})
 
 SET(_src
   wrappers.cc
+  export_point.cc
+  export_triangulation.cc
   )
 
 FILE(GLOB _header
diff --git a/source/python/export_point.cc b/source/python/export_point.cc
new file mode 100644 (file)
index 0000000..1ddfbdb
--- /dev/null
@@ -0,0 +1,40 @@
+// ---------------------------------------------------------------------
+//
+// 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/base/point.h>
+
+namespace PyDealII
+{
+
+  double get_x(const dealii::Point<2> &point)
+  {
+    return point(0);
+  }
+
+  double get_y(const dealii::Point<2> &point)
+  {
+    return point(1);
+  }
+
+  void export_point()
+  {
+    boost::python::class_<dealii::Point<2>> ("Point", "Constructor. Requires x and y.",
+                                             boost::python::init<double, double>())
+                                         .add_property("x", get_x, "Return the x component of the point.")
+                                         .add_property("y", get_y, "Return the y component of the point.");
+  }
+
+}
diff --git a/source/python/export_triangulation.cc b/source/python/export_triangulation.cc
new file mode 100644 (file)
index 0000000..c626cc6
--- /dev/null
@@ -0,0 +1,177 @@
+// ---------------------------------------------------------------------
+//
+// 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/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>
+
+namespace PyDealII
+{
+
+  unsigned int n_active_cells(const dealii::Triangulation<2> &triangulation)
+  {
+    return triangulation.n_active_cells();
+  }
+
+
+
+  void generate_hyper_cube(dealii::Triangulation<2> &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);
+  }
+
+
+
+  void generate_subdivided_hyper_cube(dealii::Triangulation<2> &triangulation,
+                                      const unsigned int        repetitions,
+                                      const double              left = 0.,
+                                      const double              right = 1.)
+  {
+    dealii::GridGenerator::subdivided_hyper_cube(triangulation, repetitions,
+                                                 left, right);
+  }
+
+
+
+  void generate_hyper_rectangle(dealii::Triangulation<2> &triangulation,
+                                const dealii::Point<2>   &p1,
+                                const dealii::Point<2>   &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);
+  }
+
+
+
+  void generate_hyper_ball(dealii::Triangulation<2> &triangulation,
+                           const dealii::Point<2>   &center = dealii::Point<2>(),
+                           const double              radius = 1.)
+  {
+    dealii::GridGenerator::hyper_ball(triangulation, center, radius);
+  }
+
+
+
+  void save(const dealii::Triangulation<2> &triangulation,
+            const std::string filename)
+  {
+    std::ofstream ofs(filename);
+    boost::archive::text_oarchive oa(ofs);
+    oa << triangulation;
+  }
+
+
+
+  void load(dealii::Triangulation<2> &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)
+
+
+  void export_triangulation()
+  {
+    boost::python::class_<dealii::Triangulation<2>> ("Triangulation")
+                                                 .def("n_active_cells", n_active_cells,
+                                                      "Return the number of active cells",
+                                                      boost::python::args("self"))
+                                                 .def("generate_hyper_cube", generate_hyper_cube,
+                                                      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,
+                                                      generate_subdivided_hyper_cube_overloads(
+                                                        boost::python::args("self", "repetitions", "left", "right"),
+                                                        "Generate a subdivided hyper_cube."))
+                                                 .def("generate_hyper_rectangle", generate_hyper_rectangle,
+                                                      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,
+                                                      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, "Serialize and save the triangulation.",
+                                                      boost::python::args("self", "filename"))
+                                                 .def("load", load, "Load and deserialize a triangulation.",
+                                                      boost::python::args("self", "filename"));
+  }
+
+}
index 6a2e1b75630b1a13902df0869f52dff27ae669cb..3c0b6c3d39c3c9c557fa9a3b8592912068f29611 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>
+
+namespace PyDealII
+{
+  void export_point();
+  void export_triangulation();
+}
 
 char const *pydealii_docstring =
   "                                                             \n"
@@ -27,36 +28,6 @@ char const *pydealii_docstring =
   "Some interesting doc.                                        \n"
   ;
 
-unsigned int n_active_cells(const dealii::Triangulation<2> &triangulation)
-{
-  return triangulation.n_active_cells();
-}
-
-void generate_hyper_cube(dealii::Triangulation<2> &triangulation, const double left=0.,
-                         const double right=0., const bool colorize=false)
-{
-  dealii::GridGenerator::hyper_cube(triangulation);
-}
-
-void save(const dealii::Triangulation<2> &triangulation,
-          const std::string filename)
-{
-  std::ofstream ofs(filename);
-  boost::archive::text_oarchive oa(ofs);
-  oa << triangulation;
-}
-
-void load(dealii::Triangulation<2> &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_MODULE(PyDealII)
 {
   boost::python::scope().attr("__doc__") = pydealii_docstring;
@@ -66,20 +37,6 @@ BOOST_PYTHON_MODULE(PyDealII)
   doc_options.enable_py_signatures();
   doc_options.disable_cpp_signatures();
 
-
-  boost::python::class_<dealii::Triangulation<2>> ("Triangulation")
-                                               .def("n_active_cells", n_active_cells,
-                                                    "Return the number of active cells",
-                                                    boost::python::args("self"))
-                                               .def("generate_hyper_cube", generate_hyper_cube,
-                                                    generate_hyper_cube_overloads(
-                                                      boost::python::args("self", "left", "right", "colorize"),
-                                                      "Generate a hyper_cube."))
-                                               .def("refine_global", &dealii::Triangulation<2>::refine_global,
-                                                    "Refine the mesh uniformly",
-                                                    boost::python::args("self", "times"))
-                                               .def("save", save, "Serialize and save the triangulation",
-                                                    boost::python::args("self", "filename"))
-                                               .def("load", load, "Load and deserialize a triangulation",
-                                                    boost::python::args("self", "filename"));
+  PyDealII::export_point();
+  PyDealII::export_triangulation();
 }

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.