]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Extend cell python wrappers and add face wrappers 9023/head
authorAlexander Grayver <agrayver@erdw.ethz.ch>
Fri, 8 Nov 2019 15:35:21 +0000 (16:35 +0100)
committerAlexander Grayver <agrayver@erdw.ethz.ch>
Thu, 14 Nov 2019 14:12:21 +0000 (15:12 +0100)
contrib/python-bindings/include/cell_accessor_wrapper.h
contrib/python-bindings/include/tria_accessor_wrapper.h [new file with mode: 0644]
contrib/python-bindings/source/CMakeLists.txt
contrib/python-bindings/source/cell_accessor_wrapper.cc
contrib/python-bindings/source/export_cell_accessor.cc
contrib/python-bindings/source/export_tria_accessor.cc [new file with mode: 0644]
contrib/python-bindings/source/tria_accessor_wrapper.cc [new file with mode: 0644]
contrib/python-bindings/source/wrappers.cc
contrib/python-bindings/tests/cell_accessor_wrapper.py
doc/news/changes/minor/20191114Grayver [new file with mode: 0644]

index 9ea48bce3e49c55a73871bbc4eb3694bd05c21d5..4fb18eb2a66661e0d93ad7805f647d198d3aafb4 100644 (file)
@@ -16,6 +16,8 @@
 #ifndef dealii_cell_accessor_wrapper_h
 #define dealii_cell_accessor_wrapper_h
 
+#include <boost/python.hpp>
+
 #include <deal.II/base/config.h>
 
 #include <deal.II/grid/tria_accessor.h>
@@ -43,6 +45,11 @@ namespace python
                         const int             level,
                         const int             index);
 
+    /**
+     * Constructor for an empty object.
+     */
+    CellAccessorWrapper();
+
     /**
      * Destructor.
      */
@@ -107,6 +114,31 @@ namespace python
      */
     int get_manifold_id() const;
 
+    /**
+     * Return the ith neighbor of a cell. If the neighbor does not exist, 
+     * i.e., if the ith face of the current object is at the boundary, 
+     * then an exception is thrown.
+     */
+    CellAccessorWrapper neighbor(const int i) const;
+
+    /**
+     * Return whether the cell is at the boundary.
+     */
+    bool at_boundary() const;
+
+    /**
+     * This is a slight variation to the at_boundary function: 
+     * for 2 dimensions it is equivalent, for three 
+     * dimensions it returns whether at least one of the 12 
+     * lines of the hexahedron is at a boundary.
+     */
+    bool has_boundary_lines() const;
+
+    /**
+     * Get faces of the underlying cell.
+     */ 
+    boost::python::list faces() const;
+
     /**
      * Exception.
      */
@@ -116,7 +148,15 @@ namespace python
                    << " does not exist. The largest vertex number "
                    << "acceptable is "<< arg2-1);
 
+    DeclException2(ExcNeighborDoesNotExist,
+                   int, int,
+                   << "Requested neighbor number " << arg1
+                   << " does not exist. The largest neighbor number "
+                   << "acceptable is "<< arg2-1);
   private:
+    template<int dim, int spacedim>
+    const CellAccessorWrapper construct_neighbor_wrapper(const int i) const;
+
     /**
      * Dimension of the underlying CellAccessor object.
      */
@@ -133,6 +173,22 @@ namespace python
     void *cell_accessor;
 
   };
+
+
+  template<int dim, int spacedim>
+  const CellAccessorWrapper CellAccessorWrapper::construct_neighbor_wrapper(const int i) const
+  {
+       CellAccessor<dim,spacedim> *cell =
+               static_cast<CellAccessor<dim,spacedim>*>(cell_accessor);
+
+       auto neighbor = cell->neighbor(i);
+
+       CellAccessorWrapper neighbor_wrapper;
+        neighbor_wrapper.dim = dim;
+        neighbor_wrapper.spacedim = spacedim;
+       neighbor_wrapper.cell_accessor = new CellAccessor<dim,spacedim>(&neighbor->get_triangulation(), neighbor->level(), neighbor->index());
+       return neighbor_wrapper;
+  }
 }
 
 DEAL_II_NAMESPACE_CLOSE
