From 7a2ce86eede1343144d8c20b9f0473d441bb3d92 Mon Sep 17 00:00:00 2001 From: wmayer Date: Sat, 3 Mar 2012 02:37:03 +0100 Subject: [PATCH] Script to create Menger sponge --- src/Mod/Mesh/App/Core/Degeneration.cpp | 24 +++++ src/Mod/Mesh/App/Core/Degeneration.h | 29 ++++++ src/Mod/Mesh/App/MeshPy.xml | 10 +++ src/Mod/Mesh/App/MeshPyImp.cpp | 30 +++++++ src/Mod/TemplatePyMod/MengerSponge.py | 120 +++++++++++++++++++++++++ 5 files changed, 213 insertions(+) create mode 100644 src/Mod/TemplatePyMod/MengerSponge.py diff --git a/src/Mod/Mesh/App/Core/Degeneration.cpp b/src/Mod/Mesh/App/Core/Degeneration.cpp index 4ec23e603..fd48104e6 100644 --- a/src/Mod/Mesh/App/Core/Degeneration.cpp +++ b/src/Mod/Mesh/App/Core/Degeneration.cpp @@ -380,6 +380,30 @@ bool MeshFixDuplicateFacets::Fixup() // ---------------------------------------------------------------------- +bool MeshEvalInternalFacets::Evaluate() +{ + _indices.clear(); + unsigned long uIndex=0; + const MeshFacetArray& rFaces = _rclMesh.GetFacets(); + + // get all facets + std::set aFaceSet; + MeshFacetArray::_TConstIterator first = rFaces.begin(); + for (MeshFacetArray::_TConstIterator it = rFaces.begin(); it != rFaces.end(); ++it, uIndex++) { + std::pair::iterator, bool> + pI = aFaceSet.insert(it); + if (!pI.second) { + // collect both elements + _indices.push_back(*pI.first - first); + _indices.push_back(uIndex); + } + } + + return _indices.empty(); +} + +// ---------------------------------------------------------------------- + bool MeshEvalDegeneratedFacets::Evaluate() { MeshFacetIterator it(_rclMesh); diff --git a/src/Mod/Mesh/App/Core/Degeneration.h b/src/Mod/Mesh/App/Core/Degeneration.h index 5a2ef182e..f2b28ffd7 100644 --- a/src/Mod/Mesh/App/Core/Degeneration.h +++ b/src/Mod/Mesh/App/Core/Degeneration.h @@ -188,6 +188,35 @@ public: bool Fixup (); }; +/** + * The MeshEvalInternalFacets class identifies internal facets of a volume mesh. + * @author Werner Mayer + */ +class MeshExport MeshEvalInternalFacets : public MeshEvaluation +{ +public: + /** + * Construction. + */ + MeshEvalInternalFacets (const MeshKernel &rclM) : MeshEvaluation( rclM ) { } + /** + * Destruction. + */ + ~MeshEvalInternalFacets () { } + /** + * Identifiy internal facets. + */ + bool Evaluate (); + /** + * Return the indices. + */ + const std::vector& GetIndices() const + { return _indices; } + +private: + std::vector _indices; +}; + /** * The MeshEvalDegeneratedFacets class searches for degenerated facets. A facet is degenerated either if its points * are collinear, i.e. they lie on a line or two points are coincident. In the latter case these points are duplicated. diff --git a/src/Mod/Mesh/App/MeshPy.xml b/src/Mod/Mesh/App/MeshPy.xml index 9e2eac8ec..6070d9fc1 100644 --- a/src/Mod/Mesh/App/MeshPy.xml +++ b/src/Mod/Mesh/App/MeshPy.xml @@ -127,6 +127,16 @@ Example: Remove a list of facet indices from the mesh + + + Builds a list of facet indices with triangles that are inside a volume mesh + + + + + Repairs the neighbourhood which might be broken + + Combine this mesh with another mesh. diff --git a/src/Mod/Mesh/App/MeshPyImp.cpp b/src/Mod/Mesh/App/MeshPyImp.cpp index 8c15864f0..ace103289 100644 --- a/src/Mod/Mesh/App/MeshPyImp.cpp +++ b/src/Mod/Mesh/App/MeshPyImp.cpp @@ -36,6 +36,7 @@ #include "MeshProperties.h" #include "Core/Algorithm.h" #include "Core/Iterator.h" +#include "Core/Degeneration.h" #include "Core/Elements.h" #include "Core/Grid.h" #include "Core/MeshKernel.h" @@ -581,6 +582,35 @@ PyObject* MeshPy::removeFacets(PyObject *args) Py_Return; } +PyObject* MeshPy::getInternalFacets(PyObject *args) +{ + if (!PyArg_ParseTuple(args, "")) + return 0; + + const MeshCore::MeshKernel& kernel = getMeshObjectPtr()->getKernel(); + MeshCore::MeshEvalInternalFacets eval(kernel); + eval.Evaluate(); + + const std::vector& indices = eval.GetIndices(); + Py::List ary(indices.size()); + Py::List::size_type pos=0; + for (std::vector::const_iterator it = indices.begin(); it != indices.end(); ++it) { + ary[pos++] = Py::Long(*it); + } + + return Py::new_reference_to(ary); +} + +PyObject* MeshPy::rebuildNeighbourHood(PyObject *args) +{ + if (!PyArg_ParseTuple(args, "")) + return 0; + + MeshCore::MeshKernel& kernel = getMeshObjectPtr()->getKernel(); + kernel.RebuildNeighbours(); + Py_Return; +} + PyObject* MeshPy::addMesh(PyObject *args) { PyObject* mesh; diff --git a/src/Mod/TemplatePyMod/MengerSponge.py b/src/Mod/TemplatePyMod/MengerSponge.py new file mode 100644 index 000000000..67beb5921 --- /dev/null +++ b/src/Mod/TemplatePyMod/MengerSponge.py @@ -0,0 +1,120 @@ +# Script to create a Menger sponge +# (c) 2012 Werner Mayer LGPL + +# The script is based on the work of daxmick at +# https://sourceforge.net/apps/phpbb/free-cad/viewtopic.php?f=3&t=2307 + +import threading +import Mesh, MeshGui +from FreeCAD import Base + +# Create a Box and Place it a coords (x,y,z) +def PlaceBox(x,y,z): + box = Mesh.createBox(1,1,1) + box.translate(x,y,z) + return box + +def Sierpinski(level,x0,y0,z0): + #print threading.current_thread().name + boxnums = pow(3,level) + thirds = boxnums / 3 + twothirds = thirds * 2 + if(level == 0): + rangerx = [x0] + rangery = [y0] + rangerz = [z0] + else: + rangerx = [ x0, x0 + thirds, x0 + twothirds ] + rangery = [ y0, y0 + thirds, y0 + twothirds ] + rangerz = [ z0, z0 + thirds, z0 + twothirds ] + block = 1 + skip=[5,11,13,14,15,17,23] + mesh=Mesh.Mesh() + for i in rangerx: + for j in rangery: + for k in rangerz: + if block not in skip: + if(level > 0): + mesh.addMesh(Sierpinski(level-1,i,j,k)) + else: + mesh.addMesh(PlaceBox(i,j,k)) + block+=1 + return mesh + +### Multi-threaded ### + +class MengerThread(threading.Thread): + def __init__(self,args): + self.args=args + self.mesh=Mesh.Mesh() + threading.Thread.__init__(self) + def run(self): + for i in self.args: + self.mesh.addMesh(Sierpinski(*i)) + +def makeMengerSponge_mt(level=3,x0=0,y0=0,z0=0): + """ + Is much slower than makeMengerSponge!!! :( + """ + if level == 0: + mesh=Sierpinski(level,x0,y0,z0) + Mesh.show(mesh) + return + + boxnums = pow(3,level) + thirds = boxnums / 3 + twothirds = thirds * 2 + + rangerx = [ x0, x0 + thirds, x0 + twothirds ] + rangery = [ y0, y0 + thirds, y0 + twothirds ] + rangerz = [ z0, z0 + thirds, z0 + twothirds ] + block = 1 + skip=[5,11,13,14,15,17,23] + + # collect the arguments for the algorithm in a list + args=[] + for i in rangerx: + for j in rangery: + for k in rangerz: + if block not in skip: + args.append((level-1,i,j,k)) + block+=1 + + numJobs = 4 + threads=[] + while numJobs > 0: + size = len(args) + count = size / numJobs + numJobs-=1 + thr=MengerThread(args[:count]) + threads.append(thr) + args=args[count:] + + print "Number of threads: %i" % (len(threads)) + for thr in threads: + thr.start() + for thr in threads: + thr.join() + + mesh=Mesh.Mesh() + for thr in threads: + mesh.addMesh(thr.mesh) + del thr.mesh + + print mesh + mesh.removeDuplicatedPoints() + mesh.removeFacets(mesh.getInternalFacets()) + mesh.rebuildNeighbourHood() + print "Mesh is solid: %s" % (mesh.isSolid()) + Mesh.show(mesh) + + +### Single-threaded ### + +def makeMengerSponge(level=3,x0=0,y0=0,z0=0): + mesh=Sierpinski(level,x0,y0,z0) + mesh.removeDuplicatedPoints() + mesh.removeFacets(mesh.getInternalFacets()) + mesh.rebuildNeighbourHood() + print "Mesh is solid: %s" % (mesh.isSolid()) + Mesh.show(mesh)