PointWrapper(boost::python::list list);
/**
- * Destructpr.
+ * Copy constructor.
+ */
+ PointWrapper(const PointWrapper &other);
+
+ /**
+ * Destructor.
*/
~PointWrapper();
+ /**
+ * Assignment operator. The dimension of the point is changed if it is
+ * different that the one of @p other.
+ */
+ PointWrapper &operator= (const PointWrapper &other);
+
/**
* Return the first component of the Point.
*/
*/
void set_z(double z);
+ /**
+ * Return the dimension of the underlying Point.
+ */
+ int get_dim() const;
+
/**
* Return a pointer that can be casted to the underlying Point.
*/
void *get_point();
private:
+ /**
+ * Delete the underlying Point and free the memory.
+ */
+ void clear();
+
+
+ /**
+ * Copy @p other PointWrapper.
+ */
+ void copy(const PointWrapper &other);
+
/**
* Dimension of the Point.
*/
+ inline
+ int PointWrapper::get_dim() const
+ {
+ return dim;
+ }
+
+
+
inline
void *PointWrapper::get_point()
{
{
// Macro to enable default arguments
- BOOST_PYTHON_MEMBER_FUNCTION_OVERLOADS(generate_hyper_cube_overloads, generate_hyper_cube, 0, 3)
+ BOOST_PYTHON_MEMBER_FUNCTION_OVERLOADS(generate_hyper_cube_overloads, generate_hyper_cube, 1, 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,
+ PointWrapper::PointWrapper(const PointWrapper &other)
+ {
+ copy(other);
+ }
+
+
+
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;
- }
+ clear();
dim = -1;
}
+ PointWrapper &PointWrapper::operator= (const PointWrapper &other)
+ {
+ clear();
+ copy(other);
+
+ return *this;
+ }
+
+
+
double PointWrapper::get_x() const
{
if (dim == 2)
AssertThrow(false,
dealii::ExcMessage("The z coordinate is only available for three-dimensional points"));
}
+
+
+
+
+ void PointWrapper::clear()
+ {
+ 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;
+ }
+ }
+
+
+
+ void PointWrapper::copy(const PointWrapper &other)
+ {
+ dim = other.dim;
+
+ if (dim == 2)
+ {
+ dealii::Point<2> *other_point = static_cast<dealii::Point<2>*>(other.point);
+ point = new dealii::Point<2> ((*other_point)[0], (*other_point)[1]);
+ }
+ else if (dim == 3)
+ {
+ dealii::Point<3> *other_point = static_cast<dealii::Point<3>*>(other.point);
+ point = new dealii::Point<3> ((*other_point)[0], (*other_point)[1], (*other_point)[2]);
+ }
+ else
+ point = nullptr;
+ }
}
void *triangulation)
{
// Cast the PointWrapper objects to Point<dim>
- std::vector<dealii::Point<dim>> points(dim);
- for (int i=0; i<dim; ++i)
+ std::vector<dealii::Point<dim>> points(dim+1);
+ for (int i=0; i<dim+1; ++i)
points[i] = *(static_cast<dealii::Point<dim>*>((wrapped_points[i]).get_point()));
dealii::Triangulation<dim> *tria =
const bool colorize,
void *triangulation)
{
+ AssertThrow(p1.get_dim() == dim,
+ dealii::ExcMessage("Dimension of p1 is not the same as the dimension of the Triangulation."));
+ AssertThrow(p2.get_dim() == dim,
+ dealii::ExcMessage("Dimension of p2 is not the same as the dimension of the Triangulation."));
// Cast the PointWrapper object to Point<dim>
dealii::Point<dim> point_1 = *(static_cast<dealii::Point<dim>*>(p1.get_point()));
dealii::Point<dim> point_2 = *(static_cast<dealii::Point<dim>*>(p2.get_point()));
const bool colorize,
void *triangulation)
{
+ AssertThrow(p1.get_dim() == dim,
+ dealii::ExcMessage("Dimension of p1 is not the same as the dimension of the Triangulation."));
+ AssertThrow(p2.get_dim() == dim,
+ dealii::ExcMessage("Dimension of p2 is not the same as the dimension of the Triangulation."));
// Cast the PointWrapper object to Point<dim>
dealii::Point<dim> point_1 = *(static_cast<dealii::Point<dim>*>(p1.get_point()));
dealii::Point<dim> point_2 = *(static_cast<dealii::Point<dim>*>(p2.get_point()));
void TriangulationWrapper::generate_simplex(boost::python::list &vertices)
{
+ AssertThrow(boost::python::len(vertices) == dim+1,
+ dealii::ExcMessage("The number of vertices should be equal to dim+1."));
// 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]);
+ std::vector<PointWrapper> wrapped_points(dim+1);
+ for (int i=0; i<dim+1; ++i)
+ {
+ wrapped_points[i] = boost::python::extract<PointWrapper>(vertices[i]);
+ AssertThrow(wrapped_points[i].get_dim() == dim,
+ dealii::ExcMessage("Point of wrong dimension."));
+ }
if (dim == 2)
internal::generate_simplex<2>(wrapped_points, triangulation);
PointWrapper &p2,
const bool colorize)
{
+ AssertThrow(boost::python::len(repetition_list) == dim,
+ dealii::ExcMessage("The list of repetitions must have the same length as the number of dimension."));
+
// Extract the repetitions from the python list
std::vector<unsigned int> repetitions(dim);
for (int i=0; i<dim; ++i)
--- /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 TestPointWrapper(unittest.TestCase):
+
+ def test_2d_point(self):
+ point = Point([0., 1.])
+ self.assertEqual(point.x, 0.)
+ self.assertEqual(point.y, 1.)
+ point.x = 1.
+ point.y = 2.
+ self.assertEqual(point.x, 1.)
+ self.assertEqual(point.y, 2.)
+
+ def test_3d_point(self):
+ point = Point([0., 1., 2.])
+ self.assertEqual(point.x, 0.)
+ self.assertEqual(point.y, 1.)
+ self.assertEqual(point.z, 2.)
+ point.x = 1.
+ point.y = 2.
+ point.z = 3.
+ self.assertEqual(point.x, 1.)
+ self.assertEqual(point.y, 2.)
+ self.assertEqual(point.z, 3.)
+
+
+if __name__ == '__main__':
+ unittest.main()
--- /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 TestTriangulationWrapper(unittest.TestCase):
+
+ def setUp(self):
+ self.dim = ['2D', '3D']
+
+ def build_hyper_cube_triangulation(self, dim):
+ triangulation_1 = Triangulation(dim)
+ triangulation_1.generate_hyper_cube()
+ return triangulation_1
+
+ def build_hyper_rectangle_triangulation(self, dim):
+ triangulation_2 = Triangulation(dim)
+ if (dim == '2D'):
+ point_1 = Point([0., 0.])
+ point_2 = Point([1., 1.])
+ else:
+ point_1 = Point([0., 0., 0.])
+ point_2 = Point([1., 1., 1.])
+ triangulation_2.generate_hyper_rectangle(point_1, point_2)
+ return triangulation_2
+
+ def test_hyper_cube(self):
+ for dim in self.dim:
+ triangulation_1 = self.build_hyper_cube_triangulation(dim)
+ n_cells = triangulation_1.n_active_cells()
+ self.assertEqual(n_cells, 1)
+
+ def test_simplex(self):
+ for dim in self.dim:
+ triangulation = Triangulation(dim)
+ if (dim == '2D'):
+ point_1 = Point([0., 0.])
+ point_2 = Point([1., 0.])
+ point_3 = Point([1., 1.])
+ vertices = [point_1, point_2, point_3]
+ else:
+ point_1 = Point([0., 0., 0.])
+ point_2 = Point([1., 0., 0.])
+ point_3 = Point([1., 1., 0.])
+ point_4 = Point([1., 1., 1.])
+ vertices = [point_1, point_2, point_3, point_4]
+ triangulation.generate_simplex(vertices)
+ n_cells = triangulation.n_active_cells()
+ self.assertEqual(n_cells, len(vertices))
+
+ def test_subdivided_hyper_cube(self):
+ for dim in self.dim:
+ triangulation = Triangulation(dim)
+ repetitions = 2
+ triangulation.generate_subdivided_hyper_cube(repetitions)
+ n_cells = triangulation.n_active_cells()
+ if (dim == '2D'):
+ self.assertEqual(n_cells, 4)
+ else:
+ self.assertEqual(n_cells, 8)
+
+ def test_hyper_rectangle(self):
+ for dim in self.dim:
+ triangulation_2 = self.build_hyper_rectangle_triangulation(dim)
+ n_cells = triangulation_2.n_active_cells()
+ self.assertEqual(n_cells, 1)
+
+ def test_subdivided_hyper_rectangle(self):
+ for dim in self.dim:
+ triangulation = Triangulation(dim)
+ if (dim == '2D'):
+ repetitions = [6, 4]
+ point_1 = Point([0., 0.])
+ point_2 = Point([1., 1.])
+ else:
+ repetitions = [2, 3, 4]
+ point_1 = Point([0., 0., 0.])
+ point_2 = Point([1., 1., 1.])
+ triangulation.generate_subdivided_hyper_rectangle(repetitions,
+ point_1,
+ point_2)
+ n_cells = triangulation.n_active_cells()
+ self.assertEqual(n_cells, 24)
+
+ def test_hyper_ball(self):
+ for dim in self.dim:
+ triangulation = Triangulation(dim)
+ if (dim == '2D'):
+ center = Point([0., 0.])
+ n_cells_ref = 5
+ else:
+ center = Point([0., 0., 0.])
+ n_cells_ref = 7
+ triangulation.generate_hyper_ball(center)
+ n_cells = triangulation.n_active_cells()
+ self.assertEqual(n_cells, n_cells_ref)
+
+ def test_shift_and_merge(self):
+ for dim in self.dim:
+ triangulation_1 = self.build_hyper_cube_triangulation(dim)
+ triangulation_2 = self.build_hyper_rectangle_triangulation(dim)
+ if (dim == '2D'):
+ triangulation_2.shift([1., 0.])
+ else:
+ triangulation_2.shift([1., 0., 0.])
+ triangulation = Triangulation(dim)
+ triangulation.merge_triangulations(triangulation_1,
+ triangulation_2)
+ n_cells = triangulation.n_active_cells()
+ self.assertEqual(n_cells, 2)
+
+ def test_refine_global(self):
+ for dim in self.dim:
+ triangulation = self.build_hyper_cube_triangulation(dim)
+ triangulation.refine_global(2)
+ n_cells = triangulation.n_active_cells()
+ if (dim == '2D'):
+ self.assertEqual(n_cells, 16)
+ else:
+ self.assertEqual(n_cells, 64)
+
+ def test_save_load(self):
+ for dim in self.dim:
+ triangulation_1 = self.build_hyper_cube_triangulation(dim)
+ triangulation_1.refine_global(1)
+ triangulation_1.save('mesh.output')
+ triangulation_2 = Triangulation(dim)
+ triangulation_2.load('mesh.output')
+ n_cells = triangulation_2.n_active_cells()
+ if (dim == '2D'):
+ self.assertEqual(n_cells, 4)
+ else:
+ self.assertEqual(n_cells, 8)
+
+
+if __name__ == '__main__':
+ unittest.main()