diff --git a/pyproject.toml b/pyproject.toml index cdbd686..76fadda 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -4,7 +4,7 @@ build-backend = "scikit_build_core.build" [project] name = "openscad_cpp_evaluator" -version = "1.20.4" +version = "1.20.5" description = "C++ OpenSCAD evaluator with Python bindings" readme = "README.md" requires-python = ">=3.12" diff --git a/src/mesh_check.cpp b/src/mesh_check.cpp index 8edb983..fc8d907 100644 --- a/src/mesh_check.cpp +++ b/src/mesh_check.cpp @@ -79,20 +79,56 @@ size_t countUnwelded(const M& m) { // The faces around `v`, as (prev, next) pairs on the opposite edge. The link // is one cycle when following those pairs from any starting corner visits // all of them; more than one walk means the surface pinches at v. +// +// `vertTris` is CSR: the triangles around vertex v are +// vertTriList[vertTriStart[v] .. vertTriStart[v+1]). +// +// Every buffer here is hoisted out of the loop and cleared per vertex. The +// straightforward spelling -- an unordered_map of vectors for the link and a +// std::set for the walk, both built fresh per vertex -- allocated on the +// order of a dozen blocks for each of a mesh's vertices, millions of them on +// a 200K-triangle model, to answer a question about roughly six neighbours. +// Same reason the edge and face counts below collect flat and sort; see +// their comment. template size_t countPinched(const M& m, - const std::vector>& vertTris) { + const std::vector& vertTriStart, + const std::vector& vertTriList) { size_t pinched = 0; - for (Vert v = 0; v < vertTris.size(); ++v) { - const auto& tris = vertTris[v]; - if (tris.size() < 2) continue; - // next[a] = b for the corner opposite v in each face, undirected so - // a reversed neighbour still links up -- winding is checked - // separately and should not show up as a pinch too. - std::unordered_map> link; - for (size_t t : tris) { + const size_t nVerts = vertTriStart.empty() ? 0 : vertTriStart.size() - 1; + + std::vector nodes; // the distinct neighbours of v + std::vector parent; // union-find over indices into `nodes` + std::vector> edges; + + // Linear search, not a map: a vertex has a handful of neighbours, and at + // that size a scan beats hashing outright -- no allocation, no hashing, + // and the whole thing stays in cache. + const auto indexOf = [&](Vert w) -> uint32_t { + for (size_t i = 0; i < nodes.size(); ++i) { + if (nodes[i] == w) return static_cast(i); + } + nodes.push_back(w); + parent.push_back(parent.size()); + return static_cast(nodes.size() - 1); + }; + const auto find = [&](size_t i) { + while (parent[i] != i) { + parent[i] = parent[parent[i]]; // path halving + i = parent[i]; + } + return i; + }; + + for (Vert v = 0; v < nVerts; ++v) { + const size_t begin = vertTriStart[v], end = vertTriStart[v + 1]; + if (end - begin < 2) continue; + nodes.clear(); + parent.clear(); + edges.clear(); + for (size_t k = begin; k < end; ++k) { Vert w[3]; - triVerts(m, t, w); + triVerts(m, vertTriList[k], w); Vert a = 0, b = 0; int found = 0; for (int i = 0; i < 3; ++i) { @@ -100,21 +136,20 @@ size_t countPinched(const M& m, (found++ == 0 ? a : b) = w[i]; } if (found < 2) continue; // degenerate; counted elsewhere - link[a].push_back(b); - link[b].push_back(a); + edges.emplace_back(indexOf(a), indexOf(b)); } - if (link.empty()) continue; - std::set visited; - std::vector stack{link.begin()->first}; - while (!stack.empty()) { - Vert cur = stack.back(); - stack.pop_back(); - if (!visited.insert(cur).second) continue; - for (Vert nxt : link[cur]) { - if (!visited.count(nxt)) stack.push_back(nxt); - } + if (nodes.empty()) continue; + // Undirected, so a reversed neighbour still links up -- winding is + // checked separately and should not show up as a pinch too. + for (const auto& e : edges) { + const size_t ra = find(e.first), rb = find(e.second); + if (ra != rb) parent[ra] = rb; + } + size_t components = 0; + for (size_t i = 0; i < parent.size(); ++i) { + if (find(i) == i) ++components; } - if (visited.size() < link.size()) ++pinched; + if (components > 1) ++pinched; } return pinched; } @@ -158,10 +193,16 @@ MeshDiagnosis checkMesh(const M& mesh) { // index below; see its comment for the measurement. std::vector> uses; // {undirected edge, traversed a->b} uses.reserve(tris * 3); - std::vector> vertTris(nVerts); std::vector> seenFaces; // sorted-corner triples, deduplicated below seenFaces.reserve(tris); + // The triangles around each vertex, CSR rather than a vector per vertex: + // one inner vector per vertex is an allocation per vertex, and they are + // only ever read back as a flat run. Counted here, filled after the scan. + std::vector vertTriStart(nVerts + 1, 0); + std::vector vertTriList; + vertTriList.reserve(tris * 3); + for (size_t t = 0; t < tris; ++t) { if (degenerate(mesh, t)) ++d.degenerateFaces; // Counted, but NOT removed from the topology. A zero-area triangle @@ -181,12 +222,26 @@ MeshDiagnosis checkMesh(const M& mesh) { seenFaces.push_back(sorted); for (int i = 0; i < 3; ++i) { - if (v[i] < nVerts) vertTris[v[i]].push_back(t); + if (v[i] < nVerts) ++vertTriStart[v[i] + 1]; const Vert a = v[i], b = v[(i + 1) % 3]; uses.emplace_back(edgeKey(undirected(a, b)), a < b); } } + for (size_t i = 0; i < nVerts; ++i) vertTriStart[i + 1] += vertTriStart[i]; + vertTriList.resize(vertTriStart[nVerts]); + { + std::vector fill(vertTriStart.begin(), vertTriStart.end() - 1); + for (size_t t = 0; t < tris; ++t) { + if (repeatsAVertex(mesh, t)) continue; // skipped above too + Vert v[3]; + triVerts(mesh, t, v); + for (int i = 0; i < 3; ++i) { + if (v[i] < nVerts) vertTriList[fill[v[i]]++] = t; + } + } + } + // A face seen N times is N-1 duplicates, which is what insert().second // counted one at a time. std::sort(seenFaces.begin(), seenFaces.end()); @@ -214,7 +269,7 @@ MeshDiagnosis checkMesh(const M& mesh) { i = j; } - d.pinchedVertices = countPinched(mesh, vertTris); + d.pinchedVertices = countPinched(mesh, vertTriStart, vertTriList); d.unweldedVertices = countUnwelded(mesh); return d; } diff --git a/tests/test_mesh_check.cpp b/tests/test_mesh_check.cpp index b406594..aef9468 100644 --- a/tests/test_mesh_check.cpp +++ b/tests/test_mesh_check.cpp @@ -76,6 +76,26 @@ TEST(MeshCheck, TwoTetrahedraSharingOneVertexPinchThere) { EXPECT_FALSE(d.manifold()); } +TEST(MeshCheck, ThreeTetrahedraOnOneVertexAreStillOnePinchedVertex) { + // The link at vertex 0 falls into THREE pieces, not two. The count is of + // pinched VERTICES, so this is still 1 -- worth pinning because the walk + // only has to notice that the link is not one piece, and an + // implementation that counts pieces instead would say 2 here and agree + // with the two-tetrahedron case by luck. + manifold::MeshGL m = tetra(); + for (int k = 0; k < 2; ++k) { + const uint32_t base = static_cast(m.vertProperties.size() / 3); + const float z = -1.0f - static_cast(k) * 2.0f; + m.vertProperties.insert(m.vertProperties.end(), + {0, 0, z, 1, 0, z, 0, 1, z}); + m.triVerts.insert(m.triVerts.end(), { + 0u, base + 1, base, 0u, base, base + 2, + 0u, base + 2, base + 1, base, base + 1, base + 2}); + } + const MeshDiagnosis d = checkMesh(m); + EXPECT_EQ(d.pinchedVertices, 1u) << d.summary(); +} + TEST(MeshCheck, AReversedFaceIsCaughtAsInconsistentWinding) { manifold::MeshGL m = tetra(); std::swap(m.triVerts[1], m.triVerts[2]); // flip one face