diff --git a/contrib/python-bindings/include/tria_accessor_wrapper.h b/contrib/python-bindings/include/tria_accessor_wrapper.h
new file mode 100644 (file)
index 0000000..76e6216
--- /dev/null
@@ -0,0 +1,128 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_tria_accessor_wrapper_h
+#define dealii_tria_accessor_wrapper_h
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/grid/tria_accessor.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace python
+{
+  class PointWrapper;
+  class TriangulationWrapper;
+
+  class TriaAccessorWrapper
+  {
+  public:
+    /**
+     * Copy constructor.
+     */
+    TriaAccessorWrapper(const TriaAccessorWrapper &other);
+
+    /**
+     * Constructor.
+     */
+    TriaAccessorWrapper(void* tria_accessor, const int structdim, const int dim, const int spacedim);
+
+    /**
+     * Destructor.
+     */
+    ~TriaAccessorWrapper();
+
+    /**
+     * Get the barycenter of the cell.
+     */
+    PointWrapper get_barycenter() const;
+
+    /**
+     * Set the ith vertex of the cell to @p point_wrapper.
+     */
+    void set_vertex(const int i, PointWrapper &point_wrapper);
+
+    /**
+     * Return the ith vertex of the cell.
+     */
+    PointWrapper get_vertex(const int i) const;
+
+    /**
+     * Set the manifold id.
+     */
+    void set_manifold_id(const int manifold_id);
+
+    /**
+     * Get the manifold id.
+     */
+    int get_manifold_id() const;
+
+    /**
+     * Set the boundary id.
+     */
+    void set_boundary_id(const int boundary_id);
+
+    /**
+     * Get the boundary id.
+     */
+    int get_boundary_id() const;
+
+    /**
+     * Set the boundary id for all objects.
+     */
+    void set_all_boundary_ids(const int boundary_id);
+
+    /**
+     * Return whether the cell is at the boundary.
+     */
+    bool at_boundary() const;
+
+    /**
+     * Exception.
+     */
+    DeclException2(ExcVertexDoesNotExist,
+                   int, int,
+                   << "Requested vertex number " << arg1
+                   << " does not exist. The largest vertex number "
+                   << "acceptable is "<< arg2-1);
+
+  private:
+    /**
+     * Dimension of the underlying TriaAccessor object
+     */
+    int structdim;
+
+    /**
+     * Dimension of the triangulation which TriaAccessor object
+     * is linked to.
+     */
+    int dim;
+
+    /**
+     * Space dimension of the underlying TriaAccessor object.
+     */
+    int spacedim;
+
+    /**
+     * Pointer that can be casted to the underlying TriaAccessor object.
+     */
+    void *tria_accessor;
+  };
+}
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index 508754b1fb39e724a7596132a0c1a3541de885c5..9a6ac24d88c8faba370e0511f872fc612b2582e3 100644 (file)
@@ -90,10 +90,12 @@ INCLUDE_DIRECTORIES(SYSTEM ${PYTHON_INCLUDE_DIRS})
 
 SET(_src
   wrappers.cc
+  export_tria_accessor.cc
   export_cell_accessor.cc
   export_point.cc
   export_triangulation.cc
   cell_accessor_wrapper.cc
+  tria_accessor_wrapper.cc
   point_wrapper.cc
   triangulation_wrapper.cc
   )
index 5151410a0f450310fb12bedf5b764f5721e28626..2eae0865385c7c447f8b5c7855566c1c2cd2061c 100644 (file)
@@ -18,6 +18,7 @@
 #include <boost/python.hpp>
 #include <point_wrapper.h>
 #include <triangulation_wrapper.h>
