fixes

* please gcc

* added documentation for mallocManaged

.github/workflows/manifold.yml

@@ -220,7 +220,7 @@ jobs:
     timeout-minutes: 30
     strategy:
       matrix:
-        variant: [none, omp, tbb, none-cuda, omp-cuda, tbb-cuda]
+        variant: [none, omp, tbb, none-cuda, omp-cuda, tbb-cuda, js]
     runs-on: ubuntu-latest
     if: github.event.pull_request.draft == false
     steps:

.gitignore

@@ -109,6 +109,9 @@ build
 buildWASM
 docs/html
 node_modules/
+bindings/wasm/examples/built
+bindings/wasm/examples/dist
+bindings/wasm/examples/public/manifold*
 __pycache__
 .vscode/c_cpp_properties.json
 *.sublime-project

.gitmodules

@@ -14,3 +14,6 @@
 [submodule "src/third_party/clipper2"]
 	path = src/third_party/clipper2
 	url = https://github.com/AngusJohnson/clipper2
+[submodule "src/third_party/quickhull"]
+	path = src/third_party/quickhull
+	url = https://github.com/akuukka/quickhull

.vscode/settings.json

@@ -98,7 +98,17 @@
     "concepts": "cpp",
     "queue": "cpp",
     "stack": "cpp",
-    "__hash_table": "cpp"
+    "__hash_table": "cpp",
+    "__std_stream": "cpp",
+    "__verbose_abort": "cpp",
+    "any": "cpp",
+    "barrier": "cpp",
+    "charconv": "cpp",
+    "coroutine": "cpp",
+    "latch": "cpp",
+    "semaphore": "cpp",
+    "span": "cpp",
+    "__sso_allocator": "cpp"
   },
   "C_Cpp.clang_format_fallbackStyle": "google",
   "editor.formatOnSave": true,

README.md

@@ -98,7 +98,7 @@ node test/manifold_test.js
 
 ### Python
 
-The CMake script will build the python binding `pymanifold` automatically. To
+The CMake script will build the python binding `manifold3d` automatically. To
 use the extension, please add `$BUILD_DIR/tools` to your `PYTHONPATH`, where
 `$BUILD_DIR` is the build directory for CMake. Examples using the python binding
 can be found in `bindings/python/examples`. To see exported samples, run:
@@ -112,8 +112,8 @@ Run the following code in the interpreter for
 python binding documentation:
 
 ```
->>> import pymanifold
->>> help(pymanifold)
+>>> import manifold3d
+>>> help(manifold3d)
 ```
 
 For more detailed documentation, please refer to the C++ API.

bindings/c/cross.cpp

@@ -106,13 +106,6 @@ ManifoldCrossSection *manifold_cross_section_intersection(
   return to_c(new (mem) CrossSection(cs));
 }
 
