revert to old, robust decimator (#1243)

.vscode/launch.json

@@ -12,7 +12,7 @@
       "program": "${workspaceFolder}/build/test/manifold_test",
       "args": [
         "--gtest_catch_exceptions=0",
-        "--gtest_filter=Samples.Scallop"
+        "--gtest_filter=Manifold.Simplify"
       ],
       "stopAtEntry": false,
       "cwd": "${workspaceFolder}/build/test",

src/face_op.cpp

@@ -80,13 +80,6 @@ PolygonsIdx ProjectPolygons(const std::vector<std::vector<int>>& polys,
   }  // for poly
   return polygons;
 }
-
-uint64_t PackEdgeFace(int edge, int face) {
-  return static_cast<uint64_t>(edge) << 32 | static_cast<uint64_t>(face);
-}
-
-int UnpackEdge(uint64_t data) { return static_cast<int>(data >> 32); }
-int UnpackFace(uint64_t data) { return static_cast<int>(data & 0xFFFFFFFF); }
 }  // namespace
 
 namespace manifold {
@@ -259,148 +252,6 @@ void Manifold::Impl::Face2Tri(const Vec<int>& faceEdge,
   CreateHalfedges(triProp, triVerts);
 }
 
-constexpr size_t kRemove = std::numeric_limits<size_t>::max();
-
-void Manifold::Impl::FlattenFaces() {
-  Vec<size_t> edgeTri(halfedge_.size());
-  const size_t numTri = NumTri();
-  const auto policy = autoPolicy(numTri);
-
-  std::vector<std::atomic<int>> vertDegree(NumVert());
-  for_each(policy, vertDegree.begin(), vertDegree.end(),
-           [](auto& v) { v.store(0); });
-
-  for_each_n(policy, countAt(0_uz), numTri,
-             [&edgeTri, &vertDegree, this](size_t tri) {
-               const auto& ref = meshRelation_.triRef[tri];
-               for (const int i : {0, 1, 2}) {
-                 const int edge = 3 * tri + i;
-                 const int pair = halfedge_[edge].pairedHalfedge;
-                 if (pair < 0) {
-                   edgeTri[edge] = kRemove;
-                   return;
-                 }
-                 if (ref.SameFace(meshRelation_.triRef[pair / 3])) {
-                   edgeTri[edge] = kRemove;
-                 } else {
-                   edgeTri[edge] = tri;
-                   ++vertDegree[halfedge_[edge].startVert];
-                 }
-               }
-             });
-
-  auto comp = [&edgeTri, this](size_t a, size_t b) {
-    if (edgeTri[a] == edgeTri[b]) return false;
-    if (edgeTri[a] == kRemove) return false;
-    if (edgeTri[b] == kRemove) return true;
-    const TriRef& refA = meshRelation_.triRef[edgeTri[a]];
-    const TriRef& refB = meshRelation_.triRef[edgeTri[b]];
-    return std::tie(refA.meshID, refA.coplanarID, refA.faceID) <
-           std::tie(refB.meshID, refB.coplanarID, refB.faceID);
-  };
-
-  Vec<size_t> newHalf2Old(halfedge_.size());
-  sequence(newHalf2Old.begin(), newHalf2Old.end());
-  stable_sort(newHalf2Old.begin(), newHalf2Old.end(), comp);
-  newHalf2Old.resize(std::lower_bound(countAt(0_uz), countAt(halfedge_.size()),
-                                      kRemove,
-                                      [&](const size_t i, const size_t val) {
-                                        return edgeTri[newHalf2Old[i]] < val;
-                                      }) -
-                     countAt(0_uz));
-
-  Vec<Halfedge> newHalfedge(newHalf2Old.size());
-  for_each_n(policy, countAt(0_uz), newHalf2Old.size(),
-             [&](size_t i) { newHalfedge[i] = halfedge_[newHalf2Old[i]]; });
-
-  Vec<int> faceEdge(1, 0);
-  for (size_t i = 1; i < newHalf2Old.size(); ++i) {
-    if (comp(newHalf2Old[i - 1], newHalf2Old[i])) {
-      faceEdge.push_back(i);
-    }
-  }
-  const size_t numFace = faceEdge.size();
-  faceEdge.push_back(newHalf2Old.size());
-
-  if (numFace < 4) {
-    MakeEmpty(Error::NoError);
-    return;  // empty
-  } else if (numFace == numTri) {
-    return;  // already flat
-  }
-
-  Vec<TriRef> halfedgeRef(halfedge_.size());
-  Vec<vec3> oldFaceNormal = std::move(faceNormal_);
-  faceNormal_.resize(numFace);
-
-  std::atomic<uint64_t> faceData(0);
-  Vec<int> faceEdge2(faceEdge.size());
-  faceEdge2[0] = 0;
-  Vec<Halfedge> newHalfedge2(newHalfedge.size());
-  for_each_n(policy, countAt(0_uz), numFace, [&](size_t inFace) {
-    std::vector<std::vector<int>> polys = AssembleHalfedges(
-        newHalfedge.cbegin() + faceEdge[inFace],
-        newHalfedge.cbegin() + faceEdge[inFace + 1], faceEdge[inFace]);
-    Vec<Halfedge> polys2;
-    for (const auto& poly : polys) {
-      int start = -1;
-      int last = -1;
-      for (const int edge : poly) {
-        // only keep vert touching 3 or more faces
-        if (vertDegree[newHalfedge[edge].startVert] < 3) continue;
-        if (start == -1) {
-          start = edge;
-        } else if (last != -1) {
-          polys2.push_back({newHalfedge[last].startVert,
-                            newHalfedge[edge].startVert, -1,
-                            newHalfedge[last].propVert});
-        }
-        last = edge;
-      }
-      polys2.push_back({newHalfedge[last].startVert,
-                        newHalfedge[start].startVert, -1,
-                        newHalfedge[last].propVert});
-    }
-
-    const int numEdge = polys2.size();
-    if (numEdge < 3) return;
-    const uint64_t inc = PackEdgeFace(numEdge, 1);
-    const uint64_t data = faceData.fetch_add(inc);
-    const int start = UnpackEdge(data);
-    const int outFace = UnpackFace(data);
-    if (numEdge > 0) {
-      std::copy(polys2.begin(), polys2.end(), newHalfedge2.begin() + start);
-    }
-    faceEdge2[outFace + 1] = start + numEdge;
-
-    const int oldFace = newHalf2Old[faceEdge[inFace]] / 3;
-    // only fill values that will get read.
-    halfedgeRef[start] = meshRelation_.triRef[oldFace];
-    faceNormal_[outFace] = oldFaceNormal[oldFace];
-  });
-
-  const int startFace = UnpackEdge(faceData);
-  const int outFace = UnpackFace(faceData);
-
-  if (startFace == 0) {
-    MakeEmpty(Error::NoError);
-    return;  // empty
-  }
-
-  newHalfedge2.resize(startFace);
-  halfedgeRef.resize(startFace);
-  faceEdge2.resize(outFace + 1);
-  halfedge_ = std::move(newHalfedge2);
-
-  // mark unreferenced verts for removal
-  for_each_n(policy, countAt(0_uz), NumVert(), [&](size_t vert) {
-    if (vertDegree[vert] < 3) vertPos_[vert] = vec3(NAN);
-  });
-
-  Face2Tri(faceEdge2, halfedgeRef, false);
-  Finish();
-}
-
 Polygons Manifold::Impl::Slice(double height) const {
   Box plane = bBox_;
   plane.min.z = plane.max.z = height;

src/impl.h

@@ -330,7 +330,6 @@ struct Manifold::Impl {
   // face_op.cpp
   void Face2Tri(const Vec<int>& faceEdge, const Vec<TriRef>& halfedgeRef,
                 bool allowConvex = false);
-  void FlattenFaces();
   Polygons Slice(double height) const;
   Polygons Project() const;
 

src/manifold.cpp

@@ -383,7 +383,8 @@ Manifold Manifold::SetTolerance(double tolerance) const {
   if (tolerance > impl->tolerance_) {
     impl->tolerance_ = tolerance;
     impl->MarkCoplanar();
-    impl->FlattenFaces();
+    impl->SimplifyTopology();
+    impl->Finish();
   } else {
     // for reducing tolerance, we need to make sure it is still at least
     // equal to epsilon.
@@ -405,7 +406,8 @@ Manifold Manifold::Simplify(double tolerance) const {
     impl->tolerance_ = tolerance;
     impl->MarkCoplanar();
   }
-  impl->FlattenFaces();
+  impl->SimplifyTopology();
+  impl->Finish();
   impl->tolerance_ = oldTolerance;
   return Manifold(impl);
 }

test/manifold_test.cpp

@@ -549,6 +549,19 @@ TEST(Manifold, SliceEmptyObject) {
   EXPECT_TRUE(empty.IsEmpty());
   CrossSection bottom = empty.Slice();
 }
+
+TEST(Manifold, Simplify) {
+  Polygons polyCircle =
+      CrossSection::Circle(1, 20).Translate({10, 0}).ToPolygons();
+  Manifold torus = Manifold::Revolve(polyCircle, 100);
+  Manifold simplified = torus.Simplify(0.4);
+  EXPECT_NEAR(torus.Volume(), simplified.Volume(), 20);
+  EXPECT_NEAR(torus.SurfaceArea(), simplified.SurfaceArea(), 10);
+
+#ifdef MANIFOLD_EXPORT
+  if (options.exportModels) ExportMesh("torus.glb", simplified.GetMeshGL(), {});
+#endif
+}
 #endif
 
 TEST(Manifold, MeshRelation) {

test/properties_test.cpp

@@ -77,7 +77,7 @@ TEST(Properties, Tolerance) {
 }
 
 TEST(Properties, ToleranceSphere) {
-  const int n = 500;
+  const int n = 1000;
   Manifold sphere = Manifold::Sphere(1, 4 * n);
   EXPECT_EQ(sphere.NumTri(), 8 * n * n);