+#include <tria_accessor_wrapper.h>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -240,6 +241,56 @@ namespace python
         static_cast<const CellAccessor<dim,spacedim>*>(cell_accessor);
       return cell->manifold_id();
     }
+
+
+
+    template <int dim, int spacedim>
+    bool at_boundary(const void *cell_accessor)
+    {
+      const CellAccessor<dim,spacedim> *cell =
+        static_cast<const CellAccessor<dim,spacedim>*>(cell_accessor);
+      return cell->at_boundary();
+    }
+
+
+
+    template <int dim, int spacedim>
+    bool has_boundary_lines(const void *cell_accessor)
+    {
+      const CellAccessor<dim,spacedim> *cell =
+        static_cast<const CellAccessor<dim,spacedim>*>(cell_accessor);
+      return cell->has_boundary_lines();
+    }
+
+
+
+    template <int dim, int spacedim>
+    const CellAccessor<dim,spacedim>* neighbor(const int i, const void *cell_accessor)
+    {
+      const CellAccessor<dim,spacedim> *cell =
+        static_cast<const CellAccessor<dim,spacedim>*>(cell_accessor);
+      return cell->neighbor(i);
+    }
+
+
+
+    template <int dim, int spacedim>
+    boost::python::list faces(const void *cell_accessor)
+    {
+      const CellAccessor<dim,spacedim> *cell =
+        static_cast<const CellAccessor<dim,spacedim>*>(cell_accessor);
+
+      boost::python::list faces_list;
+
+      auto face_iterators = cell->face_iterators();
+      for(auto &it : face_iterators)
+      {
+        TriaAccessor<dim-1,dim,spacedim> *face_accessor = new TriaAccessor<dim-1,dim,spacedim>(*it);
+        faces_list.append(TriaAccessorWrapper(face_accessor, dim - 1, dim, spacedim));
+      }
+
+      return faces_list;
+    }
   }
 
 
@@ -273,6 +324,15 @@ namespace python
 
 
 
+  CellAccessorWrapper::CellAccessorWrapper()
+    :
+    dim(0),
+    spacedim(0),
+    cell_accessor(nullptr)
+  {
+  }
+
+
 
   CellAccessorWrapper::CellAccessorWrapper(TriangulationWrapper &triangulation_wrapper,
                                            const int             level,
@@ -318,8 +378,8 @@ namespace python
           }
         else if ((dim == 2) && (spacedim == 3))
           {
-            CellAccessor<3,3> *tmp =
-              static_cast<CellAccessor<3,3>*>(cell_accessor);
+            CellAccessor<2,3> *tmp =
+              static_cast<CellAccessor<2,3>*>(cell_accessor);
             delete tmp;
           }
         else
@@ -384,6 +444,7 @@ namespace python
   }
 
 
+
   PointWrapper CellAccessorWrapper::get_barycenter() const
   {
     if ((dim == 2) && (spacedim == 2))
@@ -472,6 +533,58 @@ namespace python
     else
       return internal::get_manifold_id<3,3>(cell_accessor);
   }
+
+
+
+  bool CellAccessorWrapper::at_boundary() const
+  {
+    if ((dim == 2) && (spacedim == 2))
+      return internal::at_boundary<2,2>(cell_accessor);
+    else if ((dim== 2) && (spacedim == 3))
+      return internal::at_boundary<2,3>(cell_accessor);
+    else
+      return internal::at_boundary<3,3>(cell_accessor);
+  }
+
+
+
+  bool CellAccessorWrapper::has_boundary_lines() const
+  {
+    if ((dim == 2) && (spacedim == 2))
+      return internal::has_boundary_lines<2,2>(cell_accessor);
+    else if ((dim== 2) && (spacedim == 3))
+      return internal::has_boundary_lines<2,3>(cell_accessor);
+    else
+      return internal::has_boundary_lines<3,3>(cell_accessor);
+  }
+
+
+
+  CellAccessorWrapper CellAccessorWrapper::neighbor(const int i) const
+  {
+    AssertThrow(i<(2*dim),
+                ExcNeighborDoesNotExist(i, 2*dim));
+
+    if ((dim == 2) && (spacedim == 2))
+      return construct_neighbor_wrapper<2, 2>(i);
+    else if ((dim== 2) && (spacedim == 3))
+      return construct_neighbor_wrapper<2, 3>(i);
+    else
+      return construct_neighbor_wrapper<3, 3>(i);
+  }
+
+
+
+  boost::python::list CellAccessorWrapper::faces() const
+  {
+    if ((dim == 2) && (spacedim == 2))
+      return internal::faces<2, 2>(cell_accessor);
+    else if ((dim== 2) && (spacedim == 3))
+      return internal::faces<2, 3>(cell_accessor);
+    else
+      return internal::faces<3, 3>(cell_accessor);
+  }
+
 }
 
 DEAL_II_NAMESPACE_CLOSE