-ManifoldCrossSection *manifold_cross_section_rect_clip(void *mem,
-                                                       ManifoldCrossSection *cs,
-                                                       ManifoldRect *r) {
-  auto clipped = from_c(cs)->RectClip(*from_c(r));
-  return to_c(new (mem) CrossSection(clipped));
-}
-
 ManifoldCrossSection *manifold_cross_section_translate(void *mem,
                                                        ManifoldCrossSection *cs,
                                                        float x, float y) {

bindings/c/include/manifoldc.h

@@ -2,9 +2,8 @@
 #include <stdint.h>
 #include <types.h>
 
-#include <cstdint>
-
 #ifdef __cplusplus
+#include <cstdint>
 extern "C" {
 #endif
 
@@ -50,6 +49,15 @@ ManifoldMeshGL *manifold_level_set_seq(void *mem,
                                        float (*sdf)(float, float, float),
                                        ManifoldBox *bounds, float edge_length,
                                        float level);
+// The _context variants of manifold_level_set allow a pointer to be passed
+// back to each invocation of the sdf function pointer, for languages that
+// need additional data.
+ManifoldMeshGL *manifold_level_set_context(
+    void *mem, float (*sdf)(float, float, float, void *), ManifoldBox *bounds,
+    float edge_length, float level, void *ctx);
+ManifoldMeshGL *manifold_level_set_seq_context(
+    void *mem, float (*sdf)(float, float, float, void *), ManifoldBox *bounds,
+    float edge_length, float level, void *ctx);
 
 // Manifold Vectors
 
@@ -195,9 +203,6 @@ ManifoldCrossSection *manifold_cross_section_difference(
     void *mem, ManifoldCrossSection *a, ManifoldCrossSection *b);
 ManifoldCrossSection *manifold_cross_section_intersection(
     void *mem, ManifoldCrossSection *a, ManifoldCrossSection *b);
-ManifoldCrossSection *manifold_cross_section_rect_clip(void *mem,
-                                                       ManifoldCrossSection *cs,
-                                                       ManifoldRect *r);
 
 // CrossSection Transformation
 

bindings/c/manifoldc.cpp

@@ -33,6 +33,20 @@ ManifoldMeshGL *level_set(void *mem, float (*sdf)(float, float, float),
   auto mesh = LevelSet(fun, *from_c(bounds), edge_length, level, pol);
   return to_c(new (mem) MeshGL(mesh));
 }
+ManifoldMeshGL *level_set_context(
+    void *mem, float (*sdf_context)(float, float, float, void *),
+    ManifoldBox *bounds, float edge_length, float level, bool seq, void *ctx) {
+  // Bind function with context argument to one without
+  using namespace std::placeholders;
+  std::function<float(float, float, float)> sdf =
+      std::bind(sdf_context, _1, _2, _3, ctx);
+  std::function<float(glm::vec3)> fun = [sdf](glm::vec3 v) {
+    return (sdf(v.x, v.y, v.z));
+  };
+  auto pol = seq ? std::make_optional(ExecutionPolicy::Seq) : std::nullopt;
+  auto mesh = LevelSet(fun, *from_c(bounds), edge_length, level, pol);
+  return to_c(new (mem) MeshGL(mesh));
+}
 }  // namespace
 
 #ifdef __cplusplus
@@ -233,6 +247,18 @@ ManifoldMeshGL *manifold_level_set_seq(void *mem,
   return level_set(mem, sdf, bounds, edge_length, level, true);
 }
 
+ManifoldMeshGL *manifold_level_set_context(
+    void *mem, float (*sdf)(float, float, float, void *), ManifoldBox *bounds,
+    float edge_length, float level, void *ctx) {
+  return level_set_context(mem, sdf, bounds, edge_length, level, false, ctx);
+}
+
+ManifoldMeshGL *manifold_level_set_seq_context(
+    void *mem, float (*sdf)(float, float, float, void *), ManifoldBox *bounds,
+    float edge_length, float level, void *ctx) {
+  return level_set_context(mem, sdf, bounds, edge_length, level, true, ctx);
+}
+
 ManifoldManifold *manifold_refine(void *mem, ManifoldManifold *m, int refine) {
   auto refined = from_c(m)->Refine(refine);
   return to_c(new (mem) Manifold(refined));

bindings/python/CMakeLists.txt

@@ -15,10 +15,10 @@
 project(python)
 
 add_subdirectory(third_party)
-pybind11_add_module(pymanifold pymanifold.cpp)
-target_link_libraries(pymanifold PRIVATE manifold)
-target_compile_options(pymanifold PRIVATE ${MANIFOLD_FLAGS})
-target_compile_features(pymanifold PUBLIC cxx_std_17)
-target_include_directories(pymanifold
+pybind11_add_module(manifold3d manifold3d.cpp)
+target_link_libraries(manifold3d PRIVATE manifold)
+target_compile_options(manifold3d PRIVATE ${MANIFOLD_FLAGS})
+target_compile_features(manifold3d PUBLIC cxx_std_17)
+target_include_directories(manifold3d
     PRIVATE ${PYBIND11_DIR}/include
 )

bindings/python/examples/bricks.py

@@ -14,7 +14,7 @@
  limitations under the License.
  """
 
-from pymanifold import Manifold
+from manifold3d import Manifold
 
 # https://gist.github.com/deckar01/ef11def51de7e71d9f288c6e5819fdb7
 

bindings/python/examples/cube_with_dents.py

@@ -14,7 +14,7 @@
  limitations under the License.
  """
  
-from pymanifold import Manifold
+from manifold3d import Manifold
 from functools import reduce
 
 # https://gist.github.com/ochafik/2db96400e3c1f73558fcede990b8a355#file-cube-with-half-spheres-dents-scad

bindings/python/examples/extrude.py

@@ -1,4 +1,4 @@
-from pymanifold import CrossSection, FillRule
+from manifold3d import CrossSection, FillRule
 
 
 def run():

bindings/python/examples/maze.py

@@ -1,7 +1,7 @@
 # One of the benchmarks in in https://gist.github.com/ochafik/2db96400e3c1f73558fcede990b8a355
 # https://www.thingiverse.com/groups/openscad/forums/general/topic:34699
 
-from pymanifold import Manifold
+from manifold3d import Manifold
 import numpy as np
 
 n = 10

bindings/python/examples/run_all.py

@@ -18,6 +18,7 @@ import pathlib
 import sys
 import importlib
 import trimesh
+import numpy as np
 from time import time
 
 if __name__ == "__main__":
@@ -34,8 +35,11 @@ if __name__ == "__main__":
         model = module.run()
         mesh = model.to_mesh()
         if export_models:
-            meshOut = trimesh.Trimesh(
-                vertices=mesh.vert_pos, faces=mesh.tri_verts)
+            if mesh.vert_properties.shape[1] > 3:
+                vertices = np.hsplit(mesh.vert_properties, 3)[0]
+            else:
+                vertices = mesh.vert_properties
+            meshOut = trimesh.Trimesh(vertices=vertices, faces=mesh.tri_verts)
             trimesh.exchange.export.export_mesh(meshOut, f'{f}.glb', 'glb')
             print(f'Exported model to {f}.glb')
         t1 = time()

bindings/python/examples/union_failure.py

@@ -1,4 +1,4 @@
-from pymanifold import Manifold
+from manifold3d import Manifold
 
 def run():
     # for some reason this causes collider error

bindings/python/manifold3d.cpp

@@ -27,12 +27,19 @@ using namespace manifold;
 typedef std::tuple<float, float> Float2;
 typedef std::tuple<float, float, float> Float3;
 
-PYBIND11_MODULE(pymanifold, m) {
-  m.doc() =
-      "Python binding for the manifold library. Please check the C++ "
-      "documentation for APIs.\n"
-      "This binding will perform copying to make the API more familiar to "
-      "OpenSCAD users.";
+constexpr auto pyArrayFlags = py::array::c_style | py::array::forcecast;
+
+PYBIND11_MAKE_OPAQUE(std::vector<int>);
+PYBIND11_MAKE_OPAQUE(std::vector<uint32_t>);
+PYBIND11_MAKE_OPAQUE(std::vector<float>);
+
+template <typename T>
+std::vector<T> toVector(const py::array_t<T> &arr) {
+  return std::vector<T>(arr.data(), arr.data() + arr.size());
+}
+
+PYBIND11_MODULE(manifold3d, m) {
+  m.doc() = "Python binding for the Manifold library.";
 
   m.def("set_min_circular_angle", Quality::SetMinCircularAngle,
         "Sets an angle constraint the default number of circular segments for "
@@ -40,7 +47,7 @@ PYBIND11_MODULE(pymanifold, m) {
         "and Manifold::Revolve() constructors. The number of segments will be "
         "rounded up to the nearest factor of four."
         "\n\n"
-        "@param angle The minimum angle in degrees between consecutive "
+        ":param angle: The minimum angle in degrees between consecutive "
         "segments. The angle will increase if the the segments hit the minimum "
         "edge length.\n"
         "Default is 10 degrees.");
@@ -51,7 +58,7 @@ PYBIND11_MODULE(pymanifold, m) {
         "and Manifold::Revolve() constructors. The number of segments will be "
         "rounded up to the nearest factor of four."
         "\n\n"
-        "@param length The minimum length of segments. The length will "
+        ":param length: The minimum length of segments. The length will "
         "increase if the the segments hit the minimum angle. Default is 1.0.");
 
   m.def("set_circular_segments", Quality::SetCircularSegments,
@@ -60,15 +67,15 @@ PYBIND11_MODULE(pymanifold, m) {
         "Manifold::Revolve() constructors. Overrides the edge length and angle "
         "constraints and sets the number of segments to exactly this value."
         "\n\n"
-        "@param number Number of circular segments. Default is 0, meaning no "
+        ":param number: Number of circular segments. Default is 0, meaning no "
         "constraint is applied.");
 
-  m.def(
-      "get_circular_segments", Quality::GetCircularSegments,
-      "Determine the result of the SetMinCircularAngle(), "
-      "SetMinCircularEdgeLength(), and SetCircularSegments() defaults."
-      "\n\n"
-      "@param radius For a given radius of circle, determine how many default");
+  m.def("get_circular_segments", Quality::GetCircularSegments,
+        "Determine the result of the SetMinCircularAngle(), "
+        "SetMinCircularEdgeLength(), and SetCircularSegments() defaults."
+        "\n\n"
+        ":param radius: For a given radius of circle, determine how many "
+        "default");
 
   py::class_<Manifold>(m, "Manifold")
       .def(py::init<>())
@@ -93,8 +100,8 @@ PYBIND11_MODULE(pymanifold, m) {
           "transform",
           [](Manifold &self, py::array_t<float> &mat) {
             auto mat_view = mat.unchecked<2>();
-            if (mat_view.shape(0) != 3 || mat_view.shape(1) != 4)
-              throw std::runtime_error("Invalid matrix shape");
+            if (mat.ndim() != 2 || mat.shape(0) != 3 || mat.shape(1) != 4)
+              throw std::runtime_error("Invalid matrix shape, expected (3, 4)");
             glm::mat4x3 mat_glm;
             for (int i = 0; i < 3; i++) {
               for (int j = 0; j < 4; j++) {
@@ -126,7 +133,7 @@ PYBIND11_MODULE(pymanifold, m) {
           "translate",
           [](Manifold &self, py::array_t<float> &t) {
             auto t_view = t.unchecked<1>();
-            if (t.shape(0) != 3)
+            if (t.ndim() != 1 || t.shape(0) != 3)
               throw std::runtime_error("Invalid vector shape");
             return self.Translate(glm::vec3(t_view(0), t_view(1), t_view(2)));
           },
@@ -135,23 +142,23 @@ PYBIND11_MODULE(pymanifold, m) {
           "Transforms are combined and applied lazily."
           "\n\n"
           ":param v: The vector to add to every vertex.")
-      .def("scale",
-           static_cast<Manifold (*)(Manifold &, float)>(
-               [](Manifold &self, float scale) {
-                 return self.Scale(glm::vec3(scale));
-               }),
-           py::arg("scale"),
-           "Scale this Manifold in space. This operation can be chained. "
-           "Transforms are combined and applied lazily."
-           "\n\n"
-           ":param scale: The scalar multiplier for each component of every "
-           "vertices.")
+      .def(
+          "scale",
+          [](Manifold &self, float scale) {
+            return self.Scale(glm::vec3(scale));
+          },
+          py::arg("scale"),
+          "Scale this Manifold in space. This operation can be chained. "
+          "Transforms are combined and applied lazily."
+          "\n\n"
+          ":param scale: The scalar multiplier for each component of every "
+          "vertices.")
       .def(
           "scale",
           [](Manifold &self, py::array_t<float> &scale) {
             auto scale_view = scale.unchecked<1>();
-            if (scale_view.shape(0) != 3)
-              throw std::runtime_error("Invalid vector shape");
+            if (scale.ndim() != 1 || scale.shape(0) != 3)
+              throw std::runtime_error("Invalid vector shape, expected (3)");
             glm::vec3 v(scale_view(0), scale_view(1), scale_view(2));
             return self.Scale(v);
           },
@@ -160,12 +167,26 @@ PYBIND11_MODULE(pymanifold, m) {
           "Transforms are combined and applied lazily."
           "\n\n"
           ":param v: The vector to multiply every vertex by component.")
+      .def(
+          "mirror",
+          [](Manifold &self, py::array_t<float> &mirror) {
+            auto v_view = mirror.unchecked<1>();
+            if (mirror.ndim() != 1 || mirror.shape(0) != 3)
+              throw std::runtime_error("Invalid vector shape, expected (3)");
+            glm::vec3 v(v_view(0), v_view(1), v_view(2));
+            return self.Mirror(v);
+          },
+          py::arg("v"),
+          "Mirror this Manifold in space. This operation can be chained. "
+          "Transforms are combined and applied lazily."
+          "\n\n"
+          ":param mirror: The vector defining the axis of mirroring.")
       .def(
           "rotate",
           [](Manifold &self, py::array_t<float> &v) {
             auto v_view = v.unchecked<1>();
-            if (v_view.shape(0) != 3)
-              throw std::runtime_error("Invalid vector shape");
+            if (v.ndim() != 1 || v.shape(0) != 3)
+              throw std::runtime_error("Invalid vector shape, expected (3)");
             return self.Rotate(v_view(0), v_view(1), v_view(2));
           },
           py::arg("v"),
@@ -207,7 +228,14 @@ PYBIND11_MODULE(pymanifold, m) {
               v.z = std::get<2>(fv);
             });
           },
-          py::arg("f"))
+          py::arg("f"),
+          "This function does not change the topology, but allows the vertices "
+          "to be moved according to any arbitrary input function. It is easy "
+          "to create a function that warps a geometrically valid object into "
+          "one which overlaps, but that is not checked here, so it is up to "
+          "the user to choose their function with discretion."
+          "\n\n"
+          ":param warpFunc: A function that modifies a given vertex position.")
       .def(
           "refine", [](Manifold &self, int n) { return self.Refine(n); },
           py::arg("n"),
@@ -221,7 +249,34 @@ PYBIND11_MODULE(pymanifold, m) {
           "\n"
           ":param n: The number of pieces to split every edge into. Must be > "
           "1.")
-      .def("to_mesh", &Manifold::GetMesh)
+      .def(
+          "to_mesh",
+          [](Manifold &self, std::optional<py::array_t<uint32_t>> &normalIdx) {
+            glm::ivec3 v(0);
+            if (normalIdx.has_value()) {
+              if (normalIdx.value().ndim() != 1 ||
+                  normalIdx.value().shape(0) != 3)
+                throw std::runtime_error("Invalid vector shape, expected (3)");
+              auto normalIdx_view = normalIdx.value().unchecked<1>();
+              v = glm::ivec3(normalIdx_view(0), normalIdx_view(1),
+                             normalIdx_view(2));
+            }
+            return self.GetMeshGL(v);
+          },
+          "The most complete output of this library, returning a MeshGL that "
+          "is designed to easily push into a renderer, including all "
+          "interleaved vertex properties that may have been input. It also "
+          "includes relations to all the input meshes that form a part of "
+          "this result and the transforms applied to each."
+          "\n\n"
+          ":param normalIdx: If the original MeshGL inputs that formed this "
+          "manifold had properties corresponding to normal vectors, you can "
+          "specify which property channels these are (x, y, z), which will "
+          "cause this output MeshGL to automatically update these normals "
+          "according to the applied transforms and front/back side. Each "
+          "channel must be >= 3 and < numProp, and all original MeshGLs must "
+          "use the same channels for their normals.",
+          py::arg("normalIdx") = py::none())
       .def("num_vert", &Manifold::NumVert,
            "The number of vertices in the Manifold.")
       .def("num_edge", &Manifold::NumEdge,
@@ -272,7 +327,7 @@ PYBIND11_MODULE(pymanifold, m) {
            "Does the Manifold have any triangles?")
       .def(
           "decompose", [](Manifold &self) { return self.Decompose(); },
-          " This operation returns a vector of Manifolds that are "
+          "This operation returns a vector of Manifolds that are "
           "topologically disconnected. If everything is connected, the vector "
           "is length one, containing a copy of the original. It is the inverse "
           "operation of Compose().")
@@ -290,10 +345,10 @@ PYBIND11_MODULE(pymanifold, m) {
           py::arg("normal"), py::arg("origin_offset"),
           "Convenient version of Split() for a half-space."
           "\n\n"
-          "@param normal This vector is normal to the cutting plane and its "
+          ":param normal: This vector is normal to the cutting plane and its "
           "length does not matter. The first result is in the direction of "
           "this vector, the second result is on the opposite side.\n"
-          "@param originOffset The distance of the plane from the origin in "
+          ":param originOffset: The distance of the plane from the origin in "
           "the direction of the normal vector.")
       .def(
           "trim_by_plane",
@@ -306,10 +361,10 @@ PYBIND11_MODULE(pymanifold, m) {
           "Identical to SplitByPlane(), but calculating and returning only the "
           "first result."
           "\n\n"
-          "@param normal This vector is normal to the cutting plane and its "
+          ":param normal: This vector is normal to the cutting plane and its "
           "length does not matter. The result is in the direction of this "
           "vector from the plane.\n"
-          "@param originOffset The distance of the plane from the origin in "
+          ":param originOffset: The distance of the plane from the origin in "
           "the direction of the normal vector.")
       .def_property_readonly(
           "bounding_box",
@@ -353,7 +408,7 @@ PYBIND11_MODULE(pymanifold, m) {
           "by sharping all edges that are incident on it, allowing cones to be "
           "formed.")
       .def_static(
-          "from_mesh", [](const Mesh &mesh) { return Manifold(mesh); },
+          "from_mesh", [](const MeshGL &mesh) { return Manifold(mesh); },
           py::arg("mesh"))
       .def_static(
           "tetrahedron", []() { return Manifold::Tetrahedron(); },
@@ -378,8 +433,8 @@ PYBIND11_MODULE(pymanifold, m) {
           "cube",
           [](py::array_t<float> &size, bool center = false) {
             auto size_view = size.unchecked<1>();
-            if (size_view.shape(0) != 3)
-              throw std::runtime_error("Invalid vector shape");
+            if (size.ndim() != 1 || size.shape(0) != 3)
+              throw std::runtime_error("Invalid vector shape, expected (3)");
             return Manifold::Cube(
                 glm::vec3(size_view(0), size_view(1), size_view(2)), center);
           },
@@ -442,101 +497,116 @@ PYBIND11_MODULE(pymanifold, m) {
                   "marking sets of triangles that can be looked up after "
                   "further operations. Assign to MeshGL.runOriginalID vector");
 
-  py::class_<Mesh>(m, "Mesh")
-      .def(py::init([](py::array_t<float> &vertPos, py::array_t<int> &triVerts,
-                       py::array_t<float> &vertNormal,
-                       py::array_t<float> &halfedgeTangent) {
-             auto vertPos_view = vertPos.unchecked<2>();
-             auto triVerts_view = triVerts.unchecked<2>();
-             auto vertNormal_view = vertNormal.unchecked<2>();
-             auto halfedgeTangent_view = halfedgeTangent.unchecked<2>();
-             if (vertPos_view.shape(1) != 3)
-               throw std::runtime_error("Invalid vert_pos shape");
-             if (triVerts_view.shape(1) != 3)
-               throw std::runtime_error("Invalid tri_verts shape");
-             if (vertNormal_view.shape(0) != 0) {
-               if (vertNormal_view.shape(1) != 3)
-                 throw std::runtime_error("Invalid vert_normal shape");
-               if (vertNormal_view.shape(0) != vertPos_view.shape(0))
-                 throw std::runtime_error(
-                     "vert_normal must have the same length as vert_pos");
-             }
-             if (halfedgeTangent_view.shape(0) != 0) {
-               if (halfedgeTangent_view.shape(1) != 4)
-                 throw std::runtime_error("Invalid halfedge_tangent shape");
-               if (halfedgeTangent_view.shape(0) != triVerts_view.shape(0) * 3)
-                 throw std::runtime_error(
-                     "halfedge_tangent must be three times as long as "
-                     "tri_verts");
-             }
-             std::vector<glm::vec3> vertPos_vec(vertPos_view.shape(0));
-             std::vector<glm::ivec3> triVerts_vec(triVerts_view.shape(0));
-             std::vector<glm::vec3> vertNormal_vec(vertNormal_view.shape(0));
-             std::vector<glm::vec4> halfedgeTangent_vec(
-                 halfedgeTangent_view.shape(0));
-             for (int i = 0; i < vertPos_view.shape(0); i++)
-               for (const int j : {0, 1, 2})
-                 vertPos_vec[i][j] = vertPos_view(i, j);
-             for (int i = 0; i < triVerts_view.shape(0); i++)
-               for (const int j : {0, 1, 2})
-                 triVerts_vec[i][j] = triVerts_view(i, j);
-             for (int i = 0; i < vertNormal_view.shape(0); i++)
-               for (const int j : {0, 1, 2})
-                 vertNormal_vec[i][j] = vertNormal_view(i, j);
-             for (int i = 0; i < halfedgeTangent_view.shape(0); i++)
-               for (const int j : {0, 1, 2, 3})
-                 halfedgeTangent_vec[i][j] = halfedgeTangent_view(i, j);
-             return Mesh({vertPos_vec, triVerts_vec, vertNormal_vec,
-                          halfedgeTangent_vec});
-           }),
-           py::arg("vert_pos"), py::arg("tri_verts"), py::arg("vert_normal"),
-           py::arg("halfedge_tangent"))
-      .def_property_readonly("vert_pos",
-                             [](Mesh &self) {
-                               const int numVert = self.vertPos.size();
-                               py::array_t<float> vert_pos({numVert, 3});
-                               auto vert_pos_view =
-                                   vert_pos.mutable_unchecked<2>();
+  py::class_<MeshGL>(m, "Mesh")
+      .def(
+          // note that reshape requires mutable array_t, but this will not
+          // affect the original array passed into the function
+          py::init(
+              [](py::array_t<float, pyArrayFlags> &vertProp,
+                 py::array_t<int, pyArrayFlags> &triVerts,
+                 const std::optional<py::array_t<uint32_t, pyArrayFlags>>
+                     &mergeFromVert,
+                 const std::optional<py::array_t<uint32_t, pyArrayFlags>>
+                     &mergeToVert,
+                 const std::optional<py::array_t<uint32_t, pyArrayFlags>>
+                     &runIndex,
+                 const std::optional<py::array_t<uint32_t, pyArrayFlags>>
+                     &runOriginalID,
+                 std::optional<py::array_t<float, pyArrayFlags>> &runTransform,
+                 const std::optional<py::array_t<uint32_t, pyArrayFlags>>
+                     &faceID,
+                 const std::optional<py::array_t<float, pyArrayFlags>>
+                     &halfedgeTangent,
+                 float precision) {
+                MeshGL out;
 
-                               for (int i = 0; i < numVert; ++i)
-                                 for (const int j : {0, 1, 2})
-                                   vert_pos_view(i, j) = self.vertPos[i][j];
-                               return vert_pos;
-                             })
-      .def_property_readonly("tri_verts",
-                             [](Mesh &self) {
-                               const int numTri = self.triVerts.size();
-                               py::array_t<int> tri_verts({numTri, 3});
-                               auto tri_verts_view =
-                                   tri_verts.mutable_unchecked<2>();
+                out.numProp = vertProp.shape(1);
+                vertProp = vertProp.reshape({-1});
+                out.vertProperties = toVector<float>(vertProp);
 
-                               for (int i = 0; i < numTri; ++i)
-                                 for (const int j : {0, 1, 2})
-                                   tri_verts_view(i, j) = self.triVerts[i][j];
-                               return tri_verts;
-                             })
+                if (triVerts.ndim() != 2 || triVerts.shape(1) != 3)
+                  throw std::runtime_error(
+                      "Invalid tri_verts shape, expected (-1, 3)");
+                triVerts = triVerts.reshape({-1});
+                out.triVerts = toVector<uint32_t>(triVerts);
+
+                if (mergeFromVert.has_value())
+                  out.mergeFromVert = toVector<uint32_t>(mergeFromVert.value());
+
+                if (mergeToVert.has_value())
+                  out.mergeToVert = toVector<uint32_t>(mergeToVert.value());
+
+                if (runIndex.has_value())
+                  out.runIndex = toVector<uint32_t>(runIndex.value());
+
+                if (runOriginalID.has_value())
+                  out.runOriginalID = toVector<uint32_t>(runOriginalID.value());
+
+                if (runTransform.has_value()) {
+                  auto runTransform1 = runTransform.value();
+                  if (runTransform1.ndim() != 3 ||
+                      runTransform1.shape(1) != 4 ||
+                      runTransform1.shape(2) != 3)
+                    throw std::runtime_error(
+                        "Invalid run_transform shape, expected (-1, 4, 3)");
+                  runTransform1 = runTransform1.reshape({-1});
+                  out.runTransform = toVector<float>(runTransform1);
+                }
+
+                if (faceID.has_value())
+                  out.faceID = toVector<uint32_t>(faceID.value());
+
+                if (halfedgeTangent.has_value()) {
+                  auto halfedgeTangent1 = halfedgeTangent.value();
+                  if (halfedgeTangent1.ndim() != 3 ||
+                      halfedgeTangent1.shape(1) != 3 ||
+                      halfedgeTangent1.shape(2) != 4)
+                    throw std::runtime_error(
+                        "Invalid halfedge_tangent shape, expected (-1, 3, 4)");
+                  halfedgeTangent1 = halfedgeTangent1.reshape({-1});
+                  out.halfedgeTangent = toVector<float>(halfedgeTangent1);
+                }
+
+                return out;
+              }),
+          py::arg("vert_properties"), py::arg("tri_verts"),
+          py::arg("merge_from_vert") = py::none(),
+          py::arg("merge_to_vert") = py::none(),
+          py::arg("run_index") = py::none(),
+          py::arg("run_original_id") = py::none(),
+          py::arg("run_transform") = py::none(),
+          py::arg("face_id") = py::none(),
+          py::arg("halfedge_tangent") = py::none(), py::arg("precision") = 0)
       .def_property_readonly(
-          "vert_normal",
-          [](Mesh &self) {
-            const int numVert = self.vertNormal.size();
-            py::array_t<float> vert_normal({numVert, 3});
-            auto vert_normal_view = vert_normal.mutable_unchecked<2>();
-
-            for (int i = 0; i < numVert; ++i)
-              for (const int j : {0, 1, 2})
-                vert_normal_view(i, j) = self.vertNormal[i][j];
-            return vert_normal;
+          "vert_properties",
+          [](const MeshGL &self) {
+            return py::array(self.vertProperties.size(),
+                             self.vertProperties.data())
+                .reshape(std::array<int, 2>{-1, (int)self.numProp});
           })
-      .def_property_readonly("halfedge_tangent", [](Mesh &self) {
-        const int numEdge = self.halfedgeTangent.size();
-        py::array_t<float> halfedge_tangent({numEdge, 4});
-        auto halfedge_tangent_view = halfedge_tangent.mutable_unchecked<2>();
-
-        for (int i = 0; i < numEdge; ++i)
-          for (const int j : {0, 1, 2, 3})
-            halfedge_tangent_view(i, j) = self.halfedgeTangent[i][j];
-        return halfedge_tangent;
-      });
+      .def_property_readonly("tri_verts",
+                             [](const MeshGL &self) {
+                               return py::array(self.triVerts.size(),
+                                                self.triVerts.data())
+                                   .reshape(std::array<int, 2>{-1, 3});
+                             })
+      .def_property_readonly("run_transform",
+                             [](const MeshGL &self) {
+                               return py::array(self.runTransform.size(),
+                                                self.runTransform.data())
+                                   .reshape(std::array<int, 3>{-1, 4, 3});
+                             })
+      .def_property_readonly("halfedge_tangent",
+                             [](const MeshGL &self) {
+                               return py::array(self.halfedgeTangent.size(),
+                                                self.halfedgeTangent.data())
+                                   .reshape(std::array<int, 3>{-1, 3, 4});
+                             })
+      .def_readonly("merge_from_vert", &MeshGL::mergeFromVert)
+      .def_readonly("merge_to_vert", &MeshGL::mergeToVert)
+      .def_readonly("run_index", &MeshGL::runIndex)
+      .def_readonly("run_original_id", &MeshGL::runOriginalID)
+      .def_readonly("face_id", &MeshGL::faceID);
 
   py::enum_<CrossSection::FillRule>(m, "FillRule")
       .value("EvenOdd", CrossSection::FillRule::EvenOdd,
@@ -578,48 +648,90 @@ PYBIND11_MODULE(pymanifold, m) {
                        CrossSection::FillRule fillrule) {
              std::vector<SimplePolygon> simplePolygons(polygons.size());
              for (int i = 0; i < polygons.size(); i++) {
-               std::vector<glm::vec2> vertices(polygons[i].size());
+               simplePolygons[i] = std::vector<glm::vec2>(polygons[i].size());
                for (int j = 0; j < polygons[i].size(); j++) {
-                 vertices[j] = {std::get<0>(polygons[i][j]),
-                                std::get<1>(polygons[i][j])};
+                 simplePolygons[i][j] = {std::get<0>(polygons[i][j]),
+                                         std::get<1>(polygons[i][j])};
                }
-               simplePolygons[i] = {vertices};
              }
              return CrossSection(simplePolygons, fillrule);
            }),
            py::arg("polygons"),
            py::arg("fillrule") = CrossSection::FillRule::Positive,
-           "Construct a 2d cross-section from a set of contours.")
-      .def("area", &CrossSection::Area)
-      .def("num_vert", &CrossSection::NumVert)
-      .def("num_contour", &CrossSection::NumContour)
-      .def("is_empty", &CrossSection::IsEmpty)
-      .def("translate",
-           [](CrossSection self, Float2 v) {
-             return self.Translate({std::get<0>(v), std::get<1>(v)});
-           })
-      .def("rotate", &CrossSection::Rotate)
-      .def("scale",
-           [](CrossSection self, Float2 s) {
-             return self.Scale({std::get<0>(s), std::get<1>(s)});
-           })
-      .def("mirror",
-           [](CrossSection self, Float2 ax) {
-             return self.Mirror({std::get<0>(ax), std::get<1>(ax)});
-           })
-      .def("transform",
-           [](CrossSection self, py::array_t<float> &mat) {
-             auto mat_view = mat.unchecked<2>();
-             if (mat_view.shape(0) != 2 || mat_view.shape(1) != 3)
-               throw std::runtime_error("Invalid matrix shape");
-             glm::mat3x2 mat_glm;
-             for (int i = 0; i < 2; i++) {
-               for (int j = 0; j < 3; j++) {
-                 mat_glm[j][i] = mat_view(i, j);
-               }
-             }
-             return self.Transform(mat_glm);
-           })
+           "Create a 2d cross-section from a set of contours (complex "
+           "polygons). A boolean union operation (with Positive filling rule "
+           "by default) performed to combine overlapping polygons and ensure "
+           "the resulting CrossSection is free of intersections."
+           "\n\n"
+           ":param contours: A set of closed paths describing zero or more "
+           "complex polygons.\n"
+           ":param fillrule: The filling rule used to interpret polygon "
+           "sub-regions in contours.")
+      .def("area", &CrossSection::Area,
+           "Return the total area covered by complex polygons making up the "
+           "CrossSection.")
+      .def("num_vert", &CrossSection::NumVert,
+           "Return the number of vertices in the CrossSection.")
+      .def("num_contour", &CrossSection::NumContour,
+           "Return the number of contours (both outer and inner paths) in the "
+           "CrossSection.")
+      .def("is_empty", &CrossSection::IsEmpty,
+           "Does the CrossSection contain any contours?")
+      .def(
+          "translate",
+          [](CrossSection self, Float2 v) {
+            return self.Translate({std::get<0>(v), std::get<1>(v)});
+          },
+          "Move this CrossSection in space. This operation can be chained. "
+          "Transforms are combined and applied lazily."
+          "\n\n"
+          ":param v: The vector to add to every vertex.")
+      .def("rotate", &CrossSection::Rotate,
+           "Applies a (Z-axis) rotation to the CrossSection, in degrees. This "
+           "operation can be chained. Transforms are combined and applied "
+           "lazily."
+           "\n\n"
+           ":param degrees: degrees about the Z-axis to rotate.")
+      .def(
+          "scale",
+          [](CrossSection self, Float2 s) {
+            return self.Scale({std::get<0>(s), std::get<1>(s)});
+          },
+          "Scale this CrossSection in space. This operation can be chained. "
+          "Transforms are combined and applied lazily."
+          "\n\n"
+          ":param v: The vector to multiply every vertex by per component.")
+      .def(
+          "mirror",
+          [](CrossSection self, Float2 ax) {
+            return self.Mirror({std::get<0>(ax), std::get<1>(ax)});
+          },
+          "Mirror this CrossSection over the arbitrary axis described by the "
+          "unit form of the given vector. If the length of the vector is zero, "
+          "an empty CrossSection is returned. This operation can be chained. "
+          "Transforms are combined and applied lazily."
+          "\n\n"
+          ":param ax: the axis to be mirrored over")
+      .def(
+          "transform",
+          [](CrossSection self, py::array_t<float> &mat) {
+            auto mat_view = mat.unchecked<2>();
+            if (mat.ndim() != 2 || mat.shape(0) != 2 || mat.shape(1) != 3)
+              throw std::runtime_error("Invalid matrix shape, expected (2, 3)");
+            glm::mat3x2 mat_glm;
+            for (int i = 0; i < 2; i++) {
+              for (int j = 0; j < 3; j++) {
+                mat_glm[j][i] = mat_view(i, j);
+              }
+            }
+            return self.Transform(mat_glm);
+          },
+          "Transform this CrossSection in space. The first two columns form a "
+          "2x2 matrix transform and the last is a translation vector. This "
+          "operation can be chained. Transforms are combined and applied "
+          "lazily."
+          "\n\n"
+          ":param m: The affine transform matrix to apply to all the vertices.")
       .def(
           "warp",
           [](CrossSection self, const std::function<Float2(Float2)> &f) {
@@ -629,15 +741,55 @@ PYBIND11_MODULE(pymanifold, m) {
               v.y = std::get<1>(fv);
             });
           },
-          py::arg("f"))
-      .def("simplify", &CrossSection::Simplify)
+          py::arg("f"),
+          "Move the vertices of this CrossSection (creating a new one) "
+          "according to any arbitrary input function, followed by a union "
+          "operation (with a Positive fill rule) that ensures any introduced "
+          "intersections are not included in the result."
+          "\n\n"
+          ":param warpFunc: A function that modifies a given vertex position.")
+      .def("simplify", &CrossSection::Simplify,
+           "Remove vertices from the contours in this CrossSection that are "
+           "less than the specified distance epsilon from an imaginary line "
+           "that passes through its two adjacent vertices. Near duplicate "
+           "vertices and collinear points will be removed at lower epsilons, "
+           "with elimination of line segments becoming increasingly aggressive "
+           "with larger epsilons."
+           "\n\n"
+           "It is recommended to apply this function following Offset, in "
+           "order to clean up any spurious tiny line segments introduced that "
+           "do not improve quality in any meaningful way. This is particularly "
+           "important if further offseting operations are to be performed, "
+           "which would compound the issue.")
       .def("offset", &CrossSection::Offset, py::arg("delta"),
            py::arg("join_type"), py::arg("miter_limit") = 2.0,
-           py::arg("arc_tolerance") = 0.0)
+           py::arg("arc_tolerance") = 0.0,
+           "Inflate the contours in CrossSection by the specified delta, "
+           "handling corners according to the given JoinType."
+           "\n\n"
+           ":param delta: Positive deltas will cause the expansion of "
+           "outlining contours to expand, and retraction of inner (hole) "
+           "contours. Negative deltas will have the opposite effect.\n"
+           ":param jt: The join type specifying the treatment of contour joins "
+           "(corners).\n"
+           ":param miter_limit: The maximum distance in multiples of delta "
+           "that vertices can be offset from their original positions with "
+           "before squaring is applied, <B>when the join type is Miter</B> "
+           "(default is 2, which is the minimum allowed). See the [Clipper2 "
+           "MiterLimit](http://www.angusj.com/clipper2/Docs/Units/"
+           "Clipper.Offset/Classes/ClipperOffset/Properties/MiterLimit.htm) "
+           "page for a visual example.\n"
+           ":param circularSegments: Number of segments per 360 degrees of "
+           "<B>JoinType::Round</B> corners (roughly, the number of vertices "
+           "that will be added to each contour). Default is calculated by the "
+           "static Quality defaults according to the radius.")
       .def(py::self + py::self, "Boolean union.")
       .def(py::self - py::self, "Boolean difference.")
       .def(py::self ^ py::self, "Boolean intersection.")
-      .def("decompose", &CrossSection::Decompose)
+      .def("decompose", &CrossSection::Decompose,
+           "This operation returns a vector of CrossSections that are "
+           "topologically disconnected, each containing one outline contour "
+           "with zero or more holes.")
       .def(
           "to_polygons",
           [](CrossSection self) {
@@ -701,11 +853,24 @@ PYBIND11_MODULE(pymanifold, m) {
             return CrossSection::Square({std::get<0>(dims), std::get<1>(dims)},
                                         center);
           },
-          py::arg("dims"), py::arg("center") = false)
+          py::arg("dims"), py::arg("center") = false,
+          "Constructs a square with the given XY dimensions. By default it is "
+          "positioned in the first quadrant, touching the origin. If any "
+          "dimensions in size are negative, or if all are zero, an empty "
+          "Manifold will be returned."
+          "\n\n"
+          ":param size: The X, and Y dimensions of the square.\n"
+          ":param center: Set to true to shift the center to the origin.")
       .def_static(
           "circle",
           [](float radius, int circularSegments) {
             return CrossSection::Circle(radius, circularSegments);
           },
-          py::arg("radius"), py::arg("circularSegments") = 0);
+          py::arg("radius"), py::arg("circularSegments") = 0,
+          "Constructs a circle of a given radius."
+          "\n\n"
+          ":param radius: Radius of the circle. Must be positive.\n"
+          ":param circularSegments: Number of segments along its diameter. "
+          "Default is calculated by the static Quality defaults according to "
+          "the radius.");
 }

bindings/wasm/CMakeLists.txt

@@ -29,7 +29,7 @@ file(MAKE_DIRECTORY ${CMAKE_CURRENT_SOURCE_DIR}/examples/built)
 
 # ensure that interface files are copied over when modified
 add_custom_target(js_deps ALL
-        DEPENDS ${CMAKE_CURRENT_BINARY_DIR}/manifold.*
+        DEPENDS manifoldjs
                ${CMAKE_CURRENT_SOURCE_DIR}/manifold*.d.ts)
 
 add_custom_command(

bindings/wasm/bindings.cpp

@@ -103,7 +103,7 @@ EMSCRIPTEN_BINDINGS(whatever) {
       .function("_Warp", &cross_js::Warp)
       .function("transform", &cross_js::Transform)
       .function("_Translate", &CrossSection::Translate)
-      .function("_Rotate", &CrossSection::Rotate)
+      .function("rotate", &CrossSection::Rotate)
       .function("_Scale", &CrossSection::Scale)
       .function("_Mirror", &CrossSection::Mirror)
       .function("_Decompose", &CrossSection::Decompose)
@@ -114,8 +114,9 @@ EMSCRIPTEN_BINDINGS(whatever) {
       .function("_Bounds", &CrossSection::Bounds)
       .function("simplify", &CrossSection::Simplify)
       .function("_Offset", &cross_js::Offset)
-      .function("_RectClip", &CrossSection::RectClip)
-      .function("_ToPolygons", &CrossSection::ToPolygons);
+      .function("_ToPolygons", &CrossSection::ToPolygons)
+      .function("hull",
+                select_overload<CrossSection() const>(&CrossSection::Hull));
 
   // CrossSection Static Methods
   function("_Square", &CrossSection::Square);

bindings/wasm/bindings.js

@@ -114,10 +114,6 @@ Module.setup = function() {
     return this._Translate(vararg2vec2(vec));
   };
 
-  Module.CrossSection.prototype.rotate = function(vec) {
-    return this._Rotate(...vec);
-  };
-
   Module.CrossSection.prototype.scale = function(vec) {
     // if only one factor provided, scale both x and y with it
     if (typeof vec == 'number') {
@@ -165,14 +161,6 @@ Module.setup = function() {
         delta, joinTypeToInt(joinType), miterLimit, circularSegments);
   };
 
-  Module.CrossSection.prototype.rectClip = function(rect) {
-    const rect2 = {
-      min: {x: rect.min[0], y: rect.min[1]},
-      max: {x: rect.max[0], y: rect.max[1]},
-    };
-    return this._RectClip(rect2);
-  };
-
   Module.CrossSection.prototype.extrude = function(
       height, nDivisions = 0, twistDegrees = 0.0, scaleTop = [1.0, 1.0],
       center = false) {
@@ -182,8 +170,9 @@ Module.setup = function() {
     return (center ? man.translate([0., 0., -height / 2.]) : man);
   };
 
-  Module.CrossSection.prototype.revolve = function(circularSegments = 0) {
-    return Module._Revolve(this, circularSegments);
+  Module.CrossSection.prototype.revolve = function(
+      circularSegments = 0, revolveDegrees = 360.0) {
+    return Module._Revolve(this, circularSegments, revolveDegrees);
   };
 
   Module.CrossSection.prototype.add = function(other) {
@@ -555,11 +544,12 @@ Module.setup = function() {
     return cs.extrude(height, nDivisions, twistDegrees, scaleTop, center);
   };
 
-  Module.Manifold.revolve = function(polygons, circularSegments = 0) {
+  Module.Manifold.revolve = function(
+      polygons, circularSegments = 0, revolveDegrees = 360.0) {
     const cs = (polygons instanceof CrossSectionCtor) ?
         polygons :
         Module.CrossSection(polygons, 'Positive');
-    return cs.revolve(circularSegments);
+    return cs.revolve(circularSegments, revolveDegrees);
   };
 
   Module.Manifold.reserveIDs = function(n) {

bindings/wasm/examples/README.md

@@ -23,6 +23,12 @@ See `package.json` for other useful scripts.
 
 Note that the `emcmake` command automatically copies your WASM build into `built/`, (here, not just under the `buildWASM` directory) which is then packaged by Vite into `dist/assets/`.
 
+To debug the WASM build directly in Chrome dev tools, simply build in debug mode:
+```
+emcmake cmake -DCMAKE_BUILD_TYPE=Debug .. && emmake make
+```
+and install the [DWARF](goo.gle/wasm-debugging-extension) Chrome extension.
+
 When testing [ManifoldCAD.org](https://manifoldcad.org/) (either locally or the
 deployed version) note that it uses a service worker for faster loading. This
 means you need to open the page twice to see updates (the first time loads the

bindings/wasm/examples/make-manifold.html

@@ -74,6 +74,8 @@
   const wasm = await Module();
   wasm.setup();
 
+  const { Manifold, Mesh } = wasm;
+
   const mv = document.querySelector('model-viewer');
   const inputEl = document.querySelector('#input');
   const downloadButton = document.querySelector('#download');
@@ -154,11 +156,11 @@
       const materials = [];
       const tmpMesh = readMesh(mesh, attributes, materials);
       tmpMesh.runOriginalID = [];
-      const firstID = wasm.reserveIDs(materials.length);
+      const firstID = Manifold.reserveIDs(materials.length);
       for (let i = 0; i < materials.length; ++i) {
         tmpMesh.runOriginalID.push(firstID + i);
       }
-      const manifoldMesh = new wasm.Mesh(tmpMesh);
+      const manifoldMesh = new Mesh(tmpMesh);
       mesh.dispose();
 
       const matAttributes = [];
@@ -168,7 +170,7 @@
 
       try {
         // Test manifoldness - will throw if not.
-        const manifold = wasm.Manifold(manifoldMesh);
+        const manifold = Manifold(manifoldMesh);
         node.setMesh(writeMesh(docIn, manifold.getMesh(), matAttributes, materials));
         manifold.delete();
         anyManifold = true;

bindings/wasm/examples/model-viewer.html

@@ -61,6 +61,8 @@
   const wasm = await Module();
   wasm.setup();
 
+  const { Manifold, Mesh } = wasm;
+
   const attributes = ['POSITION', 'TEXCOORD_0'];
   const materialMap = new Map();
 
@@ -75,12 +77,12 @@
       const materials = [];
       const manifoldMesh = readMesh(node.getMesh(), attributes, materials);
       manifoldMesh.runOriginalID = [];
-      const firstID = wasm.reserveIDs(materials.length);
+      const firstID = Manifold.reserveIDs(materials.length);
       for (let i = 0; i < materials.length; ++i) {
         manifoldMesh.runOriginalID.push(firstID + i);
         ids.push(firstID + i);
       }
-      manifolds.push(wasm.Manifold(new wasm.Mesh(manifoldMesh)));
+      manifolds.push(Manifold(new Mesh(manifoldMesh)));
     }
     // pull in materials, TODO: replace with transfer() when available
     const startIdx = doc.getRoot().listMaterials().length;
@@ -92,7 +94,7 @@
       materialMap.set(id, material);
     }
 
-    return wasm.union(manifolds);
+    return Manifold.union(manifolds);
   }
 
   const space = await readGLB('/models/space.glb', attributes, materialMap);
@@ -102,7 +104,7 @@
   doc.createScene().addChild(node);
 
   function csg(operation) {
-    push2MV(wasm[operation](space, moon));
+    push2MV(Manifold[operation](space, moon));
   };
 
   document.querySelector('select').onchange = function (event) {

bindings/wasm/examples/three.html

@@ -52,6 +52,8 @@
   const wasm = await Module();
   wasm.setup();
 
+  const { Manifold, Mesh } = wasm;
+
   // we have manifold module, let's do some three.js
   const camera = new THREE.PerspectiveCamera(30, 1, 0.01, 10);
   camera.position.z = 1;
@@ -64,12 +66,12 @@
 
   const icosahedron = simplify(new THREE.IcosahedronGeometry(0.16));
 
-  const manifold_1 = wasm.cube([0.2, 0.2, 0.2], true);
-  const manifold_2 = new wasm.Manifold(geometry2mesh(icosahedron));
+  const manifold_1 = Manifold.cube([0.2, 0.2, 0.2], true);
+  const manifold_2 = new Manifold(geometry2mesh(icosahedron));
 
   const csg = function (operation) {
     mesh.geometry?.dispose();
-    mesh.geometry = mesh2geometry(wasm[operation](manifold_1, manifold_2).getMesh());
+    mesh.geometry = mesh2geometry(Manifold[operation](manifold_1, manifold_2).getMesh());
   };
 
   document.querySelector('select').onchange = function (event) {
@@ -92,7 +94,7 @@
   function geometry2mesh(geometry) {
     const vertProperties = geometry.attributes.position.array;
     const triVerts = geometry.index.array;
-    return new wasm.Mesh({ numProp: 3, vertProperties, triVerts });
+    return new Mesh({ numProp: 3, vertProperties, triVerts });
   }
 
   function mesh2geometry(mesh) {

bindings/wasm/manifold-encapsulated-types.d.ts

@@ -91,7 +91,7 @@ export class CrossSection {
    * @param circularSegments Number of segments along its diameter. Default is
    * calculated by the static Defaults.
    */
-  revolve(circularSegments?: number): Manifold;
+  revolve(circularSegments?: number, revolveDegrees?: number): Manifold;
 
   // Transformations
 
@@ -113,16 +113,12 @@ export class CrossSection {
   translate(v: Vec2): CrossSection;
 
   /**
-   * Applies an Euler angle rotation to the cross-section, first about the X
-   * axis, then Y, then Z, in degrees. We use degrees so that we can minimize
-   * rounding error, and eliminate it completely for any multiples of 90
-   * degrees. Additionally, more efficient code paths are used to update the
-   * cross-section when the transforms only rotate by multiples of 90 degrees.
-   * This operation can be chained. Transforms are combined and applied lazily.
+   * Applies a (Z-axis) rotation to the CrossSection, in degrees. This operation
+   * can be chained. Transforms are combined and applied lazily.
    *
-   * @param v [X, Y, Z] rotation in degrees.
+   * @param degrees degrees about the Z-axis to rotate.
    */
-  rotate(v: Vec2): CrossSection;
+  rotate(v: number): CrossSection;
 
   /**
    * Scale this CrossSection in space. This operation can be chained. Transforms
@@ -177,6 +173,11 @@ export class CrossSection {
       delta: number, joinType?: JoinType, miterLimit?: number,
       circularSegments?: number): CrossSection;
 
+  /**
+   * Compute the convex hull of the contours in this CrossSection.
+   */
+  hull(): CrossSection;
+
   /**
    * Remove vertices from the contours in this CrossSection that are less than
    * the specified distance epsilon from an imaginary line that passes through
@@ -245,14 +246,6 @@ export class CrossSection {
    */
   static intersection(polygons: (CrossSection|Polygons)[]): CrossSection;
 
-  /**
-   * Compute the intersection between a cross-section and an axis-aligned
-   * rectangle. This operation has much higher performance (O(n) vs
-   * >O(n^3)) than the general purpose intersection algorithm
-   * used for sets of cross-sections.
-   */
-  rectClip(rect: Rect): CrossSection;
-
   // Topological Operations
 
   /**
@@ -421,9 +414,11 @@ export class Manifold {
    * @param polygons A set of non-overlapping polygons to revolve.
    * @param circularSegments Number of segments along its diameter. Default is
    * calculated by the static Defaults.
+   * @param revolveDegrees Number of degrees to revolve. Default is 360 degrees.
    */
-  static revolve(polygons: CrossSection|Polygons, circularSegments?: number):
-      Manifold;
+  static revolve(
+      polygons: CrossSection|Polygons, circularSegments?: number,
+      revolveDegrees?: number): Manifold;
 
   // Mesh Conversion
 

flake.lock

@@ -17,16 +17,16 @@
     },
     "nixpkgs": {
       "locked": {
-        "lastModified": 1657296039,
-        "narHash": "sha256-Ghh39+aS+pw5sTP/ZO8VIKE6sBhMadDaQZtf+3yu4Vc=",
+        "lastModified": 1689008574,
+        "narHash": "sha256-VFMgyHDiqsGDkRg73alv6OdHJAqhybryWHv77bSCGIw=",
         "owner": "NixOS",
         "repo": "nixpkgs",
-        "rev": "71d7a4c037dc4f3e98d5c4a81b941933cf5bf675",
+        "rev": "4a729ce4b1fe5ec4fffc71c67c96aa5184ebb462",
         "type": "github"
       },
       "original": {
         "id": "nixpkgs",
-        "ref": "nixos-22.05",
+        "ref": "nixos-unstable",
         "type": "indirect"
       }
     },

flake.nix

@@ -1,6 +1,6 @@
 {
   inputs.flake-utils.url = "github:numtide/flake-utils";
-  inputs.nixpkgs.url = "nixpkgs/nixos-22.05";
+  inputs.nixpkgs.url = "nixpkgs/nixos-unstable";
 
   outputs = { self, nixpkgs, flake-utils }:
     flake-utils.lib.eachDefaultSystem
@@ -25,7 +25,7 @@
                   "manifold-${parallel-backend}";
               version = "beta";
               src = self;
-              nativeBuildInputs = (with pkgs; [ cmake (python39.withPackages(ps: with ps; [trimesh])) ]) ++ build-tools ++
+              nativeBuildInputs = (with pkgs; [ cmake (python39.withPackages (ps: with ps; [ trimesh ])) ]) ++ build-tools ++
                 (if cuda-support then with pkgs.cudaPackages; [ cuda_nvcc cuda_cudart cuda_cccl pkgs.addOpenGLRunpath ] else [ ]);
               cmakeFlags = [
                 "-DMANIFOLD_PYBIND=ON"
@@ -43,7 +43,7 @@
                 mkdir -p $out
                 cp src/manifold/libmanifold.a $out/
                 cp extras/perfTest $out
-                cp bindings/python/pymanifold* $out
+                cp bindings/python/manifold3d* $out
               '';
             };
           parallelBackends = [

meshIO/src/meshIO.cpp

@@ -79,6 +79,13 @@ void ExportScene(aiScene* scene, const std::string& filename) {
   //   std::cout << i << ", id = " << desc->id << ", " << desc->description
   //             << std::endl;
   // }
+  auto ext = filename.substr(filename.length() - 4, 4);
+  if (ext == ".3mf") {
+    // Workaround https://github.com/assimp/assimp/issues/3816
+    aiNode* old_root = scene->mRootNode;
+    scene->mRootNode = new aiNode();
+    scene->mRootNode->addChildren(1, &old_root);
+  }
 
   auto result = exporter.Export(scene, GetType(filename), filename);
 

src/cross_section/include/cross_section.h

@@ -77,6 +77,7 @@ class CrossSection {
                FillRule fillrule = FillRule::Positive);
   CrossSection(const Polygons& contours,
                FillRule fillrule = FillRule::Positive);
+  CrossSection(const Rect& rect);
   static CrossSection Square(const glm::vec2 dims, bool center = false);
   static CrossSection Circle(float radius, int circularSegments = 0);
   ///@}
@@ -143,7 +144,6 @@ class CrossSection {
   CrossSection& operator-=(const CrossSection&);
   CrossSection operator^(const CrossSection&) const;
   CrossSection& operator^=(const CrossSection&);
-  CrossSection RectClip(const Rect& rect) const;
   ///@}
 
   /** @name Topological
@@ -153,6 +153,15 @@ class CrossSection {
   std::vector<CrossSection> Decompose() const;
   ///@}
 
+  /** @name Convex Hulling
+   */
+  ///@{
+  CrossSection Hull() const;
+  static CrossSection Hull(const std::vector<CrossSection>& crossSections);
+  static CrossSection Hull(const SimplePolygon poly);
+  static CrossSection Hull(const Polygons polys);
+  ///@}
+  ///
   /** @name Conversion
    */
   ///@{
@@ -197,6 +206,7 @@ class Rect {
    */
   ///@{
   glm::vec2 Size() const;
+  float Area() const;
   float Scale() const;
   glm::vec2 Center() const;
   bool Contains(const glm::vec2& pt) const;

src/cross_section/src/cross_section.cpp

@@ -132,6 +132,58 @@ void decompose_hole(const C2::PolyTreeD* outline,
   }
 }
 
+bool V2Lesser(glm::vec2 a, glm::vec2 b) {
+  if (a.x == b.x) return a.y < b.y;
+  return a.x < b.x;
+}
+
+double Cross(glm::vec2 a, glm::vec2 b) { return a.x * b.y - a.y * b.x; }
+
+double SquaredDistance(glm::vec2 a, glm::vec2 b) {
+  auto d = a - b;
+  return glm::dot(d, d);
+}
+
+bool IsCw(glm::vec2 a, glm::vec2 b, glm::vec2 c) {
+  double lhs = Cross(a - c, b - c);
+  double rhs = 1e-18 * SquaredDistance(a, c) * SquaredDistance(b, c);
+  return lhs * std::abs(lhs) <= rhs;
+}
+
+void HullBacktrack(const SimplePolygon& pts, const int idx,
+                   std::vector<int>& keep, const int hold) {
+  const int stop = keep.size() - hold;
+  int i = 0;
+  while (i < stop && !IsCw(pts[idx], pts[keep[keep.size() - 1]],
+                           pts[keep[keep.size() - 2]])) {
+    keep.pop_back();
+    i++;
+  }
+}
+
+// Based on method described here:
+// https://www.hackerearth.com/practice/math/geometry/line-sweep-technique/tutorial/
+C2::PathD HullImpl(SimplePolygon& pts) {
+  int len = pts.size();
+  if (len < 3) return C2::PathD();  // not enough points to create a polygon
+  std::sort(pts.begin(), pts.end(), V2Lesser);
+  auto keep = std::vector<int>{0, 1};
+  for (int i = 2; i < len; i++) {
+    HullBacktrack(pts, i, keep, 1);
+    keep.push_back(i);
+  }
+  int nLower = keep.size();
+  for (int i = 0; i < len - 1; i++) {
+    int idx = len - 2 - i;
+    HullBacktrack(pts, idx, keep, nLower);
+    if (idx > 0) keep.push_back(idx);
+  }
+  auto path = C2::PathD(keep.size());
+  for (int i = 0; i < keep.size(); i++) {
+    path[i] = v2_to_pd(pts[keep[i]]);
+  }
+  return path;
+}
 }  // namespace
 
 namespace manifold {
@@ -202,6 +254,20 @@ CrossSection::CrossSection(const Polygons& contours, FillRule fillrule) {
   paths_ = shared_paths(C2::Union(ps, fr(fillrule), precision_));
 }
 
+/**
+ * Create a 2d cross-section from an axis-aligned rectangle (bounding box).
+ *
+ * @param rect An axis-aligned rectangular bounding box.
+ */
+CrossSection::CrossSection(const Rect& rect) {
+  C2::PathD p(4);
+  p[0] = C2::PointD(rect.min.x, rect.min.y);
+  p[1] = C2::PointD(rect.max.x, rect.min.y);
+  p[2] = C2::PointD(rect.max.x, rect.max.y);
+  p[3] = C2::PointD(rect.min.x, rect.max.y);
+  paths_ = shared_paths(C2::PathsD{p});
+}
+
 // Private
 // All access to paths_ should be done through the GetPaths() method, which
 // applies the accumulated transform_
@@ -391,18 +457,6 @@ std::vector<CrossSection> CrossSection::Decompose() const {
   return comps;
 }
 
-/**
- * Compute the intersection between a cross-section and an axis-aligned
- * rectangle. This operation has much higher performance (<B>O(n)</B> vs
- * <B>O(n<SUP>3</SUP>)</B>) than the general purpose intersection algorithm
- * used for sets of cross-sections.
- */
-CrossSection CrossSection::RectClip(const Rect& rect) const {
-  auto r = C2::RectD(rect.min.x, rect.min.y, rect.max.x, rect.max.y);
-  auto ps = C2::ExecuteRectClip(r, GetPaths(), false, precision_);
-  return CrossSection(ps);
-}
-
 /**
  * Move this CrossSection in space. This operation can be chained. Transforms
  * are combined and applied lazily.
@@ -558,6 +612,41 @@ CrossSection CrossSection::Offset(double delta, JoinType jointype,
   return CrossSection(ps);
 }
 
+CrossSection CrossSection::Hull(
+    const std::vector<CrossSection>& crossSections) {
+  int n = 0;
+  for (auto cs : crossSections) n += cs.NumVert();
+  SimplePolygon pts;
+  pts.reserve(n);
+  for (auto cs : crossSections) {
+    auto paths = cs.GetPaths();
+    for (auto path : paths) {
+      for (auto p : path) {
+        pts.push_back(v2_of_pd(p));
+      }
+    }
+  }
+  return CrossSection(C2::PathsD{HullImpl(pts)});
+}
+
+CrossSection CrossSection::Hull() const {
+  return Hull(std::vector<CrossSection>{*this});
+}
+
+CrossSection CrossSection::Hull(SimplePolygon pts) {
+  return CrossSection(C2::PathsD{HullImpl(pts)});
+}
+
+CrossSection CrossSection::Hull(const Polygons polys) {
+  SimplePolygon pts;
+  for (auto poly : polys) {
+    for (auto p : poly) {
+      pts.push_back(p);
+    }
+  }
+  return Hull(pts);
+}
+
 /**
  * Return the total area covered by complex polygons making up the
  * CrossSection.
@@ -642,6 +731,14 @@ Rect::Rect(const glm::vec2 a, const glm::vec2 b) {
  */
 glm::vec2 Rect::Size() const { return max - min; }
 
+/**
+ * Return the area of the rectangle.
+ */
+float Rect::Area() const {
+  auto sz = Size();
+  return sz.x * sz.y;
+}
+
 /**
  * Returns the absolute-largest coordinate value of any contained
  * point.
@@ -657,7 +754,7 @@ float Rect::Scale() const {
 glm::vec2 Rect::Center() const { return 0.5f * (max + min); }
 
 /**
- * Does this rectangle contain (includes equal) the given point?
+ * Does this rectangle contain (includes on border) the given point?
  */
 bool Rect::Contains(const glm::vec2& p) const {
   return glm::all(glm::greaterThanEqual(p, min)) &&
@@ -766,13 +863,6 @@ Rect Rect::Transform(const glm::mat3x2& m) const {
  * Return a CrossSection with an outline defined by this axis-aligned
  * rectangle.
  */
-CrossSection Rect::AsCrossSection() const {
-  SimplePolygon p(4);
-  p[0] = glm::vec2(min.x, max.y);
-  p[1] = glm::vec2(max.x, max.y);
-  p[2] = glm::vec2(max.x, min.y);
-  p[3] = glm::vec2(min.x, min.y);
-  return CrossSection(p);
-}
+CrossSection Rect::AsCrossSection() const { return CrossSection(*this); }
 
 }  // namespace manifold

src/manifold/CMakeLists.txt

@@ -25,7 +25,7 @@ endif()
 target_include_directories(${PROJECT_NAME} PUBLIC ${PROJECT_SOURCE_DIR}/include)
 target_link_libraries(${PROJECT_NAME}
     PUBLIC utilities cross_section
-    PRIVATE collider polygon ${MANIFOLD_INCLUDE} graphlite Clipper2
+    PRIVATE collider polygon ${MANIFOLD_INCLUDE} graphlite Clipper2 quickhull
 )
 
 target_compile_features(${PROJECT_NAME}

src/manifold/include/manifold.h

@@ -141,7 +141,8 @@ class Manifold {
                           int nDivisions = 0, float twistDegrees = 0.0f,
                           glm::vec2 scaleTop = glm::vec2(1.0f));
   static Manifold Revolve(const CrossSection& crossSection,
-                          int circularSegments = 0);
+                          int circularSegments = 0,
+                          float revolveDegrees = 360.0f);
   ///@}
 
   /** @name Topological
@@ -233,6 +234,10 @@ class Manifold {
   Manifold TrimByPlane(glm::vec3 normal, float originOffset) const;
   ///@}
 
+  Manifold Hull() const;
+  static Manifold Hull(const std::vector<Manifold>& manifolds);
+  static Manifold Hull(const std::vector<glm::vec3>& pts);
+
   /** @name Testing hooks
    *  These are just for internal testing.
    */
@@ -248,7 +253,6 @@ class Manifold {
   Manifold(std::shared_ptr<CsgNode> pNode_);
   Manifold(std::shared_ptr<Impl> pImpl_);
   static Manifold Invalid();
-
   mutable std::shared_ptr<CsgNode> pNode_;
 
   CsgLeafNode& GetCsgLeafNode() const;

src/manifold/src/boolean_result.cpp

@@ -481,7 +481,10 @@ void CreateProperties(Manifold::Impl &outR, const VecDH<TriRef> &refPQ,
                           inP.halfedge_.cptrD(), inQ.halfedge_.cptrD(),
                           outR.halfedge_.cptrD(), outR.precision_}));
 
-  std::map<std::tuple<int, int, int, int>, int> propIdx;
+  using Entry = std::pair<glm::ivec3, int>;
+  int idMissProp = outR.NumVert();
+  std::vector<std::vector<Entry>> propIdx(outR.NumVert() + 1);
+  outR.meshRelation_.properties.reserve(outR.NumVert() * numProp);
   int idx = 0;
 
   for (int tri = 0; tri < numTri; ++tri) {
@@ -490,47 +493,61 @@ void CreateProperties(Manifold::Impl &outR, const VecDH<TriRef> &refPQ,
 
     const int triPQ = refPQ[tri].tri;
     const bool PQ = refPQ[tri].meshID == 0;
+    const int oldNumProp = PQ ? numPropP : numPropQ;
     const auto &properties =
         PQ ? inP.meshRelation_.properties : inQ.meshRelation_.properties;
-    const auto &triProp = PQ ? inP.meshRelation_.triProperties[triPQ]
-                             : inQ.meshRelation_.triProperties[triPQ];
+    const glm::ivec3 &triProp = oldNumProp == 0 ? glm::ivec3(-1)
+                                : PQ ? inP.meshRelation_.triProperties[triPQ]
+                                     : inQ.meshRelation_.triProperties[triPQ];
 
     for (const int i : {0, 1, 2}) {
       const int vert = outR.halfedge_[3 * tri + i].startVert;
-      const glm::vec3 uvw = bary[3 * tri + i];
+      const glm::vec3 &uvw = bary[3 * tri + i];
 
-      auto key = std::make_tuple(PQ, vert, -1, -1);
-      int edge = -1;
-      for (const int j : {0, 1, 2}) {
-        if (uvw[j] == 1) {
-          // On a retained vert, the propVert must also match
-          std::get<2>(key) = triProp[j];
-          edge = -1;
+      glm::ivec4 key(PQ, idMissProp, -1, -1);
+      if (oldNumProp > 0) {
+        key[1] = vert;
+        int edge = -1;
+        for (const int j : {0, 1, 2}) {
+          if (uvw[j] == 1) {
+            // On a retained vert, the propVert must also match
+            key[2] = triProp[j];
+            edge = -1;
+            break;
+          }
+          if (uvw[j] == 0) edge = j;
+        }
+        if (edge >= 0) {
+          // On an edge, both propVerts must match
+          const int p0 = triProp[Next3(edge)];
+          const int p1 = triProp[Prev3(edge)];
+          key[2] = glm::min(p0, p1);
+          key[3] = glm::max(p0, p1);
+        }
+      }
+
+      auto &bin = propIdx[key.y];
+      bool bFound = false;
+      for (int k = 0; k < bin.size(); ++k) {
+        if (bin[k].first == glm::ivec3(key.x, key.z, key.w)) {
+          bFound = true;
+          outR.meshRelation_.triProperties[tri][i] = bin[k].second;
           break;
         }
-        if (uvw[j] == 0) edge = j;
       }
-      if (edge >= 0) {
-        // On an edge, both propVerts must match
-        const int p0 = triProp[Next3(edge)];
-        const int p1 = triProp[Prev3(edge)];
-        std::get<2>(key) = glm::min(p0, p1);
-        std::get<3>(key) = glm::max(p0, p1);
-      }
-
-      const auto it = propIdx.find(key);
-      if (it != propIdx.end()) {
-        outR.meshRelation_.triProperties[tri][i] = it->second;
-        continue;
-      }
-      outR.meshRelation_.triProperties[tri][i] = idx;
-      propIdx.insert({key, idx++});
+      if (bFound) continue;
+      bin.push_back(std::make_pair(glm::ivec3(key.x, key.z, key.w), idx));
+      outR.meshRelation_.triProperties[tri][i] = idx++;
 
       for (int p = 0; p < numProp; ++p) {
-        glm::vec3 oldProps;
-        for (const int j : {0, 1, 2})
-          oldProps[j] = properties[numProp * triProp[j] + p];
-        outR.meshRelation_.properties.push_back(glm::dot(uvw, oldProps));
+        if (p < oldNumProp) {
+          glm::vec3 oldProps;
+          for (const int j : {0, 1, 2})
+            oldProps[j] = properties[oldNumProp * triProp[j] + p];
+          outR.meshRelation_.properties.push_back(glm::dot(uvw, oldProps));
+        } else {
+          outR.meshRelation_.properties.push_back(0);
+        }
       }
     }
   }
@@ -551,12 +568,12 @@ Manifold::Impl Boolean3::Result(OpType op) const {
   const int c2 = op == OpType::Add ? 1 : 0;
   const int c3 = op == OpType::Intersect ? 1 : -1;
 
-  if (w03_.size() == 0) {
-    if (w30_.size() != 0 && op == OpType::Add) {
+  if (inP_.IsEmpty()) {
+    if (!inQ_.IsEmpty() && op == OpType::Add) {
       return inQ_;
     }
     return Manifold::Impl();
-  } else if (w30_.size() == 0) {
+  } else if (inQ_.IsEmpty()) {
     if (op == OpType::Intersect) {
       return Manifold::Impl();
     }

src/manifold/src/constructors.cpp

@@ -317,92 +317,122 @@ Manifold Manifold::Extrude(const CrossSection& crossSection, float height,
  * @param crossSection A set of non-overlapping polygons to revolve.
  * @param circularSegments Number of segments along its diameter. Default is
  * calculated by the static Defaults.
+ * @param revolveDegrees Number of degrees to revolve. Default is 360 degrees.
  */
 Manifold Manifold::Revolve(const CrossSection& crossSection,
-                           int circularSegments) {
-  auto polygons = crossSection.ToPolygons();
+                           int circularSegments, float revolveDegrees) {
+  Polygons polygons = crossSection.ToPolygons();
+
   if (polygons.size() == 0) {
     return Invalid();
   }
 
+  const Rect bounds = crossSection.Bounds();
+
+  if (bounds.max.x <= 0) {
+    return Invalid();
+  } else if (bounds.min.x < 0) {
+    // Take the x>=0 slice.
+    glm::vec2 min = bounds.min;
+    glm::vec2 max = bounds.max;
+    CrossSection posBoundingBox = CrossSection(
+        {{0.0, min.y}, {max.x, min.y}, {max.x, max.y}, {0.0, max.y}});
+
+    polygons = (crossSection ^ posBoundingBox).ToPolygons();
+  }
+
   float radius = 0.0f;
   for (const auto& poly : polygons) {
     for (const auto& vert : poly) {
       radius = fmax(radius, vert.x);
     }
   }
-  int nDivisions = circularSegments > 2 ? circularSegments
-                                        : Quality::GetCircularSegments(radius);
+
+  if (revolveDegrees > 360.0f) {
+    revolveDegrees = 360.0f;
+  }
+  const bool isFullRevolution = revolveDegrees == 360.0f;
+
+  const int nDivisions = circularSegments > 2
+                             ? circularSegments
+                             : Quality::GetCircularSegments(radius);
+
   auto pImpl_ = std::make_shared<Impl>();
   auto& vertPos = pImpl_->vertPos_;
   VecDH<glm::ivec3> triVertsDH;
   auto& triVerts = triVertsDH;
-  float dPhi = 360.0f / nDivisions;
+
+  std::vector<int> startPoses;
+  std::vector<int> endPoses;
+
+  const float dPhi = revolveDegrees / nDivisions;
+  // first and last slice are distinguished if not a full revolution.
+  const int nSlices = isFullRevolution ? nDivisions : nDivisions + 1;
+
   for (const auto& poly : polygons) {
-    int start = -1;
-    for (int polyVert = 0; polyVert < poly.size(); ++polyVert) {
-      if (poly[polyVert].x <= 0) {
-        start = polyVert;
-        break;
+    std::size_t nPosVerts = 0;
+    std::size_t nRevolveAxisVerts = 0;
+    for (auto& pt : poly) {
+      if (pt.x > 0) {
+        nPosVerts++;
+      } else {
+        nRevolveAxisVerts++;
       }
     }
-    if (start == -1) {  // poly all positive
-      for (int polyVert = 0; polyVert < poly.size(); ++polyVert) {
-        int startVert = vertPos.size();
-        int lastStart =
-            startVert +
-            (polyVert == 0 ? nDivisions * (poly.size() - 1) : -nDivisions);
-        for (int slice = 0; slice < nDivisions; ++slice) {
-          int lastSlice = (slice == 0 ? nDivisions : slice) - 1;
-          float phi = slice * dPhi;
-          glm::vec2 pos = poly[polyVert];
-          vertPos.push_back({pos.x * cosd(phi), pos.x * sind(phi), pos.y});
-          triVerts.push_back({startVert + slice, startVert + lastSlice,
-                              lastStart + lastSlice});
-          triVerts.push_back(
-              {lastStart + lastSlice, lastStart + slice, startVert + slice});
+
+    for (int polyVert = 0; polyVert < poly.size(); ++polyVert) {
+      const int startPosIndex = vertPos.size();
+
+      if (!isFullRevolution) startPoses.push_back(startPosIndex);
+
+      const int prevStartPosIndex =
+          startPosIndex +
+          (polyVert == 0 ? nRevolveAxisVerts + (nSlices * nPosVerts) - nSlices
+           : poly[polyVert - 1].x == 0.0 ? -1
+                                         : -nSlices);
+      glm::vec2 currPolyVertex = poly[polyVert];
+      glm::vec2 prevPolyVertex =
+          poly[polyVert == 0 ? poly.size() - 1 : polyVert - 1];
+
+      for (int slice = 0; slice < nSlices; ++slice) {
+        const float phi = slice * dPhi;
+        if (slice == 0 || currPolyVertex.x > 0) {
+          vertPos.push_back({currPolyVertex.x * cosd(phi),
+                             currPolyVertex.x * sind(phi), currPolyVertex.y});
+        }
+
+        if (isFullRevolution || slice > 0) {
+          const int lastSlice = (slice == 0 ? nDivisions : slice) - 1;
+          if (currPolyVertex.x > 0.0) {
+            triVerts.push_back(
+                {startPosIndex + slice, startPosIndex + lastSlice,
+                 // "Reuse" vertex of first slice if it lies on the revolve axis
+                 (prevPolyVertex.x == 0.0 ? prevStartPosIndex
+                                          : prevStartPosIndex + lastSlice)});
+          }
+
+          if (prevPolyVertex.x > 0.0) {
+            triVerts.push_back(
+                {prevStartPosIndex + lastSlice, prevStartPosIndex + slice,
+                 (currPolyVertex.x == 0.0 ? startPosIndex
+                                          : startPosIndex + slice)});
+          }
         }
       }
-    } else {  // poly crosses zero
-      int polyVert = start;
-      glm::vec2 pos = poly[polyVert];
-      do {
-        glm::vec2 lastPos = pos;
-        polyVert = (polyVert + 1) % poly.size();
-        pos = poly[polyVert];
-        if (pos.x > 0) {
-          if (lastPos.x <= 0) {
-            float a = pos.x / (pos.x - lastPos.x);
-            vertPos.push_back({0.0f, 0.0f, glm::mix(pos.y, lastPos.y, a)});
-          }
-          int startVert = vertPos.size();
-          for (int slice = 0; slice < nDivisions; ++slice) {
-            int lastSlice = (slice == 0 ? nDivisions : slice) - 1;
-            float phi = slice * dPhi;
-            glm::vec2 pos = poly[polyVert];
-            vertPos.push_back({pos.x * cosd(phi), pos.x * sind(phi), pos.y});
-            if (lastPos.x > 0) {
-              triVerts.push_back({startVert + slice, startVert + lastSlice,
-                                  startVert - nDivisions + lastSlice});
-              triVerts.push_back({startVert - nDivisions + lastSlice,
-                                  startVert - nDivisions + slice,
-                                  startVert + slice});
-            } else {
-              triVerts.push_back(
-                  {startVert - 1, startVert + slice, startVert + lastSlice});
-            }
-          }
-        } else if (lastPos.x > 0) {
-          int startVert = vertPos.size();
-          float a = pos.x / (pos.x - lastPos.x);
-          vertPos.push_back({0.0f, 0.0f, glm::mix(pos.y, lastPos.y, a)});
-          for (int slice = 0; slice < nDivisions; ++slice) {
-            int lastSlice = (slice == 0 ? nDivisions : slice) - 1;
-            triVerts.push_back({startVert, startVert - nDivisions + lastSlice,
-                                startVert - nDivisions + slice});
-          }
-        }
-      } while (polyVert != start);
+      if (!isFullRevolution) endPoses.push_back(vertPos.size() - 1);
+    }
+  }
+
+  // Add front and back triangles if not a full revolution.
+  if (!isFullRevolution) {
+    std::vector<glm::ivec3> frontTriangles =
+        Triangulate(polygons, pImpl_->precision_);
+    for (auto& t : frontTriangles) {
+      triVerts.push_back({startPoses[t.x], startPoses[t.y], startPoses[t.z]});
+    }
+
+    for (auto& t : frontTriangles) {
+      triVerts.push_back({endPoses[t.z], endPoses[t.y], endPoses[t.x]});
     }
   }
 

src/manifold/src/csg_tree.cpp

@@ -45,6 +45,15 @@ struct UpdateHalfedge {
   }
 };
 
+struct UpdateTriProp {
+  const int nextProp;
+
+  __host__ __device__ glm::ivec3 operator()(glm::ivec3 tri) {
+    tri += nextProp;
+    return tri;
+  }
+};
+
 struct UpdateMeshIDs {
   const int offset;
 
@@ -142,9 +151,12 @@ Manifold::Impl CsgLeafNode::Compose(
   int numVert = 0;
   int numEdge = 0;
   int numTri = 0;
+  int numPropVert = 0;
   std::vector<int> vertIndices;
   std::vector<int> edgeIndices;
   std::vector<int> triIndices;
+  std::vector<int> propVertIndices;
+  int numPropOut = 0;
   for (auto &node : nodes) {
     float nodeOldScale = node->pImpl_->bBox_.Scale();
     float nodeNewScale =
@@ -158,9 +170,15 @@ Manifold::Impl CsgLeafNode::Compose(
     vertIndices.push_back(numVert);
     edgeIndices.push_back(numEdge * 2);
     triIndices.push_back(numTri);
+    propVertIndices.push_back(numPropVert);
     numVert += node->pImpl_->NumVert();
     numEdge += node->pImpl_->NumEdge();
     numTri += node->pImpl_->NumTri();
+    const int numProp = node->pImpl_->NumProp();
+    numPropOut = glm::max(numPropOut, numProp);
+    numPropVert +=
+        numProp == 0 ? 1
+                     : node->pImpl_->meshRelation_.properties.size() / numProp;
   }
 
   Manifold::Impl combined;
@@ -170,6 +188,11 @@ Manifold::Impl CsgLeafNode::Compose(
   combined.faceNormal_.resize(numTri);
   combined.halfedgeTangent_.resize(2 * numEdge);
   combined.meshRelation_.triRef.resize(numTri);
+  if (numPropOut > 0) {
+    combined.meshRelation_.numProp = numPropOut;
+    combined.meshRelation_.properties.resize(numPropOut * numPropVert, 0);
+    combined.meshRelation_.triProperties.resize(numTri);
+  }
   auto policy = autoPolicy(numTri);
 
   // if we are already parallelizing for each node, do not perform multithreaded
@@ -182,12 +205,44 @@ Manifold::Impl CsgLeafNode::Compose(
   for_each_n_host(
       nodes.size() > 1 ? ExecutionPolicy::Par : ExecutionPolicy::Seq,
       countAt(0), nodes.size(),
-      [&nodes, &vertIndices, &edgeIndices, &triIndices, &combined,
-       policy](int i) {
+      [&nodes, &vertIndices, &edgeIndices, &triIndices, &propVertIndices,
+       numPropOut, &combined, policy](int i) {
         auto &node = nodes[i];
         copy(policy, node->pImpl_->halfedgeTangent_.begin(),
              node->pImpl_->halfedgeTangent_.end(),
              combined.halfedgeTangent_.begin() + edgeIndices[i]);
+        transform(
+            policy, node->pImpl_->halfedge_.begin(),
+            node->pImpl_->halfedge_.end(),
+            combined.halfedge_.begin() + edgeIndices[i],
+            UpdateHalfedge({vertIndices[i], edgeIndices[i], triIndices[i]}));
+
+        if (numPropOut > 0) {
+          auto start =
+              combined.meshRelation_.triProperties.begin() + triIndices[i];
+          if (node->pImpl_->NumProp() > 0) {
+            auto &triProp = node->pImpl_->meshRelation_.triProperties;
+            transform(policy, triProp.begin(), triProp.end(), start,
+                      UpdateTriProp({propVertIndices[i]}));
+
+            const int numProp = node->pImpl_->NumProp();
+            auto &oldProp = node->pImpl_->meshRelation_.properties;
+            auto &newProp = combined.meshRelation_.properties;
+            for (int p = 0; p < numProp; ++p) {
+              strided_range<VecDH<float>::IterC> oldRange(
+                  oldProp.begin() + p, oldProp.end(), numProp);
+              strided_range<VecDH<float>::Iter> newRange(
+                  newProp.begin() + numPropOut * propVertIndices[i] + p,
+                  newProp.end(), numPropOut);
+              copy(policy, oldRange.begin(), oldRange.end(), newRange.begin());
+            }
+          } else {
+            // point all triangles at single new property of zeros.
+            fill(policy, start, start + node->pImpl_->NumTri(),
+                 glm::ivec3(propVertIndices[i]));
+          }
+        }
+
         if (node->transform_ == glm::mat4x3(1.0f)) {
           copy(policy, node->pImpl_->vertPos_.begin(),
                node->pImpl_->vertPos_.end(),
@@ -195,11 +250,6 @@ Manifold::Impl CsgLeafNode::Compose(
           copy(policy, node->pImpl_->faceNormal_.begin(),
                node->pImpl_->faceNormal_.end(),
                combined.faceNormal_.begin() + triIndices[i]);
-          transform(
-              policy, node->pImpl_->halfedge_.begin(),
-              node->pImpl_->halfedge_.end(),
-              combined.halfedge_.begin() + edgeIndices[i],
-              UpdateHalfedge({vertIndices[i], edgeIndices[i], triIndices[i]}));
         } else {
           // no need to apply the transform to the node, just copy the vertices
           // and face normals and apply transform on the fly
@@ -215,11 +265,6 @@ Manifold::Impl CsgLeafNode::Compose(
           copy_n(policy, faceNormalBegin, node->pImpl_->faceNormal_.size(),
                  combined.faceNormal_.begin() + triIndices[i]);
 
-          transform(
-              policy, node->pImpl_->halfedge_.begin(),
-              node->pImpl_->halfedge_.end(),
-              combined.halfedge_.begin() + edgeIndices[i],
-              UpdateHalfedge({vertIndices[i], edgeIndices[i], triIndices[i]}));
           const bool invert = glm::determinant(glm::mat3(node->transform_)) < 0;
           for_each_n(policy,
                      zip(combined.halfedgeTangent_.begin() + edgeIndices[i],

src/manifold/src/edge_op.cpp

@@ -49,7 +49,7 @@ struct ShortEdge {
   const glm::vec3* vertPos;
   const float precision;
 
-  __host__ __device__ bool operator()(int edge) {
+  __host__ __device__ bool operator()(int edge) const {
     if (halfedge[edge].pairedHalfedge < 0) return false;
     // Flag short edges
     const glm::vec3 delta =
@@ -62,7 +62,7 @@ struct FlagEdge {
   const Halfedge* halfedge;
   const TriRef* triRef;
 
-  __host__ __device__ bool operator()(int edge) {
+  __host__ __device__ bool operator()(int edge) const {
     if (halfedge[edge].pairedHalfedge < 0) return false;
     // Flag redundant edges - those where the startVert is surrounded by only
     // two original triangles.
@@ -87,7 +87,7 @@ struct SwappableEdge {
   const glm::vec3* triNormal;
   const float precision;
 
-  __host__ __device__ bool operator()(int edge) {
+  __host__ __device__ bool operator()(int edge) const {
     if (halfedge[edge].pairedHalfedge < 0) return false;
 
     int tri = halfedge[edge].face;
@@ -110,6 +110,15 @@ struct SwappableEdge {
            Is01Longest(v[0], v[1], v[2]);
   }
 };
+
+struct SortEntry {
+  int start;
+  int end;
+  int index;
+  inline bool operator<(const SortEntry& other) const {
+    return start == other.start ? end < other.end : start < other.start;
+  }
+};
 }  // namespace
 
 namespace manifold {
@@ -136,21 +145,29 @@ namespace manifold {
 void Manifold::Impl::SimplifyTopology() {
   if (!halfedge_.size()) return;
 
+  int nbEdges = halfedge_.size();
+  int numFlagged = 0;
+  VecDH<uint8_t> bflags(nbEdges);
+  uint8_t* bflagsPtr = bflags.ptrD();
+
+  VecDH<SortEntry> entries(nbEdges);
+  SortEntry* entriesPtr = entries.ptrD();
+
   auto policy = autoPolicy(halfedge_.size());
+  for_each_n(policy, countAt(0), nbEdges, [=] __host__ __device__(int i) {
+    entriesPtr[i].start = halfedge_[i].startVert;
+    entriesPtr[i].end = halfedge_[i].endVert;
+    entriesPtr[i].index = i;
+  });
 
-  VecDH<Halfedge> halfedge(halfedge_);
-  VecDH<int> idx(halfedge_.size());
-  sequence(policy, idx.begin(), idx.end());
-  sort_by_key(policy, halfedge.begin(), halfedge.end(), idx.begin());
-
-  VecDH<int> flaggedEdges(halfedge_.size());
-
-  int numFlagged =
-      copy_if<decltype(flaggedEdges.begin())>(
-          policy, idx.begin(), idx.end() - 1, countAt(0), flaggedEdges.begin(),
-          DuplicateEdge({halfedge.cptrD()})) -
-      flaggedEdges.begin();
-  flaggedEdges.resize(numFlagged);
+  sort(policy, entries.begin(), entries.end());
+  for (int i = 0; i < nbEdges - 1; ++i) {
+    if (entries[i].start == entries[i + 1].start &&
+        entries[i].end == entries[i + 1].end) {
+      DedupeEdge(entries[i].index);
+      numFlagged++;
+    }
+  }
 
 #ifdef MANIFOLD_DEBUG
   if (ManifoldParams().verbose && numFlagged > 0) {
@@ -159,15 +176,18 @@ void Manifold::Impl::SimplifyTopology() {
   }
 #endif
 
-  for (const int edge : flaggedEdges) DedupeEdge(edge);
-
-  flaggedEdges.resize(halfedge_.size());
-  numFlagged =
-      copy_if<decltype(flaggedEdges.begin())>(
-          policy, countAt(0), countAt(halfedge_.size()), flaggedEdges.begin(),
-          ShortEdge({halfedge_.cptrD(), vertPos_.cptrD(), precision_})) -
-      flaggedEdges.begin();
-  flaggedEdges.resize(numFlagged);
+  {
+    numFlagged = 0;
+    ShortEdge se{halfedge_.cptrD(), vertPos_.cptrD(), precision_};
+    for_each_n(policy, countAt(0), nbEdges,
+               [=] __host__ __device__(int i) { bflagsPtr[i] = se(i); });
+    for (int i = 0; i < nbEdges; ++i) {
+      if (bflags[i]) {
+        CollapseEdge(i);
+        numFlagged++;
+      }
+    }
+  }
 
 #ifdef MANIFOLD_DEBUG
   if (ManifoldParams().verbose && numFlagged > 0) {
@@ -176,15 +196,18 @@ void Manifold::Impl::SimplifyTopology() {
   }
 #endif
 
-  for (const int edge : flaggedEdges) CollapseEdge(edge);
-
-  flaggedEdges.resize(halfedge_.size());
-  numFlagged =
-      copy_if<decltype(flaggedEdges.begin())>(
-          policy, countAt(0), countAt(halfedge_.size()), flaggedEdges.begin(),
-          FlagEdge({halfedge_.cptrD(), meshRelation_.triRef.cptrD()})) -
-      flaggedEdges.begin();
-  flaggedEdges.resize(numFlagged);
+  {
+    numFlagged = 0;
+    FlagEdge se{halfedge_.cptrD(), meshRelation_.triRef.cptrD()};
+    for_each_n(policy, countAt(0), nbEdges,
+               [=] __host__ __device__(int i) { bflagsPtr[i] = se(i); });
+    for (int i = 0; i < nbEdges; ++i) {
+      if (bflags[i]) {
+        CollapseEdge(i);
+        numFlagged++;
+      }
+    }
+  }
 
 #ifdef MANIFOLD_DEBUG
   if (ManifoldParams().verbose && numFlagged > 0) {
@@ -193,26 +216,25 @@ void Manifold::Impl::SimplifyTopology() {
   }
 #endif
 
-  for (const int edge : flaggedEdges) CollapseEdge(edge);
-
-  flaggedEdges.resize(halfedge_.size());
-  numFlagged =
-      copy_if<decltype(flaggedEdges.begin())>(
-          policy, countAt(0), countAt(halfedge_.size()), flaggedEdges.begin(),
-          SwappableEdge({halfedge_.cptrD(), vertPos_.cptrD(),
-                         faceNormal_.cptrD(), precision_})) -
-      flaggedEdges.begin();
-  flaggedEdges.resize(numFlagged);
+  {
+    numFlagged = 0;
+    SwappableEdge se{halfedge_.cptrD(), vertPos_.cptrD(), faceNormal_.cptrD(),
+                     precision_};
+    for_each_n(policy, countAt(0), nbEdges,
+               [=] __host__ __device__(int i) { bflagsPtr[i] = se(i); });
+    for (int i = 0; i < nbEdges; ++i) {
+      if (bflags[i]) {
+        RecursiveEdgeSwap(i);
+        numFlagged++;
+      }
+    }
+  }
 
 #ifdef MANIFOLD_DEBUG
   if (ManifoldParams().verbose && numFlagged > 0) {
     std::cout << "found " << numFlagged << " edges to swap" << std::endl;
   }
 #endif
-
-  for (const int edge : flaggedEdges) {
-    RecursiveEdgeSwap(edge);
-  }
 }
 
 void Manifold::Impl::DedupeEdge(const int edge) {
@@ -280,11 +302,13 @@ void Manifold::Impl::PairUp(int edge0, int edge1) {
 // (edgeEdge.endVert must == startEdge.endVert), updating each edge to point
 // to vert instead.
 void Manifold::Impl::UpdateVert(int vert, int startEdge, int endEdge) {
-  while (startEdge != endEdge) {
-    halfedge_[startEdge].endVert = vert;
-    startEdge = NextHalfedge(startEdge);
-    halfedge_[startEdge].startVert = vert;
-    startEdge = halfedge_[startEdge].pairedHalfedge;
+  int current = startEdge;
+  while (current != endEdge) {
+    halfedge_[current].endVert = vert;
+    current = NextHalfedge(current);
+    halfedge_[current].startVert = vert;
+    current = halfedge_[current].pairedHalfedge;
+    ASSERT(current != startEdge, logicErr, "infinite loop in decimator!");
   }
 }
 
@@ -561,4 +585,27 @@ void Manifold::Impl::RecursiveEdgeSwap(const int edge) {
   RecursiveEdgeSwap(halfedge_[tri0edge[1]].pairedHalfedge);
   RecursiveEdgeSwap(halfedge_[tri1edge[0]].pairedHalfedge);
 }
+
+void Manifold::Impl::SplitPinchedVerts() {
+  std::vector<bool> vertProcessed(NumVert(), false);
+  std::vector<bool> halfedgeProcessed(halfedge_.size(), false);
+  for (int i = 0; i < halfedge_.size(); ++i) {
+    if (halfedgeProcessed[i]) continue;
+    int vert = halfedge_[i].startVert;
+    if (vertProcessed[vert]) {
+      vertPos_.push_back(vertPos_[vert]);
+      vert = NumVert() - 1;
+    } else {
+      vertProcessed[vert] = true;
+    }
+    int current = i;
+    do {
+      halfedgeProcessed[current] = true;
+      halfedge_[current].startVert = vert;
+      current = halfedge_[current].pairedHalfedge;
+      halfedge_[current].endVert = vert;
+      current = NextHalfedge(current);
+    } while (current != i);
+  }
+}
 }  // namespace manifold

src/manifold/src/impl.cpp

@@ -29,6 +29,8 @@
 namespace {
 using namespace manifold;
 
+constexpr uint64_t kRemove = std::numeric_limits<uint64_t>::max();
+
 __host__ __device__ void AtomicAddVec3(glm::vec3& target,
                                        const glm::vec3& add) {
   for (int i : {0, 1, 2}) {
@@ -488,10 +490,25 @@ Manifold::Impl::Impl(const MeshGL& meshGL,
 Manifold::Impl::Impl(const Mesh& mesh, const MeshRelationD& relation,
                      const std::vector<float>& propertyTolerance,
                      bool hasFaceIDs)
-    : vertPos_(mesh.vertPos),
-      halfedgeTangent_(mesh.halfedgeTangent),
-      meshRelation_(relation) {
-  VecDH<glm::ivec3> triVerts = mesh.triVerts;
+    : vertPos_(mesh.vertPos), halfedgeTangent_(mesh.halfedgeTangent) {
+  meshRelation_ = {relation.originalID, relation.numProp, relation.properties,
+                   relation.meshIDtransform};
+
+  VecDH<glm::ivec3> triVerts;
+  for (int i = 0; i < mesh.triVerts.size(); ++i) {
+    const glm::ivec3 tri = mesh.triVerts[i];
+    // Remove topological degenerates
+    if (tri[0] != tri[1] && tri[1] != tri[2] && tri[2] != tri[0]) {
+      triVerts.push_back(tri);
+      if (relation.triRef.size() > 0) {
+        meshRelation_.triRef.push_back(relation.triRef[i]);
+      }
+      if (relation.triProperties.size() > 0) {
+        meshRelation_.triProperties.push_back(relation.triProperties[i]);
+      }
+    }
+  }
+
   if (!IsIndexInBounds(triVerts)) {
     MarkFailure(Error::VertexOutOfBounds);
     return;
@@ -510,6 +527,9 @@ Manifold::Impl::Impl(const Mesh& mesh, const MeshRelationD& relation,
     MarkFailure(Error::NotManifold);
     return;
   }
+
+  SplitPinchedVerts();
+
   CalculateNormals();
 
   InitializeOriginal();
@@ -661,12 +681,12 @@ void Manifold::Impl::CreateFaces(const std::vector<float>& propertyTolerance) {
  */
 void Manifold::Impl::CreateHalfedges(const VecDH<glm::ivec3>& triVerts) {
   const int numTri = triVerts.size();
-  const int numEdge = 3 * numTri / 2;
+  const int numHalfedge = 3 * numTri;
   // drop the old value first to avoid copy
   halfedge_.resize(0);
-  halfedge_.resize(2 * numEdge);
-  VecDH<uint64_t> edge(2 * numEdge);
-  VecDH<int> ids(2 * numEdge);
+  halfedge_.resize(numHalfedge);
+  VecDH<uint64_t> edge(numHalfedge);
+  VecDH<int> ids(numHalfedge);
   auto policy = autoPolicy(numTri);
   sequence(policy, ids.begin(), ids.end());
   for_each_n(policy, zip(countAt(0), triVerts.begin()), numTri,
@@ -680,8 +700,8 @@ void Manifold::Impl::CreateHalfedges(const VecDH<glm::ivec3>& triVerts) {
   // Once sorted, the first half of the range is the forward halfedges, which
   // correspond to their backward pair at the same offset in the second half
   // of the range.
-  for_each_n(policy, countAt(0), numEdge,
-             LinkHalfedges({halfedge_.ptrD(), ids.ptrD(), numEdge}));
+  for_each_n(policy, countAt(0), numHalfedge / 2,
+             LinkHalfedges({halfedge_.ptrD(), ids.ptrD(), numHalfedge / 2}));
 }
 
 /**

src/manifold/src/impl.h

@@ -37,10 +37,10 @@ struct Manifold::Impl {
     /// The originalID of this Manifold if it is an original; -1 otherwise.
     int originalID = -1;
     int numProp = 0;
-    VecDH<TriRef> triRef;
-    VecDH<glm::ivec3> triProperties;
     VecDH<float> properties;
     std::map<int, Relation> meshIDtransform;
+    VecDH<TriRef> triRef;
+    VecDH<glm::ivec3> triProperties;
   };
 
   Box bBox_;
@@ -127,6 +127,7 @@ struct Manifold::Impl {
   void UpdateVert(int vert, int startEdge, int endEdge);
   void FormLoop(int current, int end);
   void CollapseTri(const glm::ivec3& triEdge);
+  void SplitPinchedVerts();
 
   // smoothing.cu
   void CreateTangents(const std::vector<Smoothness>&);

src/manifold/src/manifold.cpp

@@ -12,6 +12,7 @@
 // See the License for the specific language governing permissions and
 // limitations under the License.
 
+#include <QuickHull.hpp>
 #include <algorithm>
 #include <map>
 #include <numeric>
@@ -284,7 +285,9 @@ MeshGL Manifold::GetMeshGL(glm::ivec3 normalIdx) const {
 
   // Duplicate verts with different props
   std::vector<int> vert2idx(impl.NumVert(), -1);
-  std::map<std::pair<int, int>, int> vertPropPair;
+  std::vector<std::vector<glm::ivec2>> vertPropPair(impl.NumVert());
+  out.vertProperties.reserve(numVert * out.numProp);
+
   for (int run = 0; run < out.runOriginalID.size(); ++run) {
     for (int tri = out.runIndex[run] / 3; tri < out.runIndex[run + 1] / 3;
          ++tri) {
@@ -294,14 +297,19 @@ MeshGL Manifold::GetMeshGL(glm::ivec3 normalIdx) const {
         const int prop = triProp[i];
         const int vert = out.triVerts[3 * tri + i];
 
-        const auto it = vertPropPair.find({vert, prop});
-        if (it != vertPropPair.end()) {
-          out.triVerts[3 * tri + i] = it->second;
-          continue;
+        auto& bin = vertPropPair[vert];
+        bool bFound = false;
+        for (int k = 0; k < bin.size(); ++k) {
+          if (bin[k].x == prop) {
+            bFound = true;
+            out.triVerts[3 * tri + i] = bin[k].y;
+            break;
+          }
         }
+        if (bFound) continue;
         const int idx = out.vertProperties.size() / out.numProp;
-        vertPropPair.insert({{vert, prop}, idx});
         out.triVerts[3 * tri + i] = idx;
+        bin.push_back({prop, idx});
 
         for (int p : {0, 1, 2}) {
           out.vertProperties.push_back(impl.vertPos_[vert][p]);
@@ -310,6 +318,7 @@ MeshGL Manifold::GetMeshGL(glm::ivec3 normalIdx) const {
           out.vertProperties.push_back(
               impl.meshRelation_.properties[prop * numProp + p]);
         }
+
         if (updateNormals) {
           glm::vec3 normal;
           const int start = out.vertProperties.size() - out.numProp;
@@ -767,4 +776,50 @@ Manifold Manifold::TrimByPlane(glm::vec3 normal, float originOffset) const {
 }
 
 ExecutionParams& ManifoldParams() { return params; }
+
+/**
+ * Compute the convex hull of a set of points. If the given points are fewer
+ * than 4, or they are all coplanar, an empty Manifold will be returned.
+ *
+ * @param pts A vector of 3-dimensional points over which to compute a convex
+ * hull.
+ */
+Manifold Manifold::Hull(const std::vector<glm::vec3>& pts) {
+  const int numVert = pts.size();
+  if (numVert < 4) return Manifold();
+
+  std::vector<quickhull::Vector3<double>> vertices(numVert);
+  for (int i = 0; i < numVert; i++) {
+    vertices[i] = {pts[i].x, pts[i].y, pts[i].z};
+  }
+
+  quickhull::QuickHull<double> qh;
+  // bools: correct triangle winding, and use original indices
+  auto hull = qh.getConvexHull(vertices, false, true);
+  const auto& triangles = hull.getIndexBuffer();
+  const int numTris = triangles.size() / 3;
+
+  Mesh mesh;
+  mesh.vertPos = pts;
+  mesh.triVerts.reserve(numTris);
+  for (int i = 0; i < numTris; i++) {
+    const int j = i * 3;
+    mesh.triVerts.push_back({triangles[j], triangles[j + 1], triangles[j + 2]});
+  }
+  return Manifold(mesh);
+}
+
+/**
+ * Compute the convex hull of this manifold.
+ */
+Manifold Manifold::Hull() const { return Hull(GetMesh().vertPos); }
+
+/**
+ * Compute the convex hull enveloping a set of manifolds.
+ *
+ * @param manifolds A vector of manifolds over which to compute a convex hull.
+ */
+Manifold Manifold::Hull(const std::vector<Manifold>& manifolds) {
+  return Compose(manifolds).Hull();
+}
 }  // namespace manifold

src/manifold/src/properties.cpp

@@ -197,9 +197,8 @@ struct CheckHalfedges {
 
   __host__ __device__ bool operator()(int edge) {
     const Halfedge halfedge = halfedges[edge];
-    if (halfedge.startVert == -1 && halfedge.endVert == -1 &&
-        halfedge.pairedHalfedge == -1)
-      return true;
+    if (halfedge.startVert == -1 && halfedge.endVert == -1) return true;
+    if (halfedge.pairedHalfedge == -1) return false;
 
     const Halfedge paired = halfedges[halfedge.pairedHalfedge];
     bool good = true;

src/manifold/src/sort.cpp

@@ -312,8 +312,9 @@ void Manifold::Impl::SortVerts() {
 
   VecDH<int> vertNew2Old(numVert);
   sequence(policy, vertNew2Old.begin(), vertNew2Old.end());
+
   sort_by_key(policy, vertMorton.begin(), vertMorton.end(),
-              zip(vertPos_.begin(), vertNew2Old.begin()));
+              vertNew2Old.begin());
 
   ReindexVerts(vertNew2Old, numVert);
 
@@ -323,10 +324,12 @@ void Manifold::Impl::SortVerts() {
       lower_bound<decltype(vertMorton.begin())>(policy, vertMorton.begin(),
                                                 vertMorton.end(), kNoCode) -
       vertMorton.begin();
-  vertPos_.resize(newNumVert);
+
+  vertNew2Old.resize(newNumVert);
+  Permute(vertPos_, vertNew2Old);
+
   if (vertNormal_.size() == numVert) {
     Permute(vertNormal_, vertNew2Old);
-    vertNormal_.resize(newNumVert);
   }
 }
 
@@ -395,7 +398,7 @@ void Manifold::Impl::SortFaces(VecDH<Box>& faceBox,
   sequence(policy, faceNew2Old.begin(), faceNew2Old.end());
 
   sort_by_key(policy, faceMorton.begin(), faceMorton.end(),
-              zip(faceBox.begin(), faceNew2Old.begin()));
+              faceNew2Old.begin());
 
   // Tris were flagged for removal with pairedHalfedge = -1 and assigned kNoCode
   // to sort them to the end, which allows them to be removed.
@@ -403,10 +406,10 @@ void Manifold::Impl::SortFaces(VecDH<Box>& faceBox,
       find<decltype(faceMorton.begin())>(policy, faceMorton.begin(),
                                          faceMorton.end(), kNoCode) -
       faceMorton.begin();
-  faceBox.resize(newNumTri);
   faceMorton.resize(newNumTri);
   faceNew2Old.resize(newNumTri);
 
+  Permute(faceBox, faceNew2Old);
   GatherFaces(faceNew2Old);
 }
 

src/third_party/CMakeLists.txt

@@ -6,5 +6,8 @@ set(CLIPPER2_EXAMPLES OFF)
 set(CLIPPER2_TESTS OFF)
 set(CLIPPER2_USINGZ "OFF" CACHE STRING "Preempt cache default of USINGZ (we only use 2d)")
 add_definitions(-D_HAS_EXCEPTIONS=0) # disable exceptions for STL
-
 add_subdirectory(clipper2/CPP)
+
+add_library(quickhull quickhull/QuickHull.cpp)
+target_include_directories(quickhull PUBLIC quickhull)
+target_compile_features(quickhull PUBLIC cxx_std_17)

src/utilities/include/vec_dh.h

@@ -166,8 +166,9 @@ class ManagedVec {
   void shrink_to_fit() {
     T *newBuffer = nullptr;
     if (size_ > 0) {
-      mallocManaged(&newBuffer, size_ * sizeof(T));
-      prefetch(newBuffer, size_ * sizeof(T), onHost);
+      int n_bytes = size_ * sizeof(T);
+      mallocManaged(&newBuffer, n_bytes);
+      prefetch(newBuffer, n_bytes, onHost);
       uninitialized_copy(autoPolicy(size_), ptr_, ptr_ + size_, newBuffer);
     }
     freeManaged(ptr_);
@@ -227,19 +228,9 @@ class ManagedVec {
 
   const T &back() const { return *(ptr_ + size_ - 1); }
 
-  T &operator[](size_t i) {
-    // ASSERT(i < size_, logicErr,
-    //        "Index out of bounds: " + std::to_string(i) +
-    //            ", size: " + std::to_string(size_));
-    return *(ptr_ + i);
-  }
+  T &operator[](size_t i) { return *(ptr_ + i); }
 
-  const T &operator[](size_t i) const {
-    // ASSERT(i < size_, logicErr,
-    //        "Index out of bounds: " + std::to_string(i) +
-    //            ", size: " + std::to_string(size_));
-    return *(ptr_ + i);
-  }
+  const T &operator[](size_t i) const { return *(ptr_ + i); }
 
   bool empty() const { return size_ == 0; }
 
@@ -252,6 +243,11 @@ class ManagedVec {
   static constexpr int DEVICE_MAX_BYTES = 1 << 16;
 
   static void mallocManaged(T **ptr, size_t bytes) {
+    // only exists to please the compiler
+    // see https://github.com/elalish/manifold/pull/496 for background
+    if (bytes >= (1ull << 63)) {
+      throw std::bad_alloc();
+    }
 #ifdef MANIFOLD_USE_CUDA
     if (CudaEnabled())
       cudaMallocManaged(reinterpret_cast<void **>(ptr), bytes);
@@ -399,15 +395,9 @@ class VecDH {
 
   const T *ptrH() const { return cptrH(); }
 
-  T &operator[](int i) {
-    impl_.prefetch_to(true);
-    return impl_[i];
-  }
+  T &operator[](int i) { return impl_[i]; }
 
-  const T &operator[](int i) const {
-    impl_.prefetch_to(true);
-    return impl_[i];
-  }
+  const T &operator[](int i) const { return impl_[i]; }
 
   T &back() { return impl_.back(); }
 

test/boolean_test.cpp

@@ -530,6 +530,44 @@ TEST(Boolean, Cubes) {
 #endif
 }
 
+TEST(Boolean, NoRetainedVerts) {
+  Manifold cube = Manifold::Cube(glm::vec3(1), true);
+  Manifold oct = Manifold::Sphere(1, 4);
+  EXPECT_NEAR(cube.GetProperties().volume, 1, 0.001);
+  EXPECT_NEAR(oct.GetProperties().volume, 1.333, 0.001);
+  EXPECT_NEAR((cube ^ oct).GetProperties().volume, 0.833, 0.001);
+}
+
+TEST(Boolean, PropertiesNoIntersection) {
+  MeshGL cubeUV = CubeUV();
+  Manifold m0(cubeUV);
+  Manifold m1 = m0.Translate(glm::vec3(1.5));
+  Manifold result = m0 + m1;
+  EXPECT_EQ(result.NumProp(), 2);
+  RelatedGL(result, {cubeUV});
+}
+
+TEST(Boolean, MixedProperties) {
+  MeshGL cubeUV = CubeUV();
+  Manifold m0(cubeUV);
+  Manifold m1 = Manifold::Cube();
+  Manifold result = m0 + m1.Translate(glm::vec3(0.5));
+  EXPECT_EQ(result.NumProp(), 2);
+  RelatedGL(result, {cubeUV, m1.GetMeshGL()});
+}
+
+TEST(Boolean, MixedNumProp) {
+  MeshGL cubeUV = CubeUV();
+  Manifold m0(cubeUV);
+  Manifold m1 = Manifold::Cube();
+  Manifold result =
+      m0 + m1.SetProperties(1, [](float* prop, glm::vec3 p, const float* n) {
+               prop[0] = 1;
+             }).Translate(glm::vec3(0.5));
+  EXPECT_EQ(result.NumProp(), 2);
+  RelatedGL(result, {cubeUV, m1.GetMeshGL()});
+}
+
 TEST(Boolean, Subtract) {
   Mesh firstMesh;
   firstMesh.vertPos = {{0, 0, 0},           {1540, 0, 0},

test/cross_section_test.cpp

@@ -16,6 +16,8 @@
 
 #include <gtest/gtest.h>
 
+#include <vector>
+
 #include "glm/geometric.hpp"
 #include "manifold.h"
 #include "polygon.h"
@@ -72,13 +74,17 @@ TEST(CrossSection, Empty) {
   EXPECT_TRUE(e.IsEmpty());
 }
 
-TEST(CrossSection, RectClip) {
-  auto sq = CrossSection::Square({10, 10});
-  auto rect = Rect({0, 0}, {10, 5});
-  auto clipped = sq.RectClip(rect);
+TEST(CrossSection, Rect) {
+  float w = 10;
+  float h = 5;
+  auto rect = Rect({0, 0}, {w, h});
+  auto cross = rect.AsCrossSection();
+  auto area = rect.Area();
 
-  EXPECT_EQ(sq.Area() / 2, clipped.Area());
+  EXPECT_FLOAT_EQ(area, w * h);
+  EXPECT_FLOAT_EQ(area, cross.Area());
   EXPECT_TRUE(rect.Contains({5, 5}));
+  EXPECT_TRUE(rect.Contains(cross.Bounds()));
   EXPECT_TRUE(rect.Contains(Rect()));
   EXPECT_TRUE(rect.DoesOverlap(Rect({5, 5}, {15, 15})));
   EXPECT_TRUE(Rect().IsEmpty());
@@ -141,3 +147,47 @@ TEST(CrossSection, Decompose) {
   Identical(Manifold::Extrude(ab, 1.).GetMesh(),
             Manifold::Extrude(recomp, 1.).GetMesh());
 }
+
+TEST(CrossSection, FillRule) {
+  SimplePolygon polygon = {
+      {-7, 13},   //
+      {-7, 12},   //
+      {-5, 9},    //
+      {-5, 8.1},  //
+      {-4.8, 8},  //
+  };
+
+  CrossSection positive(polygon);
+  EXPECT_NEAR(positive.Area(), 0.683, 0.001);
+
+  CrossSection negative(polygon, CrossSection::FillRule::Negative);
+  EXPECT_NEAR(negative.Area(), 0.193, 0.001);
+
+  CrossSection evenOdd(polygon, CrossSection::FillRule::EvenOdd);
+  EXPECT_NEAR(evenOdd.Area(), 0.875, 0.001);
+
+  CrossSection nonZero(polygon, CrossSection::FillRule::NonZero);
+  EXPECT_NEAR(nonZero.Area(), 0.875, 0.001);
+}
+
+TEST(CrossSection, Hull) {
+  auto circ = CrossSection::Circle(10, 360);
+  auto circs = std::vector<CrossSection>{circ, circ.Translate({0, 30}),
+                                         circ.Translate({30, 0})};
+  auto circ_tri = CrossSection::Hull(circs);
+  auto centres = SimplePolygon{{0, 0}, {0, 30}, {30, 0}, {15, 5}};
+  auto tri = CrossSection::Hull(centres);
+
+#ifdef MANIFOLD_EXPORT
+  if (options.exportModels) {
+    auto circ_tri_ex = Manifold::Extrude(circ_tri, 10);
+    ExportMesh("cross_section_hull_circ_tri.glb", circ_tri_ex.GetMesh(), {});
+  }
+#endif
+
+  auto circ_area = circ.Area();
+  EXPECT_FLOAT_EQ(circ_area, (circ - circ.Scale({0.8, 0.8})).Hull().Area());
+  EXPECT_FLOAT_EQ(
+      circ_area * 2.5,
+      (CrossSection::BatchBoolean(circs, OpType::Add) - tri).Area());
+}

test/manifold_test.cpp

@@ -197,13 +197,28 @@ TEST(Manifold, ExtrudeCone) {
   EXPECT_FLOAT_EQ(donut.GetProperties().volume, 4.0f);
 }
 
+Polygons RotatePolygons(Polygons polys, const int index) {
+  Polygons rotatedPolys;
+  for (auto& polygon : polys) {
+    auto rotatedPolygon = polygon;
+    std::rotate(rotatedPolygon.begin(), rotatedPolygon.begin() + index,
+                rotatedPolygon.end());
+    rotatedPolys.push_back(rotatedPolygon);
+  }
+  return rotatedPolys;
+}
+
 TEST(Manifold, Revolve) {
   Polygons polys = SquareHole();
-  Manifold vug = Manifold::Revolve(polys, 48);
-  EXPECT_EQ(vug.Genus(), -1);
-  auto prop = vug.GetProperties();
-  EXPECT_NEAR(prop.volume, 14.0f * glm::pi<float>(), 0.2f);
-  EXPECT_NEAR(prop.surfaceArea, 30.0f * glm::pi<float>(), 0.2f);
+  Manifold vug;
+  for (int i = 0; i < polys[0].size(); i++) {
+    Polygons rotatedPolys = RotatePolygons(polys, i);
+    vug = Manifold::Revolve(rotatedPolys, 48);
+    EXPECT_EQ(vug.Genus(), -1);
+    auto prop = vug.GetProperties();
+    EXPECT_NEAR(prop.volume, 14.0f * glm::pi<float>(), 0.2f);
+    EXPECT_NEAR(prop.surfaceArea, 30.0f * glm::pi<float>(), 0.2f);
+  }
 }
 
 TEST(Manifold, Revolve2) {
@@ -215,6 +230,38 @@ TEST(Manifold, Revolve2) {
   EXPECT_NEAR(prop.surfaceArea, 96.0f * glm::pi<float>(), 1.0f);
 }
 
+TEST(Manifold, PartialRevolveOnYAxis) {
+  Polygons polys = SquareHole(2.0f);
+  Polygons offsetPolys = SquareHole(10.0f);
+
+  Manifold revolute;
+  for (int i = 0; i < polys[0].size(); i++) {
+    Polygons rotatedPolys = RotatePolygons(polys, i);
+    revolute = Manifold::Revolve(rotatedPolys, 48, 180);
+    EXPECT_EQ(revolute.Genus(), 1);
+    auto prop = revolute.GetProperties();
+    EXPECT_NEAR(prop.volume, 24.0f * glm::pi<float>(), 1.0f);
+    EXPECT_NEAR(
+        prop.surfaceArea,
+        48.0f * glm::pi<float>() + 4.0f * 4.0f * 2.0f - 2.0f * 2.0f * 2.0f,
+        1.0f);
+  }
+}
+
+TEST(Manifold, PartialRevolveOffset) {
+  Polygons polys = SquareHole(10.0f);
+
+  Manifold revolute;
+  for (int i = 0; i < polys[0].size(); i++) {
+    Polygons rotatedPolys = RotatePolygons(polys, i);
+    revolute = Manifold::Revolve(rotatedPolys, 48, 180);
+    auto prop = revolute.GetProperties();
+    EXPECT_EQ(revolute.Genus(), 1);
+    EXPECT_NEAR(prop.surfaceArea, 777.0f, 1.0f);
+    EXPECT_NEAR(prop.volume, 376.0f, 1.0f);
+  }
+}
+
 TEST(Manifold, Warp) {
   CrossSection square = CrossSection::Square({1, 1});
   Manifold shape = Manifold::Extrude(square, 2, 10).Warp([](glm::vec3& v) {
@@ -615,3 +662,98 @@ TEST(Manifold, MultiCompose) {
                          part.Mirror({1, 0, 0}).Translate({10, 10, 0})});
   EXPECT_FLOAT_EQ(finalAssembly.GetProperties().volume, 4000);
 }
+
+TEST(Manifold, MergeDegenerates) {
+  MeshGL cube = Manifold::Cube(glm::vec3(1), true).GetMeshGL();
+  MeshGL squash;
+  squash.vertProperties = cube.vertProperties;
+  squash.triVerts = cube.triVerts;
+  // Move one vert to the position of its neighbor and remove one triangle
+  // linking them to break the manifold.
+  squash.vertProperties[squash.vertProperties.size() - 1] *= -1;
+  squash.triVerts.resize(squash.triVerts.size() - 3);
+  // Rotate the degenerate triangle to the middle to catch more problems.
+  std::rotate(squash.triVerts.begin(), squash.triVerts.begin() + 3 * 5,
+              squash.triVerts.end());
+  // Merge should remove the now duplicate vertex.
+  EXPECT_TRUE(squash.Merge());
+  // Manifold should remove the triangle with two references to the same vert.
+  Manifold squashed = Manifold(squash);
+  EXPECT_FALSE(squashed.IsEmpty());
+  EXPECT_EQ(squashed.Status(), Manifold::Error::NoError);
+}
+
+TEST(Manifold, PinchedVert) {
+  Mesh shape;
+  shape.vertPos = {{0, 0, 0},         //
+                   {1, 1, 0},         //
+                   {1, -1, 0},        //
+                   {-0.00001, 0, 0},  //
+                   {-1, -1, -0},      //
+                   {-1, 1, 0},        //
+                   {0, 0, 2},         //
+                   {0, 0, -2}};
+  shape.triVerts = {{0, 2, 6},  //
+                    {2, 1, 6},  //
+                    {1, 0, 6},  //
+                    {4, 3, 6},  //
+                    {3, 5, 6},  //
+                    {5, 4, 6},  //
+                    {2, 0, 4},  //
+                    {0, 3, 4},  //
+                    {3, 0, 1},  //
+                    {3, 1, 5},  //
+                    {7, 2, 4},  //
+                    {7, 4, 5},  //
+                    {7, 5, 1},  //
+                    {7, 1, 2}};
+  Manifold touch(shape);
+  EXPECT_FALSE(touch.IsEmpty());
+  EXPECT_EQ(touch.Status(), Manifold::Error::NoError);
+  EXPECT_EQ(touch.Genus(), 0);
+}
+
+TEST(Manifold, TictacHull) {
+  const float tictacRad = 100;
+  const float tictacHeight = 500;
+  const int tictacSeg = 1000;
+  const float tictacMid = tictacHeight - 2 * tictacRad;
+  const auto sphere = Manifold::Sphere(tictacRad, tictacSeg);
+  const std::vector<Manifold> spheres{sphere,
+                                      sphere.Translate({0, 0, tictacMid})};
+  const auto tictac = Manifold::Hull(spheres);
+
+#ifdef MANIFOLD_EXPORT
+  if (options.exportModels) {
+    ExportMesh("tictac_hull.glb", tictac.GetMesh(), {});
+  }
+#endif
+
+  EXPECT_EQ(sphere.NumVert() + tictacSeg, tictac.NumVert());
+}
+
+TEST(Manifold, HollowHull) {
+  auto sphere = Manifold::Sphere(100, 360);
+  auto hollow = sphere - sphere.Scale({0.8, 0.8, 0.8});
+  const float sphere_vol = sphere.GetProperties().volume;
+  EXPECT_FLOAT_EQ(hollow.Hull().GetProperties().volume, sphere_vol);
+}
+
+TEST(Manifold, CubeHull) {
+  std::vector<glm::vec3> cubePts = {
+      {0, 0, 0},       {1, 0, 0},   {0, 1, 0},      {0, 0, 1},  // corners
+      {1, 1, 0},       {0, 1, 1},   {1, 0, 1},      {1, 1, 1},  // corners
+      {0.5, 0.5, 0.5}, {0.5, 0, 0}, {0.5, 0.7, 0.2}  // internal points
+  };
+  auto cube = Manifold::Hull(cubePts);
+  EXPECT_FLOAT_EQ(cube.GetProperties().volume, 1);
+}
+
+TEST(Manifold, EmptyHull) {
+  const std::vector<glm::vec3> tooFew{{0, 0, 0}, {1, 0, 0}, {0, 1, 0}};
+  EXPECT_TRUE(Manifold::Hull(tooFew).IsEmpty());
+
+  const std::vector<glm::vec3> coplanar{
+      {0, 0, 0}, {1, 0, 0}, {0, 1, 0}, {1, 1, 0}};
+  EXPECT_TRUE(Manifold::Hull(coplanar).IsEmpty());
+}

test/test.h

@@ -41,6 +41,7 @@ MeshGL CubeSTL();
 MeshGL WithIndexColors(const Mesh& in);
 MeshGL WithPositionColors(const Manifold& in);
 MeshGL WithNormals(const Manifold& in);
+MeshGL CubeUV();
 float GetMaxProperty(const MeshGL& mesh, int channel);
 float GetMinProperty(const MeshGL& mesh, int channel);
 void Identical(const Mesh& mesh1, const Mesh& mesh2);

test/test_main.cpp

@@ -249,6 +249,31 @@ MeshGL WithNormals(const Manifold& in) {
   return out;
 }
 
+MeshGL CubeUV() {
+  MeshGL mgl;
+  mgl.numProp = 5;
+  mgl.vertProperties = {0.5,  -0.5, 0.5,  0.5,  0.66,  //
+                        -0.5, -0.5, 0.5,  0.25, 0.66,  //
+                        0.5,  0.5,  0.5,  0.5,  0.33,  //
+                        -0.5, 0.5,  0.5,  0.25, 0.33,  //
+                        -0.5, -0.5, -0.5, 1.0,  0.66,  //
+                        0.5,  -0.5, -0.5, 0.75, 0.66,  //
+                        -0.5, 0.5,  -0.5, 1.0,  0.33,  //
+                        0.5,  0.5,  -0.5, 0.75, 0.33,  //
+                        -0.5, -0.5, -0.5, 0.0,  0.66,  //
+                        -0.5, 0.5,  -0.5, 0.0,  0.33,  //
+                        -0.5, 0.5,  -0.5, 0.25, 0.0,   //
+                        0.5,  0.5,  -0.5, 0.5,  0.0,   //
+                        -0.5, -0.5, -0.5, 0.25, 1.0,   //
+                        0.5,  -0.5, -0.5, 0.5,  1.0};
+  mgl.triVerts = {3, 1, 0, 3, 0, 2, 7,  5,  4, 7,  4, 6, 2, 0, 5,  2, 5,  7,
+                  9, 8, 1, 9, 1, 3, 11, 10, 3, 11, 3, 2, 0, 1, 12, 0, 12, 13};
+  mgl.mergeFromVert = {8, 12, 13, 9, 10, 11};
+  mgl.mergeToVert = {4, 4, 5, 6, 6, 7};
+  mgl.runOriginalID.push_back(Manifold::ReserveIDs(1));
+  return mgl;
+}
+
 float GetMaxProperty(const MeshGL& mesh, int channel) {
   float max = -std::numeric_limits<float>::infinity();
   const int numVert = mesh.NumVert();
@@ -344,7 +369,7 @@ void RelatedGL(const Manifold& out, const std::vector<MeshGL>& originals,
           ASSERT_NEAR(glm::length(normal), 1, 0.0001);
           ASSERT_GT(glm::dot(normal, outNormal), 0);
         } else {
-          for (int p = 3; p < output.numProp; ++p) {
+          for (int p = 3; p < inMesh.numProp; ++p) {
             const float propOut =
                 output.vertProperties[vert * output.numProp + p];