diff --git a/meshpy/triangle.py b/meshpy/triangle.py new file mode 100644 index 0000000000000000000000000000000000000000..32d5b3131560302b5c913314da8e86c782d4800e --- /dev/null +++ b/meshpy/triangle.py @@ -0,0 +1,134 @@ +from meshpy.common import MeshInfoBase +import meshpy._triangle as internals + + + + +class MeshInfo(internals.MeshInfo, MeshInfoBase): + _constituents = [ + "points", "point_attributes", "point_markers", + "elements", "element_attributes", "element_volumes", + "neighbors", + "segments", "segment_markers", + "holes", + "regions", + "edges", "edge_markers", + "normals", + ] + + def __getstate__(self): + def dump_array(array): + try: + return [ + [array[i,j] for j in range(array.unit)] + for i in range(len(array))] + except RuntimeError: + # not allocated + return None + + return self.number_of_point_attributes, \ + self.number_of_element_attributes, \ + [(name, dump_array(getattr(self, name))) for name in self._constituents] + + def __setstate__(self, (p_attr_count, e_attr_count, state)): + self.number_of_point_attributes = p_attr_count + self.number_of_element_attributes = e_attr_count + for name, array in state: + if name not in self._constituents: + raise RuntimeError, "Unknown constituent during unpickling" + + dest_array = getattr(self, name) + + if array is None: + dest_array.deallocate() + else: + if len(dest_array) != len(array): + dest_array.resize(len(array)) + + for i,tup in enumerate(array): + for j,v in enumerate(tup): + dest_array[i, j] = v + + def set_segments(self, segments, segment_markers=None): + self.segments.resize(len(segments)) + + for i, seg in enumerate(segments): + self.segments[i] = seg + + if segment_markers is not None: + assert len(segment_markers) == len(segments) + for i, mark in enumerate(segment_markers): + self.SegmentMarkers[i] = mark + + + + + + + + +def build(mesh_info, verbose=False, refinement_func=None): + opts = "pzqj" + if verbose: + opts += "VV" + else: + opts += "Q" + if refinement_func is not None: + opts += "u" + + mesh = MeshInfo() + internals.triangulate(opts, mesh_info, mesh, MeshInfo(), refinement_func) + return mesh + + + + +def refine(input_p, verbose=False, refinement_func=None): + opts = "razj" + if input_p.Segments.size() != 0: + opts += "p" + if verbose: + opts += "VV" + else: + opts += "Q" + if refinement_func is not None: + opts += "u" + output_p = MeshInfo() + internals.triangulate(opts, input_p, output_p, MeshInfo(), refinement_func) + return output_p + + + + +def dump_parameters(par): + def dump_array(name, array): + subs = array.unit() + print "array %s: %d elements, %d values per element" % (name, array.size(), subs) + if array.size() == 0 or array.unit() == 0: + return + try: + array.get(0) + except RuntimeError: + print " not allocated" + return + + for i in range(array.size()): + print " %d:" % i, + for j in range(subs): + print array.get_sub(i,j), + print + + for name, arr in par._constituents: + dump_array(name, getattr(par, name)) + + + + +def write_gnuplot_mesh(filename, out_p): + gp_file = file(filename, "w") + + for points in out_p.elements: + for pt in points: + gp_file.write("%f %f 0\n" % tuple(out_p.points[pt])) + gp_file.write("%f %f 0\n\n" % tuple(out_p.points[points[0]])) + diff --git a/src/foreign_array.hpp b/src/foreign_array.hpp index 0f72be5de90e71c4fbb3f3af620a163e5f2acd0d..7e8820c8a402c02a2a278489dcf093068d86801c 100644 --- a/src/foreign_array.hpp +++ b/src/foreign_array.hpp @@ -6,6 +6,7 @@ #include #include +#include @@ -68,19 +69,25 @@ class tSizeChangeNotifier template -class tForeignArray : public tSizeChangeNotifier, public tSizeChangeNotificationReceiver, +class tReadOnlyForeignArray : public tSizeChangeNotifier, public tSizeChangeNotificationReceiver, public boost::noncopyable { + protected: ElementT *&Contents; int &NumberOf; unsigned Unit; tSizeChangeNotifier *SlaveTo; std::string Name; + bool DeleteOnDestruction; public: - tForeignArray(const std::string &name, - ElementT *&cts, int &number_of, unsigned unit = 1, tSizeChangeNotifier *slave_to = NULL) - : Contents(cts), NumberOf(number_of), Unit(unit), SlaveTo(slave_to), Name(name) + typedef ElementT value_type; + + tReadOnlyForeignArray(const std::string &name, + ElementT *&cts, int &number_of, unsigned unit=1, tSizeChangeNotifier *slave_to=NULL, + bool delete_on_destruction=false) + : Contents(cts), NumberOf(number_of), Unit(unit), SlaveTo(slave_to), Name(name), + DeleteOnDestruction(delete_on_destruction) { Contents = NULL; if (SlaveTo) @@ -92,11 +99,14 @@ class tForeignArray : public tSizeChangeNotifier, public tSizeChangeNotification setSize(0); } - ~tForeignArray() + ~tReadOnlyForeignArray() { if (SlaveTo) SlaveTo->unregisterForNotification(this); + if (DeleteOnDestruction) + deallocate(); + if (!SlaveTo) NumberOf = 0; } @@ -114,7 +124,7 @@ class tForeignArray : public tSizeChangeNotifier, public tSizeChangeNotification void deallocate() { if (Contents != NULL) - free(Contents); + delete[] Contents; Contents = NULL; } @@ -146,7 +156,7 @@ class tForeignArray : public tSizeChangeNotifier, public tSizeChangeNotification NumberOf = size; if (Contents != NULL) - free(Contents); + deallocate(); if (size == 0 || Unit == 0) Contents = NULL; @@ -169,50 +179,66 @@ class tForeignArray : public tSizeChangeNotifier, public tSizeChangeNotification } } - void set(unsigned index, ElementT value) + ElementT &get(unsigned index) { if (index >= NumberOf * Unit) throw std::runtime_error("index out of bounds"); if (Contents == NULL) throw std::runtime_error("Array unallocated"); - Contents[ index ] = value; + return Contents[ index ]; } - void setSub(unsigned index, unsigned sub_index, ElementT value) + ElementT &getSub(unsigned index, unsigned sub_index) { - set(index * Unit + sub_index, value); + return get(index * Unit + sub_index); } +}; + + + + + +template +class tForeignArray : public tReadOnlyForeignArray +{ + typedef tReadOnlyForeignArray super; - ElementT get(unsigned index) + public: + tForeignArray(const std::string &name, + ElementT *&cts, int &number_of, unsigned unit=1, + tSizeChangeNotifier *slave_to=NULL, + bool delete_on_destruction=false) + : super(name, cts, number_of, unit, slave_to, delete_on_destruction) { - if (index >= NumberOf * Unit) + } + + void set(unsigned index, ElementT value) + { + if (index >= this->NumberOf * this->Unit) throw std::runtime_error("index out of bounds"); - if (Contents == NULL) + if (this->Contents == NULL) throw std::runtime_error("Array unallocated"); - return Contents[ index ]; + this->Contents[ index ] = value; } - ElementT getSub(unsigned index, unsigned sub_index) + void setSub(unsigned index, unsigned sub_index, ElementT value) { - return get(index * Unit + sub_index); + set(index * this->Unit + sub_index, value); } tForeignArray &operator=(tForeignArray const &src) { - if (SlaveTo) - assert(src.size() == SlaveTo->size()); + if (this->SlaveTo) + assert(src.size() == this->SlaveTo->size()); else setSize(src.size()); setUnit(src.Unit); if (src.Contents) - memcpy(Contents, src.Contents, sizeof(ElementT) * Unit * src.size()); + memcpy(this->Contents, src.Contents, sizeof(ElementT) * this->Unit * src.size()); else - { - free(Contents); - Contents = NULL; - } + this->deallocate(); return *this; } diff --git a/src/foreign_array_wrap.hpp b/src/foreign_array_wrap.hpp new file mode 100644 index 0000000000000000000000000000000000000000..4d5860667e900050f4b34af0959a37789dc4dcd4 --- /dev/null +++ b/src/foreign_array_wrap.hpp @@ -0,0 +1,175 @@ +#ifndef _HEADER_SEEN_FOREIGN_ARRAY_WRAP +#define _HEADER_SEEN_FOREIGN_ARRAY_WRAP + + + + +#include "foreign_array.hpp" +#include + + + + +using namespace boost::python; + + + + +#define PYTHON_ERROR(TYPE, REASON) \ +{ \ + PyErr_SetString(PyExc_##TYPE, REASON); \ + throw error_already_set(); \ +} + + + + +namespace { + /* This wrap helper works as long as the value_type is a plain old data (POD) + * type. + * + * In exchange for this, it nicely wraps the "unit" abstraction provided by + * foreign arrays. + */ + template + struct tPODForeignArrayWrapHelper + { + typedef typename FA::value_type value_type; + + static object getitem(FA &self, long idx) + { + if (idx < 0) idx += self.size(); + if (idx >= (long) self.size()) PYTHON_ERROR(IndexError, "index out of bounds"); + + if (self.unit() > 1) + { + list l; + for (unsigned i = 0; i(idx[0]); + long i_sub = extract(idx[1]); + + if (i_main < 0) idx += self.size(); + if (i_main >= (long) self.size()) PYTHON_ERROR(IndexError, "index out of bounds"); + if (i_sub < 0) idx += self.unit(); + if (i_sub >= (long) self.unit()) PYTHON_ERROR(IndexError, "subindex out of bounds"); + + return object(self.getSub(i_main, i_sub)); + } + + static void setitem(FA &self, long idx, object value) + { + if (idx < 0) idx += self.size(); + if (idx >= (long) self.size()) PYTHON_ERROR(IndexError, "index out of bounds"); + + if (self.unit() > 1) + { + if ((long) self.unit() != len(value)) + PYTHON_ERROR(ValueError, "value must be a sequence of length self.unit"); + + for (long i = 0; i(value[i])); + } + else + self.set(idx, extract(value)); + } + + static void setitem(FA &self, tuple idx, const value_type &v) + { + if (len(idx) != 2) + PYTHON_ERROR(IndexError, "expected index tuple of length 2"); + long i_main = extract(idx[0]); + long i_sub = extract(idx[1]); + + if (i_main < 0) idx += self.size(); + if (i_main >= (long) self.size()) PYTHON_ERROR(IndexError, "index out of bounds"); + if (i_sub < 0) idx += self.unit(); + if (i_sub >= (long) self.unit()) PYTHON_ERROR(IndexError, "subindex out of bounds"); + + self.setSub(i_main, i_sub, v); + } + }; + + + + + /* This wrap helper works for more complicated data structures, for which we + * just ship out internal references--boost::python takes care of life support + * for us. + * + * In exchange for this, it does not allow setting entries or support the unit API. + */ + template + struct tStructureForeignArrayWrapHelper + { + typedef typename FA::value_type value_type; + + static value_type &getitem(FA &self, long idx) + { + if (idx < 0) idx += self.size(); + if (idx >= (long) self.size()) PYTHON_ERROR(IndexError, "index out of bounds"); + + return self.get(idx); + } + }; +} + + + + +template +void exposePODForeignArray(const std::string &name) +{ + typedef tForeignArray cl; + typedef tPODForeignArrayWrapHelper w_cl; + typedef typename cl::value_type value_type; + + boost::python::class_(name.c_str(), no_init) + .def("__len__", &cl::size) + .def("resize", &cl::setSize) + .def("setup", &cl::setup) + .add_property("unit", &cl::unit) + .def("__getitem__", (object (*)(cl &, long)) &w_cl::getitem) + .def("__getitem__", (object (*)(cl &, tuple)) &w_cl::getitem) + .def("__setitem__", (void (*)(cl &, long, object)) &w_cl::setitem) + .def("__setitem__", (void (*)(cl &, tuple, const value_type &)) &w_cl::setitem) + .def("deallocate", &cl::deallocate) + ; +} + + + + + +template +void exposeStructureForeignArray(const std::string &name) +{ + typedef tForeignArray cl; + typedef tStructureForeignArrayWrapHelper w_cl; + typedef typename cl::value_type value_type; + + boost::python::class_(name.c_str(), no_init) + .def("__len__", &cl::size) + .def("resize", &cl::setSize) + .def("setup", &cl::setup) + .add_property("unit", &cl::unit) + .def("__getitem__", &w_cl::getitem, return_internal_reference<>()) + .def("deallocate", &cl::deallocate) + ; +} + + + + + +#endif diff --git a/src/tetgen.cpp b/src/tetgen.cpp index 348a256569c9319e6c970e3fff9bbd50fff39f6b..1b0c805dc6147dff4173e1b3c4a9e7445b9112d0 100644 --- a/src/tetgen.cpp +++ b/src/tetgen.cpp @@ -47,6 +47,34 @@ void terminatetetgen(int x) // Begin of class 'tetgenio' implementation // +tetgenio::polygon::polygon() +{ + vertexlist = (int *) NULL; + numberofvertices = 0; +} + +tetgenio::polygon::~polygon() +{ + if (vertexlist) + delete [] vertexlist; +} + +tetgenio::facet::facet() +{ + polygonlist = (polygon *) NULL; + numberofpolygons = 0; + holelist = (REAL *) NULL; + numberofholes = 0; +} + +tetgenio::facet::~facet() +{ + if (polygonlist) + delete[] polygonlist; + if (holelist) + delete[] holelist; +} + /////////////////////////////////////////////////////////////////////////////// // // // initialize() Initialize all variables of 'tetgenio'. // @@ -181,19 +209,9 @@ void tetgenio::deinitialize() } if (facetlist != (facet *) NULL) { - for (i = 0; i < numberoffacets; i++) { - f = &facetlist[i]; - for (j = 0; j < f->numberofpolygons; j++) { - p = &f->polygonlist[j]; - delete [] p->vertexlist; - } - delete [] f->polygonlist; - if (f->holelist != (REAL *) NULL) { - delete [] f->holelist; - } - } delete [] facetlist; } + if (facetmarkerlist != (int *) NULL) { delete [] facetmarkerlist; } diff --git a/src/tetgen.h b/src/tetgen.h index 55b8f3a6041aadef49e4737fca5fe9371d0bb85e..d0aa42e52219eebc1125dd0f3498ec5b6b25dc97 100644 --- a/src/tetgen.h +++ b/src/tetgen.h @@ -135,6 +135,7 @@ #include // String lib: strcpy(), strcat(), strcmp(), ... #include // Math lib: sin(), sqrt(), pow(), ... #include +#include /////////////////////////////////////////////////////////////////////////////// // // @@ -190,10 +191,13 @@ class tetgenio { // 'vertexlist' is a list of vertex indices (integers), its length is // indicated by 'numberofvertices'. The vertex indices are odered in // either counterclockwise or clockwise way. - typedef struct { + struct polygon : public boost::noncopyable { int *vertexlist; int numberofvertices; - } polygon; + + polygon(); + ~polygon(); + }; static void init(polygon* p) { p->vertexlist = (int *) NULL; @@ -204,12 +208,15 @@ class tetgenio { // to represent a planar straight line graph (PSLG) in two dimension. // A PSLG contains a list of polygons. It also may conatin holes in it, // indicated by a list of hole points (their coordinates). - typedef struct { + struct facet : public boost::noncopyable { polygon *polygonlist; int numberofpolygons; REAL *holelist; int numberofholes; - } facet; + + facet(); + ~facet(); + }; static void init(facet* f) { f->polygonlist = (polygon *) NULL; diff --git a/src/triangle.c b/src/triangle.c index f7a57004aa3dd897c1c9b9f4c24ca8955150b2e1..eba297377a5df2102232eff91b32f8ae0dd36e8b 100644 --- a/src/triangle.c +++ b/src/triangle.c @@ -532,13 +532,6 @@ struct osub { int ssorient; /* Ranges from 0 to 1. */ }; -/* The vertex data structure. Each vertex is actually an array of REALs. */ -/* The number of REALs is unknown until runtime. An integer boundary */ -/* marker, and sometimes a pointer to a triangle, is appended after the */ -/* REALs. */ - -typedef REAL *vertex; - /* A queue used to store encroached subsegments. Each subsegment's vertices */ /* are stored so that we can check whether a subsegment is still the same. */ diff --git a/src/triangle.h b/src/triangle.h index 9df1f39ea45df28eeb315f4a7c6a59af4eae3151..3f08ddec29cec9f45a545b278cdb99928cff3714 100644 --- a/src/triangle.h +++ b/src/triangle.h @@ -1,3 +1,13 @@ +#ifndef _HEADER_SEEN_TRIANGLE_H +#define _HEADER_SEEN_TRIANGLE_H + + + + +#ifdef __cplusplus +extern "C" +{ +#endif /*****************************************************************************/ /* */ /* (triangle.h) */ @@ -248,6 +258,23 @@ /* */ /*****************************************************************************/ +#ifdef SINGLE +#define REAL float +#else /* not SINGLE */ +#define REAL double +#endif /* not SINGLE */ + +#define VOID void + +/* The vertex data structure. Each vertex is actually an array of REALs. */ +/* The number of REALs is unknown until runtime. An integer boundary */ +/* marker, and sometimes a pointer to a triangle, is appended after the */ +/* REALs. */ + +typedef REAL *vertex; + +int triunsuitable(vertex triorg, vertex tridest, vertex triapex, REAL area); + struct triangulateio { REAL *pointlist; /* In / out */ REAL *pointattributelist; /* In / out */ @@ -287,3 +314,11 @@ void trifree(VOID *memptr); void triangulate(); void trifree(); #endif /* not ANSI_DECLARATORS */ +#ifdef __cplusplus +} +#endif + + + + +#endif diff --git a/src/wrap_tetgen.cpp b/src/wrap_tetgen.cpp index 83334442e729c430266744161139db9896382b2c..18b6c80c67d6ef4604581b1f7c70f6f0dd2a7b51 100644 --- a/src/wrap_tetgen.cpp +++ b/src/wrap_tetgen.cpp @@ -3,7 +3,7 @@ #include #include #include -#include "foreign_array.hpp" +#include "foreign_array_wrap.hpp" @@ -14,7 +14,7 @@ using namespace std; -struct tMeshDescriptor : public tetgenio, public boost::noncopyable +struct tMeshInfo : public tetgenio, public boost::noncopyable { public: tForeignArray Points; // in/out @@ -30,8 +30,11 @@ struct tMeshDescriptor : public tetgenio, public boost::noncopyable tForeignArray Facets; + tForeignArray Holes; + tForeignArray Regions; + public: - tMeshDescriptor() + tMeshInfo() : Points("points", pointlist, numberofpoints, 3), PointAttributes("point_attributes", pointattributelist, numberofpoints, 0, &Points), AdditionalPoints("additional_points", addpointlist, numberofaddpoints, 3), @@ -46,16 +49,18 @@ struct tMeshDescriptor : public tetgenio, public boost::noncopyable Neighbors("neighbors", neighborlist, numberoftetrahedra, 4, &Elements), - Facets("facets", facetlist, numberoffacets) + Facets("facets", facetlist, numberoffacets), /* Segments("Segments", segmentlist, numberofsegments, 2), SegmentMarkers("SegmentMarkers", segmentmarkerlist, numberofsegments, 1, &Segments), + */ Holes("Holes", holelist, numberofholes, 3), - Regions("Regions", regionlist, numberofregions, 4), + Regions("Regions", regionlist, numberofregions, 5) + /* Edges("Edges", edgelist, numberofedges, 2), EdgeMarkers("EdgeMarkers", edgemarkerlist, numberofedges, 1, &Edges), Normals("Normals", normlist, numberofedges, 2, &Edges) @@ -135,45 +140,59 @@ tTriangulationParameters *copyTriangulationParameters(const tTriangulationParame -void tetrahedralizeWrapper(tetgenbehavior &bhv, tMeshDescriptor &in, tMeshDescriptor &out) +void tetrahedralizeWrapper(tetgenbehavior &bhv, tMeshInfo &in, tMeshInfo &out) { tetrahedralize(&bhv, &in, &out); } -template -void exposeForeignArray(T, const string &name) + +template +struct manage_new_internal_reference + : with_custodian_and_ward_postcall<0, owner_arg, Base> +{ + typedef manage_new_object result_converter; +}; + + + + +tForeignArray *facet_get_polygons(tetgenio::facet &self) +{ + return new tForeignArray("polygons", + self.polygonlist, self.numberofpolygons); +} + + + + +tForeignArray *facet_get_holes(tetgenio::facet &self) +{ + return new tForeignArray("holes", self.holelist, self.numberofholes); +} + + + + + +tForeignArray *polygon_get_vertices(tetgenio::polygon &self) { - typedef tForeignArray cl; - - class_ - (name.c_str(), no_init) - .def("__len__", &cl::size) - .def("size", &cl::size) - .def("set_size", &cl::setSize) - .def("setup", &cl::setup) - .def("unit", &cl::unit) - .def("set", &cl::set) - .def("set_sub", &cl::setSub) - .def("get", &cl::get) - .def("get_sub", &cl::getSub) - .def("deallocate", &cl::deallocate) - ; + return new tForeignArray("vertices", self.vertexlist, self.numberofvertices); } -BOOST_PYTHON_MODULE(internals) +BOOST_PYTHON_MODULE(_tetgen) { def("tetrahedralize", tetrahedralizeWrapper); { - typedef tMeshDescriptor cl; + typedef tMeshInfo cl; class_ - ("MeshDescriptor", init<>()) + ("MeshInfo", init<>()) .def_readonly("points", &cl::Points) .def_readonly("point_attributes", &cl::PointAttributes) .def_readonly("point_markers", &cl::PointMarkers) @@ -188,11 +207,15 @@ BOOST_PYTHON_MODULE(internals) /* .def_readonly("Segments", &cl::Segments) .def_readonly("SegmentMarkers", &cl::SegmentMarkers) + */ + + .def_readonly("facets", &cl::Facets) - .def_readonly("Holes", &cl::Holes) + .def_readonly("holes", &cl::Holes) - .def_readonly("Regions", &cl::Regions) + .def_readonly("regions", &cl::Regions) + /* .def_readonly("Edges", &cl::Edges) .def_readonly("EdgeMarkers", &cl::EdgeMarkers) @@ -213,13 +236,83 @@ BOOST_PYTHON_MODULE(internals) //.enable_pickling() ; } + { typedef tetgenio::facet cl; - class_ - ("Facet", no_init) + class_("Facet", no_init) + .def("get_polygons", facet_get_polygons, manage_new_internal_reference<>()) + .def("get_holes", facet_get_holes, manage_new_internal_reference<>()) + ; + } + + { + typedef tetgenio::polygon cl; + class_("Polygon", no_init) + .def("get_vertices", polygon_get_vertices, manage_new_internal_reference<>()) + ; + } + + { +#define DEF_RW_MEMBER(NAME) \ + def_readwrite(#NAME, &cl::NAME) + + typedef tetgenbehavior cl; + class_("Options", init<>()) + .DEF_RW_MEMBER(plc) + .DEF_RW_MEMBER(refine) + .DEF_RW_MEMBER(quality) + .DEF_RW_MEMBER(smooth) + .DEF_RW_MEMBER(metric) + .DEF_RW_MEMBER(bgmesh) + .DEF_RW_MEMBER(varvolume) + .DEF_RW_MEMBER(fixedvolume) + .DEF_RW_MEMBER(insertaddpoints) + .DEF_RW_MEMBER(regionattrib) + .DEF_RW_MEMBER(offcenter) + .DEF_RW_MEMBER(conformdel) + .DEF_RW_MEMBER(diagnose) + .DEF_RW_MEMBER(zeroindex) + .DEF_RW_MEMBER(order) + .DEF_RW_MEMBER(facesout) + .DEF_RW_MEMBER(edgesout) + .DEF_RW_MEMBER(neighout) + .DEF_RW_MEMBER(meditview) + .DEF_RW_MEMBER(gidview) + .DEF_RW_MEMBER(geomview) + .DEF_RW_MEMBER(nobound) + .DEF_RW_MEMBER(nonodewritten) + .DEF_RW_MEMBER(noelewritten) + .DEF_RW_MEMBER(nofacewritten) + .DEF_RW_MEMBER(noiterationnum) + .DEF_RW_MEMBER(nomerge) + .DEF_RW_MEMBER(nobisect) + .DEF_RW_MEMBER(noflip) + .DEF_RW_MEMBER(nojettison) + .DEF_RW_MEMBER(steiner) + .DEF_RW_MEMBER(fliprepair) + .DEF_RW_MEMBER(docheck) + .DEF_RW_MEMBER(quiet) + .DEF_RW_MEMBER(verbose) + .DEF_RW_MEMBER(tol) + .DEF_RW_MEMBER(useshelles) + .DEF_RW_MEMBER(minratio) + .DEF_RW_MEMBER(goodratio) + .DEF_RW_MEMBER(minangle) + .DEF_RW_MEMBER(goodangle) + .DEF_RW_MEMBER(maxvolume) + .DEF_RW_MEMBER(maxdihedral) + .DEF_RW_MEMBER(alpha1) + .DEF_RW_MEMBER(alpha2) + .DEF_RW_MEMBER(alpha3) + .DEF_RW_MEMBER(epsilon) + .DEF_RW_MEMBER(epsilon2) + + .def("parse_switches", (bool (tetgenbehavior::*)(char *)) &cl::parse_commandline) ; } - exposeForeignArray(REAL(), "RealArray"); - exposeForeignArray(int(), "IntArray"); + exposePODForeignArray("RealArray"); + exposePODForeignArray("IntArray"); + exposeStructureForeignArray("FacetArray"); + exposeStructureForeignArray("PolygonArray"); } diff --git a/src/wrap_triangle.cpp b/src/wrap_triangle.cpp index 29b1fa78a2af80493d4938a5ba4a15da7fea6838..977d529d672c6818f7f13957a5d52c1b43e12ede 100644 --- a/src/wrap_triangle.cpp +++ b/src/wrap_triangle.cpp @@ -2,7 +2,7 @@ #include #include #include -#include "foreign_array.hpp" +#include "foreign_array_wrap.hpp" @@ -13,7 +13,7 @@ using namespace std; -struct tMeshDescriptor : public triangulateio, public boost::noncopyable +struct tMeshInfo : public triangulateio, public boost::noncopyable { public: tForeignArray Points; // in/out @@ -37,29 +37,29 @@ struct tMeshDescriptor : public triangulateio, public boost::noncopyable tForeignArray Normals; // out only public: - tMeshDescriptor() - : Points("points", pointlist, numberofpoints, 2), - PointAttributes("point_attributes", pointattributelist, numberofpoints, 0, &Points), - PointMarkers("point_markers", pointmarkerlist, numberofpoints, 1, &Points), + tMeshInfo() + : Points("points", pointlist, numberofpoints, 2, NULL, true), + PointAttributes("point_attributes", pointattributelist, numberofpoints, 0, &Points, true), + PointMarkers("point_markers", pointmarkerlist, numberofpoints, 1, &Points, true), - Elements("elements", trianglelist, numberoftriangles, 3), + Elements("elements", trianglelist, numberoftriangles, 3, NULL, true), ElementAttributes("element_attributes", triangleattributelist, - numberoftriangles, 0, &Elements), + numberoftriangles, 0, &Elements, true), ElementVolumes("element_volumes", trianglearealist, - numberoftriangles, 1, &Elements), + numberoftriangles, 1, &Elements, true), Neighbors("neighbors", neighborlist, - numberoftriangles, 3, &Elements), + numberoftriangles, 3, &Elements, true), - Segments("segments", segmentlist, numberofsegments, 2), - SegmentMarkers("segment_markers", segmentmarkerlist, numberofsegments, 1, &Segments), + Segments("segments", segmentlist, numberofsegments, 2, NULL, true), + SegmentMarkers("segment_markers", segmentmarkerlist, numberofsegments, 1, &Segments, true), - Holes("holes", holelist, numberofholes, 2), + Holes("holes", holelist, numberofholes, 2, NULL, true), - Regions("regions", regionlist, numberofregions, 4), + Regions("regions", regionlist, numberofregions, 4, NULL, true), - Edges("edges", edgelist, numberofedges, 2), - EdgeMarkers("edge_markers", edgemarkerlist, numberofedges, 1, &Edges), - Normals("normals", normlist, numberofedges, 2, &Edges) + Edges("edges", edgelist, numberofedges, 2, NULL, true), + EdgeMarkers("edge_markers", edgemarkerlist, numberofedges, 1, &Edges, true), + Normals("normals", normlist, numberofedges, 2, &Edges, true) { numberofpointattributes = 0; numberofcorners = 3; @@ -88,7 +88,7 @@ struct tMeshDescriptor : public triangulateio, public boost::noncopyable numberoftriangleattributes = attrs; } - tMeshDescriptor &operator=(const tMeshDescriptor &src) + tMeshInfo &operator=(const tMeshInfo &src) { numberofpointattributes = src.numberofpointattributes ; numberofcorners = src.numberofcorners; @@ -121,9 +121,9 @@ struct tMeshDescriptor : public triangulateio, public boost::noncopyable -tMeshDescriptor *copyMesh(const tMeshDescriptor &src) +tMeshInfo *copyMesh(const tMeshInfo &src) { - auto_ptr copy(new tMeshDescriptor); + auto_ptr copy(new tMeshInfo); *copy = src; return copy.release(); } @@ -185,9 +185,9 @@ int triunsuitable(vertex triorg, vertex tridest, vertex triapex, REAL area) -void triangulateWrapper(char *options, tMeshDescriptor &in, - tMeshDescriptor &out, - tMeshDescriptor &voronoi, +void triangulateWrapper(char *options, tMeshInfo &in, + tMeshInfo &out, + tMeshInfo &voronoi, object refinement_func) { RefinementFunction = refinement_func; @@ -204,38 +204,14 @@ void triangulateWrapper(char *options, tMeshDescriptor &in, -template -void exposeForeignArray(T, const string &name) -{ - typedef tForeignArray cl; - - class_ - (name.c_str(), no_init) - .def("__len__", &cl::size) - .def("size", &cl::size) - .def("set_size", &cl::setSize) - .def("setup", &cl::setup) - .def("unit", &cl::unit) - .def("set", &cl::set) - .def("set_sub", &cl::setSub) - .def("get", &cl::get) - .def("get_sub", &cl::getSub) - .def("deallocate", &cl::deallocate) - ; -} - - - - - -BOOST_PYTHON_MODULE(internals) +BOOST_PYTHON_MODULE(_triangle) { def("triangulate", triangulateWrapper); { - typedef tMeshDescriptor cl; + typedef tMeshInfo cl; class_ - ("MeshDescriptor", init<>()) + ("MeshInfo", init<>()) .def_readonly("points", &cl::Points) .def_readonly("point_attributes", &cl::PointAttributes) .def_readonly("point_markers", &cl::PointMarkers) @@ -245,17 +221,17 @@ BOOST_PYTHON_MODULE(internals) .def_readonly("element_volumes", &cl::ElementVolumes) .def_readonly("neighbors", &cl::Neighbors) - .def_readonly("Segments", &cl::Segments) - .def_readonly("SegmentMarkers", &cl::SegmentMarkers) + .def_readonly("segments", &cl::Segments) + .def_readonly("segment_markers", &cl::SegmentMarkers) - .def_readonly("Holes", &cl::Holes) + .def_readonly("holes", &cl::Holes) - .def_readonly("Regions", &cl::Regions) + .def_readonly("regions", &cl::Regions) - .def_readonly("Edges", &cl::Edges) - .def_readonly("EdgeMarkers", &cl::EdgeMarkers) + .def_readonly("edges", &cl::Edges) + .def_readonly("edge_markers", &cl::EdgeMarkers) - .def_readonly("Normals", &cl::Normals) + .def_readonly("normals", &cl::Normals) .add_property("number_of_point_attributes", &cl::numberOfPointAttributes, @@ -269,8 +245,8 @@ BOOST_PYTHON_MODULE(internals) ; } - exposeForeignArray(REAL(), "RealArray"); - exposeForeignArray(int(), "IntArray"); + exposePODForeignArray("RealArray"); + exposePODForeignArray("IntArray"); class_, tVertex, boost::noncopyable>("Vertex", no_init) .add_property("x", &tVertex::x) diff --git a/test/test_triangle.py b/test/test_triangle.py new file mode 100644 index 0000000000000000000000000000000000000000..805d1a34f9b0e80df3f4d7f8e3f661a1af05b4b4 --- /dev/null +++ b/test/test_triangle.py @@ -0,0 +1,38 @@ +import meshpy.triangle as triangle +import math +import pickle + +segments = 50 + +points = [ (1,0),(1,1),(-1,1),(-1,-1),(1,-1),(1,0) ] + +for i in range( 0, segments + 1 ): + angle = i * 2 * math.pi / segments + points.append( ( 0.5 * math.cos( angle ), 0.5 * math.sin( angle ) ) ) + +def round_trip_connect(start, end): + result = [] + for i in range(start, end): + result.append((i, i+1)) + result.append((end, start)) + return result + +def needs_refinement( vert_origin, vert_destination, vert_apex, area ): + bary_x = (vert_origin.x + vert_destination.x + vert_apex.x) / 3 + bary_y = (vert_origin.y + vert_destination.y + vert_apex.y) / 3 + + dist_center = math.sqrt( bary_x**2 + bary_y**2 ) + max_area = math.fabs( 0.002 * (dist_center-0.5) ) + 0.0001 + return area > max_area + +info = triangle.MeshInfo() +info.set_points(points) +info.set_holes([(0,0)]) +info.set_segments(round_trip_connect(0, len(points)-1)) + +mesh = triangle.triangulate(info, refinement_func=needs_refinement) + +pickled = pickle.dumps(mesh) +mesh_2 = pickle.loads(pickled) + +triangle.write_gnuplot_mesh("triangles.dat", mesh_2)