index 17084ee42a8ffca46d3d48fa23523beb179ab22d..809ea64c5e7127e2c2afe09ee777a43ff9415acd 100644 (file)
@@ -78,6 +78,35 @@ namespace python
 
 
 
+   const char at_boundary_docstring [] =
+    " Return whether the cell is at the boundary                        \n"
+    ; 
+
+
+
+   const char has_boundary_lines_docstring [] =
+    " This is a slight variation to the at_boundary function:           \n"
+    " for two dimensions, it is equivalent, for three                   \n"
+    " dimensions it returns whether at least one of the 12              \n" 
+    " lines of the hexahedron is at a boundary.                         \n"
+    ;
+
+
+
+   const char neighbor_docstring [] =
+    " Return the ith neighbor of a cell. If the neighbor does not exist,\n" 
+    " i.e., if the ith face of the current object is at the boundary,   \n" 
+    " then an exception is thrown.                                      \n"
+    ; 
+
+
+
+   const char faces_docstring [] =
+    " Return list of cell's faces.                                      \n"
+    ; 
+
+
+
   void export_cell_accessor()
   {
     boost::python::class_<CellAccessorWrapper>("CellAccessor",
@@ -102,7 +131,19 @@ namespace python
          boost::python::args("self", "i", "point_wrapper"))
     .def("get_vertex", &CellAccessorWrapper::get_vertex,
          get_vertex_docstring,
-         boost::python::args("self", "i"));
+         boost::python::args("self", "i"))
+    .def("at_boundary", &CellAccessorWrapper::at_boundary,
+         at_boundary_docstring,
+         boost::python::args("self"))
+    .def("has_boundary_lines", &CellAccessorWrapper::has_boundary_lines,
+         has_boundary_lines_docstring,
+         boost::python::args("self"))
+    .def("neighbor", &CellAccessorWrapper::neighbor,
+         neighbor_docstring,
+         boost::python::args("self", "i"))
+    .def("faces", &CellAccessorWrapper::faces,
+         faces_docstring,
+         boost::python::args("self"));
   }
 }
 
diff --git a/contrib/python-bindings/source/export_tria_accessor.cc b/contrib/python-bindings/source/export_tria_accessor.cc
new file mode 100644 (file)
index 0000000..fcb7c68
--- /dev/null
@@ -0,0 +1,81 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 - 2017 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#include <tria_accessor_wrapper.h>
+#include <triangulation_wrapper.h>
+#include <boost/python.hpp>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace python
+{
+  const char manifold_id_docstring [] =
+    "Get/Set the manifold_id of the face                                \n"
+    ;
+
+  const char boundary_id_docstring [] =
+    "Get/Set the boundary_id of the face                                \n"
+    ;
+
+  const char set_all_boundary_ids_docstring [] =
+    "Do as set_boundary_id() but also set the boundary indicators       \n"
+    "of the objects that bound the current object.                      \n"  
+    ;
+
+  const char barycenter_docstring [] =
+    "Return the barycenter of the current face                          \n"
+    ;
+
+  const char set_vertex_docstring [] =
+    " Set the ith vertex of the face to point_wrapper                   \n"
+    ;
+
+  const char get_vertex_docstring [] =
+    " Get the ith vertex of the face                                    \n"
+    ;
+
+   const char at_boundary_docstring [] =
+    " Return whether the face is at the boundary                        \n"
+    ; 
+
+  void export_tria_accessor()
+  {
+    boost::python::class_<TriaAccessorWrapper>("TriaAccessor",
+                                               boost::python::init<void*, const int, const int, const int> ())
+    .add_property("boundary_id", &TriaAccessorWrapper::get_boundary_id,
+                  &TriaAccessorWrapper::set_boundary_id,
+                  boundary_id_docstring)
+    .add_property("manifold_id", &TriaAccessorWrapper::get_manifold_id,
+                  &TriaAccessorWrapper::set_manifold_id,
+                  manifold_id_docstring)
+    .def("barycenter", &TriaAccessorWrapper::get_barycenter,
+         barycenter_docstring,
+         boost::python::args("self"))
+    .def("set_vertex", &TriaAccessorWrapper::set_vertex,
+         set_vertex_docstring,
+         boost::python::args("self", "i", "point_wrapper"))
+    .def("get_vertex", &TriaAccessorWrapper::get_vertex,
+         get_vertex_docstring,
+         boost::python::args("self", "i"))
+    .def("at_boundary", &TriaAccessorWrapper::at_boundary,
+         at_boundary_docstring,
+         boost::python::args("self"))
+    .def("set_all_boundary_ids", &TriaAccessorWrapper::set_all_boundary_ids,
+         set_all_boundary_ids_docstring,
+         boost::python::args("self", "boundary_id"));
+  }
+}
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/contrib/python-bindings/source/tria_accessor_wrapper.cc b/contrib/python-bindings/source/tria_accessor_wrapper.cc
new file mode 100644 (file)
index 0000000..6104853
--- /dev/null
@@ -0,0 +1,321 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 - 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#include <tria_accessor_wrapper.h>
+
+#include <boost/python.hpp>
+#include <point_wrapper.h>
+#include <triangulation_wrapper.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace python
+{
+  namespace internal
+  {
+    template <int structdim, int dim, int spacedim>
+    PointWrapper get_barycenter(const void *tria_accessor)
+    {
+      const TriaAccessor<structdim,dim,spacedim> *accessor =
+        static_cast<const TriaAccessor<structdim,dim,spacedim>*>(tria_accessor);
+      Point<spacedim> barycenter = accessor->barycenter();
+      boost::python::list barycenter_list;
+      for (int i=0; i<dim; ++i)
+        barycenter_list.append(barycenter[i]);
+
+      return PointWrapper(barycenter_list);
+    }
+
+
+
+    template <int structdim, int dim, int spacedim>
+    void set_boundary_id(const int boundary_id, void *tria_accessor)
+    {
+      TriaAccessor<structdim,dim,spacedim> *accessor =
+        static_cast<TriaAccessor<structdim, dim,spacedim>*>(tria_accessor);
+      accessor->set_boundary_id(boundary_id);
+    }
+
+
+
+    template <int structdim, int dim, int spacedim>
+    int get_boundary_id(const void *tria_accessor)
+    {
+      const TriaAccessor<structdim,dim,spacedim> *accessor =
+        static_cast<const TriaAccessor<structdim,dim,spacedim>*>(tria_accessor);
+
+      return accessor->boundary_id();
+    }
+
+
+
+    template <int structdim, int dim, int spacedim>
+    void set_all_boundary_ids(const int boundary_id, void *tria_accessor)
+    {
+      TriaAccessor<structdim,dim,spacedim> *accessor =
+        static_cast<TriaAccessor<structdim,dim,spacedim>*>(tria_accessor);
+      accessor->set_all_boundary_ids(boundary_id);
+    }
+
+
+    template <int structdim, int dim, int spacedim>
+    void set_vertex(const int     i,
+                    PointWrapper &point_wrapper,
+                    void         *tria_accessor)
+    {
+      TriaAccessor<structdim,dim,spacedim> *accessor =
+        static_cast<TriaAccessor<structdim,dim,spacedim>*>(tria_accessor);
+      Point<spacedim> *point =
+        static_cast<Point<spacedim>*>(point_wrapper.get_point());
+
+      accessor->vertex(i) = *point;
+    }
+
+
+
+    template <int structdim, int dim, int spacedim>
+    PointWrapper get_vertex(const int i, const void *tria_accessor)
+    {
+      const TriaAccessor<structdim,dim,spacedim> *accessor =
+        static_cast<const TriaAccessor<structdim,dim,spacedim>*>(tria_accessor);
+      Point<spacedim> vertex = accessor->vertex(i);
+
+      boost::python::list coordinates;
+      for (int i=0; i<spacedim; ++i)
+        coordinates.append(vertex[i]);
+
+      return PointWrapper(coordinates);
+    }
+
+
+
+    template <int structdim, int dim, int spacedim>
+    void set_manifold_id(const int manifold_id, void *tria_accessor)
+    {
+      TriaAccessor<structdim,dim,spacedim> *accessor =
+        static_cast<TriaAccessor<structdim,dim,spacedim>*>(tria_accessor);
+      accessor->set_manifold_id(manifold_id);
+    }
+
+
+
+    template <int structdim, int dim, int spacedim>
+    int get_manifold_id(const void *tria_accessor)
+    {
+      const TriaAccessor<structdim,dim,spacedim> *accessor =
+        static_cast<const TriaAccessor<structdim,dim,spacedim>*>(tria_accessor);
+      return accessor->manifold_id();
+    }
+
+
+
+    template <int structdim, int dim, int spacedim>
+    bool at_boundary(const void *tria_accessor)
+    {
+      const TriaAccessor<structdim,dim,spacedim> *accessor =
+        static_cast<const TriaAccessor<structdim,dim,spacedim>*>(tria_accessor);
+      return accessor->at_boundary();
+    }
+  }
+
+
+
+  TriaAccessorWrapper::TriaAccessorWrapper(const TriaAccessorWrapper &other)
+    :
+    structdim(other.structdim),
+    dim(other.dim),
+    spacedim(other.spacedim)
+  {
+    if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+      {
+        TriaAccessor<1,2,2> *other_accessor =
+          static_cast<TriaAccessor<1,2,2>*>(other.tria_accessor);
+        tria_accessor = new TriaAccessor<1,2,2>(*other_accessor);
+      }
+    else if ((dim == 2) && (spacedim == 3) && (structdim == 1))
+      {
+        TriaAccessor<1,2,3> *other_accessor =
+          static_cast<TriaAccessor<1,2,3>*>(other.tria_accessor);
+        tria_accessor = new TriaAccessor<1,2,3>(*other_accessor);
+      }
+    else if ((dim == 3) && (spacedim == 3) && (structdim == 2))
+      {
+        TriaAccessor<2,3,3> *other_accessor =
+          static_cast<TriaAccessor<2,3,3>*>(other.tria_accessor);
+        tria_accessor = new TriaAccessor<2,3,3>(*other_accessor);
+      }
+    else
+      AssertThrow(false, ExcMessage("Wrong structdim-dim-spacedim combination."));
+  }
+
+  TriaAccessorWrapper::TriaAccessorWrapper(void* tria_accessor, 
+                                           const int structdim, 
+                                           const int dim, 
+                                           const int spacedim):
+    structdim(structdim),
+    dim(dim),
+    spacedim(spacedim),
+    tria_accessor(tria_accessor)
+  {}
+
+
+  TriaAccessorWrapper::~TriaAccessorWrapper()
+  {
+    if (dim != -1)
+      {
+        if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+          {
+            // We cannot call delete on a void pointer so cast the void pointer back
+            // first.
+            TriaAccessor<1,2,2> *tmp =
+              static_cast<TriaAccessor<1,2,2>*>(tria_accessor);
+            delete tmp;
+          }
+        else if ((dim == 2) && (spacedim == 3) && (structdim == 1))
+          {
+            TriaAccessor<1,2,3> *tmp =
+              static_cast<TriaAccessor<1,2,3>*>(tria_accessor);
+            delete tmp;
+          }
+        else
+          {
+            TriaAccessor<2,3,3> *tmp =
+              static_cast<TriaAccessor<2,3,3>*>(tria_accessor);
+            delete tmp;
+          }
+
+        dim = -1;
+        spacedim = -1;
+        structdim = -1;
+        tria_accessor = nullptr;
+      }
+  }
+
+
+
+  PointWrapper TriaAccessorWrapper::get_barycenter() const
+  {
+    if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+      return internal::get_barycenter<1,2,2>(tria_accessor);
+    else if ((dim == 2) && (spacedim == 3) && (structdim == 1))
+      return internal::get_barycenter<1,2,3>(tria_accessor);
+    else
+      return internal::get_barycenter<2,3,3>(tria_accessor);
+  }
+
+
+
+  void TriaAccessorWrapper::set_boundary_id(const int boundary_id)
+  {
+    if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+      internal::set_boundary_id<1,2,2>(boundary_id, tria_accessor);
+    else if ((dim == 2) && (spacedim == 3) && (structdim == 1))
+      internal::set_boundary_id<1,2,3>(boundary_id, tria_accessor);
+    else
+      internal::set_boundary_id<2,3,3>(boundary_id, tria_accessor);
+  }
+
+
+
+  void TriaAccessorWrapper::set_all_boundary_ids(const int boundary_id)
+  {
+    if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+      internal::set_all_boundary_ids<1,2,2>(boundary_id, tria_accessor);
+    else if ((dim == 2) && (spacedim == 3) && (structdim == 1))
+      internal::set_all_boundary_ids<1,2,3>(boundary_id, tria_accessor);
+    else
+      internal::set_all_boundary_ids<2,3,3>(boundary_id, tria_accessor);
+  }
+
+
+
+  int TriaAccessorWrapper::get_boundary_id() const
+  {
+    if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+      return internal::get_boundary_id<1,2,2>(tria_accessor);
+    else if ((dim == 2) && (spacedim == 3) && (structdim == 1))
+      return internal::get_boundary_id<1,2,3>(tria_accessor);
+    else
+      return internal::get_boundary_id<2,3,3>(tria_accessor);
+  }
+
+
+
+  void TriaAccessorWrapper::set_vertex(const int     i,
+                                       PointWrapper &point_wrapper)
+  {
+    AssertThrow(i<std::pow(2,dim),
+                ExcVertexDoesNotExist(i, std::pow(2,dim)));
+    if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+      internal::set_vertex<1,2,2>(i, point_wrapper, tria_accessor);
+    else if ((dim == 2) && (spacedim == 3) && (structdim == 1))
+      internal::set_vertex<1,2,3>(i, point_wrapper, tria_accessor);
+    else
+      internal::set_vertex<2,3,3>(i, point_wrapper, tria_accessor);
+  }
+
+
+
+  PointWrapper TriaAccessorWrapper::get_vertex(const int i) const
+  {
+    AssertThrow(i<std::pow(2,dim),
+                ExcVertexDoesNotExist(i, std::pow(2,dim)));
+    if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+      return internal::get_vertex<1,2,2>(i, tria_accessor);
+    else if ((dim == 2) && (spacedim == 3) && (structdim == 1))
+      return internal::get_vertex<1,2,3>(i, tria_accessor);
+    else
+      return internal::get_vertex<2,3,3>(i, tria_accessor);
+  }
+
+
+
+  void TriaAccessorWrapper::set_manifold_id(const int manifold_id)
+  {
+    if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+      internal::set_manifold_id<1,2,2>(manifold_id, tria_accessor);
+    else if ((dim == 2) && (spacedim == 3) && (structdim == 1))
+      internal::set_manifold_id<1,2,3>(manifold_id, tria_accessor);
+    else
+      internal::set_manifold_id<2,3,3>(manifold_id, tria_accessor);
+  }
+
+
+
+  int TriaAccessorWrapper::get_manifold_id() const
+  {
+    if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+      return internal::get_manifold_id<1,2,2>(tria_accessor);
+    else if ((dim== 2) && (spacedim == 3) && (structdim == 1))
+      return internal::get_manifold_id<1,2,3>(tria_accessor);
+    else
+      return internal::get_manifold_id<2,3,3>(tria_accessor);
+  }
+
+
+
+  bool TriaAccessorWrapper::at_boundary() const
+  {
+    if ((dim == 2) && (spacedim == 2) && (structdim == 1))
+      return internal::at_boundary<1,2,2>(tria_accessor);
+    else if ((dim== 2) && (spacedim == 3) && (structdim == 1))
+      return internal::at_boundary<1,2,3>(tria_accessor);
+    else
+      return internal::at_boundary<2,3,3>(tria_accessor);
+  }
+
+}
+
+DEAL_II_NAMESPACE_CLOSE
index 889cbea6133fa48de522ff6cb9aa7f849fb892a4..52bba15d0494d1ae4839e3bb4407f83a1f3b1226 100644 (file)
@@ -22,6 +22,7 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace python
 {
+  void export_tria_accessor();
   void export_cell_accessor();
   void export_point();
   void export_triangulation();
@@ -54,6 +55,7 @@ BOOST_PYTHON_MODULE(Debug)
   // message is printed.
   dealii::deal_II_exceptions::disable_abort_on_exception();
 
+  dealii::python::export_tria_accessor();
   dealii::python::export_cell_accessor();
   dealii::python::export_point();
   dealii::python::export_triangulation();
@@ -70,6 +72,7 @@ BOOST_PYTHON_MODULE(Release)
   doc_options.enable_py_signatures();
   doc_options.disable_cpp_signatures();
 
+  dealii::python::export_tria_accessor();
   dealii::python::export_cell_accessor();
   dealii::python::export_point();
   dealii::python::export_triangulation();
index fd28b63d93ce9d1a9e0ea8db49a6fe6fa1d7d604..662208ccc4353a7e82d729beee2d97ffe547aa91 100644 (file)
@@ -16,7 +16,6 @@
 import unittest
 from PyDealII.Debug import *
 
-
 class TestCellAccessorWrapper(unittest.TestCase):
 
     def setUp(self):
@@ -24,6 +23,25 @@ class TestCellAccessorWrapper(unittest.TestCase):
         self.triangulation.generate_hyper_cube()
         self.triangulation.refine_global(1)
 
+    def test_faces(self):
+        for cell in self.triangulation.active_cells():
+            faces = cell.faces()
+            self.assertEqual(len(faces), 4)
+
+    def test_at_boundary(self):
+        for cell in self.triangulation.active_cells():
+            self.assertEqual(cell.at_boundary(), cell.has_boundary_lines()) 
+
+    def test_faces(self):
+        n_neighbors = 0
+        for cell in self.triangulation.active_cells():
+            faces = cell.faces()
+            for i in range(len(faces)):
+                if not faces[i].at_boundary():
+                    neighbor = cell.neighbor(i)
+                    n_neighbors += 1
+        self.assertEqual(n_neighbors, 8)
+
     def test_material_id(self):
         material_id = 0
         for cell in self.triangulation.active_cells():
diff --git a/doc/news/changes/minor/20191114Grayver b/doc/news/changes/minor/20191114Grayver
new file mode 100644 (file)
index 0000000..0593e08
--- /dev/null
@@ -0,0 +1,4 @@
+Added: extend python bindings to allow for working with cell's
+neighbors, iterate over faces and modify boundary ids. 
+<br>
+(Alexander Grayver, 2019/11/14)

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.