v2.3 (#685)

bumped version

.github/workflows/build_wheels.yml

@@ -1,6 +1,9 @@
 name: build_wheels
 
-on: [release, workflow_dispatch]
+on: 
+  release:
+    types: [published]
+  workflow_dispatch:
 
 jobs:
   build_wheels:
@@ -38,7 +41,7 @@ jobs:
       with:
         path: dist/*.tar.gz
 
-  pypi-publish-beta:
+  pypi-publish:
     name: Upload release to PyPI
     needs: [build_wheels, make_sdist]
     runs-on: ubuntu-latest

.github/workflows/check_format.yml

@@ -17,5 +17,8 @@ jobs:
         source: '.'
         exclude: '*/third_party'
         extensions: 'h,cpp,js,ts,html'
-        clangFormatVersion: 12
-
+        clangFormatVersion: 11
+    - uses: psf/black@stable
+      with:
+        options: "--check --verbose"
+        src: "./bindings/python/examples"

.github/workflows/manifold.yml

@@ -22,8 +22,8 @@ jobs:
     - name: Install dependencies
       run: |
         sudo apt-get -y update
-        DEBIAN_FRONTEND=noninteractive sudo apt install -y libglm-dev libgtest-dev libassimp-dev git libtbb-dev pkg-config libpython3-dev python3 python3-distutils python3-pip lcov
-        pip install trimesh
+        DEBIAN_FRONTEND=noninteractive sudo apt install -y libglm-dev libgtest-dev libassimp-dev git libtbb-dev libthrust-dev pkg-config libpython3-dev python3 python3-distutils python3-pip lcov
+        pip install trimesh pytest
     - uses: actions/checkout@v3
       with:
         submodules: recursive
@@ -42,6 +42,7 @@ jobs:
         ./manifold_test
         cd ../../
         python3 bindings/python/examples/run_all.py -e
+        python3 -m pytest
     - name: Coverage Report
       # only do code coverage for default sequential backend, it seems that TBB
       # backend will cause failure
@@ -105,8 +106,8 @@ jobs:
         # setup emscripten
         git clone https://github.com/emscripten-core/emsdk.git
         cd emsdk
-        ./emsdk install latest
-        ./emsdk activate latest
+        ./emsdk install 3.1.45
+        ./emsdk activate 3.1.45
     - uses: jwlawson/actions-setup-cmake@v1.12
     - name: Build WASM
       run: |

.github/workflows/publish_npm.yml

@@ -0,0 +1,34 @@
+name: publish_npm
+
+on: 
+  release:
+    types: [published]
+  workflow_dispatch:
+
+jobs:
+  publish_npm:
+    runs-on: ubuntu-22.04
+
+    steps:
+    - uses: actions/checkout@v3
+
+    - name: Setup Node
+      uses: actions/setup-node@v2
+      with:
+        node-version: '20.x'
+        registry-url: 'https://registry.npmjs.org'
+
+    - name: Download built examples
+      uses: dawidd6/action-download-artifact@v2
+      with:
+        workflow: manifold.yml
+        workflow_conclusion: completed
+        name: wasm
+        path: ./bindings/wasm/
+
+    - name: Publish to npm
+      run: |
+        cd ./bindings/wasm/
+        npm publish
+      env:
+        NODE_AUTH_TOKEN: ${{secrets.NPM_TOKEN}}

CMakeLists.txt

@@ -21,12 +21,25 @@ set(MANIFOLD_PAR "NONE" CACHE STRING "Parallel backend, either \"TBB\" or \"NONE
 set(MANIFOLD_FLAGS "" CACHE STRING "Manifold compiler flags")
 option(MANIFOLD_EXPORT "Build mesh export (via assimp) utility library" OFF)
 option(MANIFOLD_TEST "Enable testing suite" ON)
+# fuzztest is a rather large dependency
+option(MANIFOLD_FUZZ "Enable fuzzing tests" OFF)
 option(MANIFOLD_DEBUG "Enable debug tracing/timing" OFF)
-
 option(MANIFOLD_PYBIND "Build python bindings" ON)
 option(MANIFOLD_CBIND "Build C (FFI) bindings" OFF)
 option(MANIFOLD_JSBIND "Build js binding" ${EMSCRIPTEN})
 
+if(MANIFOLD_FUZZ)
+  if(NOT CMAKE_CXX_COMPILER_ID STREQUAL "Clang")
+    message(FATAL_ERROR "fuzztest only supports clang" )
+  endif()
+  # we should enable debug checks
+  set(MANIFOLD_DEBUG ON)
+  # enable fuzztest fuzzing mode
+  set(FUZZTEST_FUZZING_MODE ON)
+  # address sanitizer required
+  set(CMAKE_CXX_FLAGS "$CMAKE_CXX_FLAGS -fsanitize=address")
+endif()
+
 option(TRACY_ENABLE "Use tracy profiling" OFF) 
 option(TRACY_MEMORY_USAGE "Track memory allocation with tracy (expensive)" OFF)
 option(BUILD_TEST_CGAL "Build CGAL performance comparisons" OFF)

README.md

@@ -1,35 +1,36 @@
 [![codecov](https://codecov.io/github/elalish/manifold/branch/master/graph/badge.svg?token=IIA8G5HVS7)](https://codecov.io/github/elalish/manifold)
+[![PyPI version](https://badge.fury.io/py/manifold3d.svg)](https://badge.fury.io/py/manifold3d)
+[![npm version](https://badge.fury.io/js/manifold-3d.svg)](https://badge.fury.io/js/manifold-3d)
 [![twitter](https://img.shields.io/twitter/follow/manifoldcad?style=social&logo=twitter)](https://twitter.com/intent/follow?screen_name=manifoldcad)
 
 ## Users
 
 [OpenSCAD](https://openscad.org/), [IFCjs](https://ifcjs.github.io/info/), [Grid.Space](https://grid.space/), and [OCADml](https://github.com/OCADml/OManifold) have all integrated our Manifold geometry kernel! Why? Because its reliability is guaranteed and it's 1,000 times faster than other libraries. See our [usage](https://github.com/elalish/manifold/discussions/340) and [performance](https://github.com/elalish/manifold/discussions/383) discussions for all the latest and to add your own projects & analyses.
 
-For example, here is a log-log plot of Manifold's performance vs. earlier OpenSCAD geometry backends:
+## Manifold Frontend Sandboxes
 
-<img src="https://elalish.github.io/manifold/samples/models/perfSpheres.png" width="350px"/>
+### [ManifoldCAD.org](https://manifoldcad.org)
 
-## [ManifoldCAD.org](https://manifoldcad.org)
+If you like OpenSCAD / JSCAD, you might also like ManifoldCAD - our own solid modelling web app where you script in JS/TS. This uses our npm package, [manifold-3d](https://www.npmjs.com/package/manifold-3d), built via WASM. It's not quite as fast as our raw C++, but it's hard to beat for interoperability.
 
-If you like OpenSCAD / JSCAD, you might also like ManifoldCAD - our own solid modelling web app. Our WASM is not multithreaded yet, but it's still quite fast and a good way to test out our Manifold library.
+*Note for Firefox users: If you find the editor is stuck on **Loading...**, setting
+`dom.workers.modules.enabled: true` in your `about:config`, as mentioned in
+[issue#328](https://github.com/elalish/manifold/issues/328#issuecomment-1473847102)
+may solve the problem.*
+
+### [Python Colab Example](https://colab.research.google.com/drive/1VxrFYHPSHZgUbl9TeWzCeovlpXrPQ5J5?usp=sharing)
+
+If you prefer Python to JS/TS, make your own copy of the example notebook above. It demonstrates interop between our [`manifold3d`](https://pypi.org/project/manifold3d/) PyPI library and the popular [`trimesh`](https://pypi.org/project/trimesh/) library, including showing the interactive model right in the notebook and saving 3D model output.
 
 ![A metallic Menger sponge](https://elalish.github.io/manifold/samples/models/mengerSponge3.webp "A metallic Menger sponge")
 
-### Note for Firefox users
-
-If you find the editor is stuck on **Loading...**, setting
-`dom.workers.modules.enabled: true` in your `about:config`, as mentioned in the
-discussion of the
-[issue#328](https://github.com/elalish/manifold/issues/328#issuecomment-1473847102)
-of this repository may solve the problem.
-
 # Manifold
 
 [**API Documentation**](https://elalish.github.io/manifold/docs/html/modules.html) | [**Algorithm Documentation**](https://github.com/elalish/manifold/wiki/Manifold-Library) | [**Blog Posts**](https://elalish.blogspot.com/search/label/Manifold) | [**Web Examples**](https://elalish.github.io/manifold/model-viewer.html)
 
 [Manifold](https://github.com/elalish/manifold) is a geometry library dedicated to creating and operating on manifold triangle meshes. A [manifold mesh](https://github.com/elalish/manifold/wiki/Manifold-Library#manifoldness) is a mesh that represents a solid object, and so is very important in manufacturing, CAD, structural analysis, etc. Further information can be found on the [wiki](https://github.com/elalish/manifold/wiki/Manifold-Library).
 
-This is a modern C++ library that Github's CI verifies builds and runs on a variety of platforms. Additionally, we build bindings for JavaScript ([manifold-3d](https://www.npmjs.com/package/manifold-3d) on npm), Python, and C to make this library more portable and easy to use.
+This is a modern C++ library that Github's CI verifies builds and runs on a variety of platforms. Additionally, we build bindings for JavaScript ([manifold-3d](https://www.npmjs.com/package/manifold-3d) on npm), Python ([manifold3d](https://pypi.org/project/manifold3d/)), and C to make this library more portable and easy to use.
 
 System Dependencies (note that we will automatically download the dependency if there is no such package on the system):
 - [`GLM`](https://github.com/g-truc/glm/): A compact header-only vector library.
@@ -87,6 +88,10 @@ Offline building:
 The build instructions used by our CI are in [manifold.yml](https://github.com/elalish/manifold/blob/master/.github/workflows/manifold.yml), which is a good source to check if something goes wrong and for instructions specific to other platforms, like Windows.
 
 ### WASM
+
+> Note that we have only tested emscripten version 3.1.45. It is known that
+  3.1.48 has some issues compiling manifold.
+
 To build the JS WASM library, first install NodeJS and set up emscripten:
 
 (on Mac):
@@ -138,10 +143,30 @@ For more detailed documentation, please refer to the C++ API.
 
 Contributions are welcome! A lower barrier contribution is to simply make a PR that adds a test, especially if it repros an issue you've found. Simply name it prepended with DISABLED_, so that it passes the CI. That will be a very strong signal to me to fix your issue. However, if you know how to fix it yourself, then including the fix in your PR would be much appreciated!
 
-## Profiling
+### Formatting
+
+There is a formatting script `format.sh` that automatically formats everything.
+It requires clang-format 11 and black formatter for python.
+
+If you have clang-format installed but without clang-11, you can specify the
+clang-format executable by setting the `CLANG_FORMAT` environment variable.
+
+### Profiling
 
 There is now basic support for the [Tracy profiler](https://github.com/wolfpld/tracy) for our tests.
 To enable tracing, compile with `-DTRACY_ENABLE=on` cmake option, and run the test with Tracy server running.
+To enable memory profiling in addition to tracing, compile with `-DTRACY_MEMORY_USAGE=ON` in addition to `-DTRACY_ENABLE=ON`.
+
+### Fuzzing Support
+
+We use https://github.com/google/fuzztest for fuzzing the triangulator.
+
+To enable fuzzing, make sure that you are using clang compiler (`-DCMAKE_CXX_COMPILER=clang -DCMAKE_C_COMPILER=clang`), running Linux, and enable fuzzing support by setting `-DMANIFOLD_FUZZ=ON`.
+
+To run the fuzzer and minimize testcase, do
+```
+../minimizer.sh ./test/polygon_fuzz --fuzz=PolygonFuzz.TriangulationNoCrash
+```
 
 ## About the author
 

bindings/c/include/manifoldc.h

@@ -111,6 +111,12 @@ ManifoldManifold *manifold_trim_by_plane(void *mem, ManifoldManifold *m,
                                          float normal_x, float normal_y,
                                          float normal_z, float offset);
 
+// 3D to 2D
+
+ManifoldCrossSection *manifold_slice(void *mem, ManifoldManifold *m,
+                                     float height);
+ManifoldCrossSection *manifold_project(void *mem, ManifoldManifold *m);
+
 // Convex Hulls
 
 ManifoldManifold *manifold_hull(void *mem, ManifoldManifold *m);

bindings/c/manifoldc.cpp

@@ -202,6 +202,17 @@ ManifoldManifold *manifold_trim_by_plane(void *mem, ManifoldManifold *m,
   return to_c(new (mem) Manifold(trimmed));
 }
 
+ManifoldCrossSection *manifold_slice(void *mem, ManifoldManifold *m,
+                                     float height) {
+  auto poly = from_c(m)->Slice(height);
+  return to_c(new (mem) CrossSection(poly));
+}
+
+ManifoldCrossSection *manifold_project(void *mem, ManifoldManifold *m) {
+  auto poly = from_c(m)->Project();
+  return to_c(new (mem) CrossSection(poly));
+}
+
 ManifoldManifold *manifold_hull(void *mem, ManifoldManifold *m) {
   auto hulled = from_c(m)->Hull();
   return to_c(new (mem) Manifold(hulled));

bindings/python/CMakeLists.txt

@@ -19,11 +19,22 @@ nanobind_add_module(
   manifold3d
   NB_STATIC STABLE_ABI LTO
   manifold3d.cpp)
-target_link_libraries(manifold3d PRIVATE manifold)
+target_link_libraries(manifold3d PRIVATE manifold sdf polygon)
 target_compile_options(manifold3d PRIVATE ${MANIFOLD_FLAGS} -DMODULE_NAME=manifold3d)
 target_compile_features(manifold3d PUBLIC cxx_std_17)
 set_target_properties(manifold3d PROPERTIES OUTPUT_NAME "manifold3d")
 
+message(Python_EXECUTABLE = ${Python_EXECUTABLE})
+add_custom_target(
+  autogen_docstrings
+  ${Python_EXECUTABLE} ${CMAKE_CURRENT_SOURCE_DIR}/gen_docs.py
+  WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}
+  BYPRODUCTS autogen_docstrings.inl
+)
+target_include_directories(manifold3d PRIVATE ${CMAKE_CURRENT_BINARY_DIR})
+add_dependencies(autogen_docstrings manifold sdf polygon)
+add_dependencies(manifold3d autogen_docstrings)
+
 if(SKBUILD)
 install(
   TARGETS manifold3d

bindings/python/README.md

@@ -0,0 +1,40 @@
+# Python Bindings
+
+## Autogenerated Doc-Strings
+
+Doc-strings for the python API wrapper are generated by gen_docs.py
+The script is run automatically during the build to keep python
+doc strings fully up-to-date with c++ sources.
+
+It scrapes documentation comments and c++ function signatures
+from the manifold c++ sources, in order to generate a c++ header file 
+exposing the comments as string variables named by function names.
+This allows python bindings to re-use the c++ comments directly.
+
+Some snake-casing of params is applied for python use case.
+
+---
+
+When modifying the Manifold C++ sources, you may need to update
+gen_docs.py. For example, top-level free functions are white-listed,
+so if you add a new one, you will need to add it in gen_docs.py.
+
+Similarly, the list of source files to parse is also white listed,
+so if you define functions in new files that need python wrappers,
+you will also need to up gen_docs.py.
+
+To verify that python docs are correct after changes, you can
+run the following commends from the manifold repo root:
+```
+pip install .
+python -c 'import manifold3d; help(manifold3d)'
+```
+
+Alternateively you could generate stubs with roughly the same info
+```
+pip install nanobind-stubgen
+pip install .
+nanobind-stubgen manifold3d
+```
+It will emit some warnings and write a file `manifold3d.pyi`
+which will show all the function signatures and docstrings.

bindings/python/examples/all_apis.py

@@ -0,0 +1,101 @@
+from manifold3d import *
+import numpy as np
+
+
+def all_root_level():
+    set_min_circular_angle(10)
+    set_min_circular_edge_length(1)
+    set_circular_segments(22)
+    n = get_circular_segments(1)
+    assert n == 22
+    poly = [[0, 0], [1, 0], [1, 1]]
+    tris = triangulate([poly])
+    tris = triangulate([np.array(poly)])
+
+
+def all_cross_section():
+    poly = [[0, 0], [1, 0], [1, 1]]
+    c = CrossSection([np.array(poly)])
+    c = CrossSection([poly])
+    c = CrossSection() + c
+    a = c.area()
+    c = CrossSection.batch_hull([c, c.translate((1, 0))])
+    b = c.bounds()
+    c = CrossSection.circle(1)
+    cs = c.decompose()
+    m = c.extrude(1)
+    c = c.hull()
+    c = CrossSection.hull_points(poly)
+    c = CrossSection.hull_points(np.array(poly))
+    e = c.is_empty()
+    c = c.mirror((0, 1))
+    n = c.num_contour()
+    n = c.num_vert()
+    c = c.offset(1, JoinType.Round)
+    m = c.revolve()
+    c = c.rotate(90)
+    c = c.scale((2, 2))
+    c = c.simplify()
+    c = CrossSection.square((1, 1))
+    p = c.to_polygons()
+    c = c.transform([[1, 0, 0], [0, 1, 0]])
+    c = c.translate((1, 1))
+    c = c.warp(lambda p: (p[0] + 1, p[1] / 2))
+    c = c.warp_batch(lambda ps: ps * [1, 0.5] + [1, 0])
+
+
+def all_manifold():
+    mesh = Manifold.sphere(1).to_mesh()
+    m = Manifold(mesh)
+    m = Manifold() + m
+    m = m.as_original()
+    m = Manifold.batch_hull([m, m.translate((0, 0, 1))])
+    b = m.bounding_box()
+    m = m.calculate_curvature(4, 5)
+    m = Manifold.compose([m, m.translate((5, 0, 0))])
+    m = Manifold.cube((1, 1, 1))
+    m = Manifold.cylinder(1, 1)
+    ms = m.decompose()
+    g = m.genus()
+    a = m.surface_area()
+    v = m.volume()
+    m = m.hull()
+    m = m.hull_points(mesh.vert_properties)
+    e = m.is_empty()
+    m = m.mirror((0, 0, 1))
+    n = m.num_edge()
+    n = m.num_prop()
+    n = m.num_prop_vert()
+    n = m.num_tri()
+    n = m.num_vert()
+    i = m.original_id()
+    p = m.precision()
+    c = m.project()
+    m = m.refine(2)
+    i = Manifold.reserve_ids(1)
+    m = m.scale((1, 2, 3))
+    m = m.set_properties(3, lambda pos, prop: pos)
+    c = m.slice(0.5)
+    m = Manifold.smooth(mesh, [0], [0.5])
+    m = Manifold.sphere(1)
+    m, n = m.split(m.translate((1, 0, 0)))
+    m, n = m.split_by_plane((0, 0, 1), 0)
+    e = m.status()
+    m = Manifold.tetrahedron()
+    mesh = m.to_mesh()
+    m = m.transform([[1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 1, 0]])
+    m = m.translate((0, 0, 0))
+    m = m.trim_by_plane((0, 0, 1), 0)
+    m = m.warp(lambda p: (p[0] + 1, p[1] / 2, p[2] * 2))
+    m = m.warp_batch(lambda ps: ps * [1, 0.5, 2] + [1, 0, 0])
+
+
+def run():
+    all_root_level()
+    all_cross_section()
+    all_manifold()
+    return Manifold()
+
+
+if __name__ == "__main__":
+    run()

bindings/python/examples/bricks.py

@@ -20,60 +20,82 @@ from manifold3d import Manifold
 
 INCHES = 25.4
 
-brick_depth = (3 + 5/8) * INCHES
-brick_height = (2 + 1/4) * INCHES
-brick_length = (7 + 5/8) * INCHES
+brick_depth = (3 + 5 / 8) * INCHES
+brick_height = (2 + 1 / 4) * INCHES
+brick_length = (7 + 5 / 8) * INCHES
 
-mortar_gap = (3/8) * INCHES
+mortar_gap = (3 / 8) * INCHES
 
 
 def brick():
-    return Manifold.cube(brick_length, brick_depth, brick_height)
+    return Manifold.cube([brick_length, brick_depth, brick_height])
 
 
 def halfbrick():
-    return Manifold.cube((brick_length - mortar_gap)/2, brick_depth, brick_height)
+    return Manifold.cube([(brick_length - mortar_gap) / 2, brick_depth, brick_height])
 
 
 def row(length):
-    bricks = [brick().translate((brick_length + mortar_gap)*x, 0, 0)
-              for x in range(length)]
-    return Manifold(bricks)
+    bricks = [
+        brick().translate([(brick_length + mortar_gap) * x, 0, 0])
+        for x in range(length)
+    ]
+    return sum(bricks, Manifold())
 
 
 def wall(length, height, alternate=0):
-    bricks = [row(length).translate(
-        ((z+alternate) % 2)*(brick_depth+mortar_gap), 0,
-        (brick_height+mortar_gap)*z) for z in range(height)]
-    return Manifold(bricks)
+    bricks = [
+        row(length).translate(
+            [
+                ((z + alternate) % 2) * (brick_depth + mortar_gap),
+                0,
+                (brick_height + mortar_gap) * z,
+            ]
+        )
+        for z in range(height)
+    ]
+    return sum(bricks, Manifold())
 
 
 def walls(length, width, height):
-    return Manifold([
-        wall(length, height),
-        wall(width, height, 1).rotate(0, 0, 90).translate(brick_depth, 0, 0),
-        wall(length, height, 1).translate(
-            0, (width)*(brick_length+mortar_gap), 0),
-        wall(width, height).rotate(0, 0, 90)
-        .translate((length+0.5)*(brick_length+mortar_gap) - mortar_gap, 0, 0),
-    ])
+    return sum(
+        [
+            wall(length, height),
+            wall(width, height, 1).rotate([0, 0, 90]).translate([brick_depth, 0, 0]),
+            wall(length, height, 1).translate(
+                [0, (width) * (brick_length + mortar_gap), 0]
+            ),
+            wall(width, height)
+            .rotate([0, 0, 90])
+            .translate(
+                [(length + 0.5) * (brick_length + mortar_gap) - mortar_gap, 0, 0]
+            ),
+        ],
+        Manifold(),
+    )
 
 
 def floor(length, width):
     results = [walls(length, width, 1)]
     if length > 1 and width > 1:
-        results.append(floor(length-1, width-1).translate(brick_depth +
-                       mortar_gap, brick_depth+mortar_gap, 0))
+        results.append(
+            floor(length - 1, width - 1).translate(
+                [brick_depth + mortar_gap, brick_depth + mortar_gap, 0]
+            )
+        )
     if length == 1 and width > 1:
-        results.append(row(width-1).rotate(0, 0, 90))
+        results.append(row(width - 1).rotate(0, 0, 90))
     if width == 1 and length > 1:
-        results.append(row(length-1).translate(2*(brick_depth +
-                       mortar_gap), brick_depth+mortar_gap, 0))
-    results.append(halfbrick().translate(
-        brick_depth+mortar_gap, brick_depth+mortar_gap, 0))
-    return Manifold(results)
+        results.append(
+            row(length - 1).translate(
+                [2 * (brick_depth + mortar_gap), brick_depth + mortar_gap, 0]
+            )
+        )
+    results.append(
+        halfbrick().translate([brick_depth + mortar_gap, brick_depth + mortar_gap, 0])
+    )
+    return sum(results, Manifold())
 
 
 def run(width=10, length=10, height=10):
     return walls(length, width, height) + floor(length, width)
-

bindings/python/examples/cube_with_dents.py

@@ -13,18 +13,21 @@
  See the License for the specific language governing permissions and
  limitations under the License.
  """
- 
+
 from manifold3d import Manifold
 from functools import reduce
 
 # https://gist.github.com/ochafik/2db96400e3c1f73558fcede990b8a355#file-cube-with-half-spheres-dents-scad
 
+
 def run(n=5, overlap=True):
-    a = Manifold.cube(n, n, 0.5).translate(-0.5, -0.5, -0.5)
+    a = Manifold.cube([n, n, 0.5]).translate([-0.5, -0.5, -0.5])
 
-    spheres = [Manifold.sphere(0.45 if overlap else 0.55, 50).translate(i, j, 0)
-               for i in range(n) for j in range(n)]
-    spheres = reduce(lambda a, b: a + b, spheres)
-
-    return a - spheres
+    spheres = [
+        Manifold.sphere(0.45 if overlap else 0.55, 50).translate([i, j, 0])
+        for i in range(n)
+        for j in range(n)
+    ]
+    # spheres = reduce(lambda a, b: a + b, spheres)
 
+    return a - sum(spheres, Manifold())

bindings/python/examples/extrude.py

@@ -10,26 +10,34 @@ def run():
     cross_section = CrossSection(polygons_points)
     polygons = cross_section.to_polygons()
     polygon = polygons[0]
-    if set(polygon) != set(polygon_points):
-        raise Exception(f"polygon={polygon} differs from polygon_points={polygon_points}")
+    if set([tuple(p) for p in polygon]) != set(polygon_points):
+        raise Exception(
+            f"polygon={polygon} differs from polygon_points={polygon_points}"
+        )
 
     # extrude a polygon to create a manifold
     extruded_polygon = cross_section.extrude(10.0)
     eps = 0.001
-    observed_volume = extruded_polygon.get_volume()
+    observed_volume = extruded_polygon.volume()
     expected_volume = 10.0
     if abs(observed_volume - expected_volume) > eps:
-        raise Exception(f"observed_volume={observed_volume} differs from expected_volume={expected_volume}")
-    observed_surface_area = extruded_polygon.get_surface_area()
+        raise Exception(
+            f"observed_volume={observed_volume} differs from expected_volume={expected_volume}"
+        )
+    observed_surface_area = extruded_polygon.surface_area()
     expected_surface_area = 42.0
     if abs(observed_surface_area - expected_surface_area) > eps:
-        raise Exception(f"observed_surface_area={observed_surface_area} differs from expected_surface_area={expected_surface_area}")
+        raise Exception(
+            f"observed_surface_area={observed_surface_area} differs from expected_surface_area={expected_surface_area}"
+        )
 
     # get bounding box from manifold
-    observed_bbox = extruded_polygon.bounding_box
+    observed_bbox = extruded_polygon.bounding_box()
     expected_bbox = (0.0, 0.0, 0.0, 1.0, 1.0, 10.0)
     if observed_bbox != expected_bbox:
-        raise Exception(f"observed_bbox={observed_bbox} differs from expected_bbox={expected_bbox}")
+        raise Exception(
+            f"observed_bbox={observed_bbox} differs from expected_bbox={expected_bbox}"
+        )
 
     return extruded_polygon
 

bindings/python/examples/gyroid_module.py

@@ -0,0 +1,59 @@
+"""
+ Copyright 2023 The Manifold Authors.
+
+ Licensed under the Apache License, Version 2.0 (the "License");
+ you may not use this file except in compliance with the License.
+ You may obtain a copy of the License at
+
+      https://www.apache.org/licenses/LICENSE-2.0
+
+ Unless required by applicable law or agreed to in writing, software
+ distributed under the License is distributed on an "AS IS" BASIS,
+ WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ See the License for the specific language governing permissions and
+ limitations under the License.
+ """
+
+import math
+import numpy as np
+from manifold3d import Mesh, Manifold
+
+
+def gyroid(x, y, z):
+    xi = x - math.pi / 4.0
+    yi = y - math.pi / 4.0
+    zi = z - math.pi / 4.0
+    return (
+        math.cos(xi) * math.sin(yi)
+        + math.cos(yi) * math.sin(zi)
+        + math.cos(zi) * math.sin(xi)
+    )
+
+
+def gyroid_levelset(level, period, size, n):
+    return Manifold(
+        Mesh.level_set(
+            gyroid,
+            [-period, -period, -period, period, period, period],
+            period / n,
+            level,
+        )
+    ).scale([size / period] * 3)
+
+
+def rhombic_dodecahedron(size):
+    box = Manifold.cube(size * math.sqrt(2.0) * np.array([1, 1, 2]), True)
+    result = box.rotate([90, 45, 0]) ^ box.rotate([90, 45, 90])
+    return result ^ box.rotate([0, 0, 45])
+
+
+def gyroid_module(size=20, n=15):
+    period = math.pi * 2.0
+    result = (
+        gyroid_levelset(-0.4, period, size, n) ^ rhombic_dodecahedron(size)
+    ) - gyroid_levelset(0.4, period, size, n)
+    return result.rotate([-45, 0, 90]).translate([0, 0, size / math.sqrt(2.0)])
+
+
+def run(size=20, n=15):
+    return gyroid_module(size, n)

bindings/python/examples/run_all.py

@@ -24,10 +24,9 @@ from time import time
 if __name__ == "__main__":
     current_file = pathlib.Path(__file__)
     current_dir = current_file.parent
-    files = [f.parts[-1][:-3]
-             for f in current_dir.glob('*.py') if f != current_file]
+    files = [f.parts[-1][:-3] for f in current_dir.glob("*.py") if f != current_file]
 
-    export_models = len(sys.argv) == 2 and sys.argv[-1] == '-e'
+    export_models = len(sys.argv) == 2 and sys.argv[-1] == "-e"
 
     for f in files:
         module = importlib.import_module(f)
@@ -41,9 +40,10 @@ if __name__ == "__main__":
             else:
                 vertices = mesh.vert_properties
                 colors = None
-            meshOut = trimesh.Trimesh(vertices=vertices, faces=mesh.tri_verts,
-                                      vertex_colors=colors)
-            trimesh.exchange.export.export_mesh(meshOut, f'{f}.glb', 'glb')
-            print(f'Exported model to {f}.glb')
+            meshOut = trimesh.Trimesh(
+                vertices=vertices, faces=mesh.tri_verts, vertex_colors=colors
+            )
+            trimesh.exchange.export.export_mesh(meshOut, f"{f}.glb", "glb")
+            print(f"Exported model to {f}.glb")
         t1 = time()
-        print(f'Took {(t1-t0)*1000:.1f}ms for {f}')
+        print(f"Took {(t1-t0)*1000:.1f}ms for {f}")

bindings/python/examples/scallop.py

@@ -0,0 +1,81 @@
+"""
+ Copyright 2023 The Manifold Authors.
+
+ Licensed under the Apache License, Version 2.0 (the "License");
+ you may not use this file except in compliance with the License.
+ You may obtain a copy of the License at
+
+      https://www.apache.org/licenses/LICENSE-2.0
+
+ Unless required by applicable law or agreed to in writing, software
+ distributed under the License is distributed on an "AS IS" BASIS,
+ WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ See the License for the specific language governing permissions and
+ limitations under the License.
+ """
+
+from manifold3d import Manifold, Mesh
+import numpy as np
+
+
+def run():
+    # A smoothed manifold demonstrating selective edge sharpening with
+    # smooth() and refine(), see more details at:
+    # https://elalish.blogspot.com/2022/03/smoothing-triangle-meshes.html
+
+    height = 10
+    radius = 30
+    offset = 20
+    wiggles = 12
+    sharpness = 0.8
+    n = 50
+
+    triangles = []
+    positions = [[-offset, 0, height], [-offset, 0, -height]]
+    sharpenedEdges = []
+
+    delta = np.pi / wiggles
+    for i in range(2 * wiggles):
+        theta = (i - wiggles) * delta
+        amp = 0.5 * height * max(np.cos(0.8 * theta), 0)
+
+        positions.append(
+            [
+                radius * np.cos(theta),
+                radius * np.sin(theta),
+                amp * (1 if i % 2 == 0 else -1),
+            ]
+        )
+
+        j = i + 1
+        if j == 2 * wiggles:
+            j = 0
+
+        smoothness = 1 - sharpness * np.cos((theta + delta / 2) / 2)
+        halfedge = 3 * len(triangles) + 1
+        sharpenedEdges.append((halfedge, smoothness))
+        triangles.append([0, 2 + i, 2 + j])
+
+        halfedge = 3 * len(triangles) + 1
+        sharpenedEdges.append((halfedge, smoothness))
+        triangles.append([1, 2 + j, 2 + i])
+
+    scallop = Mesh(
+        tri_verts=np.array(triangles, np.int32),
+        vert_properties=np.array(positions, np.float32),
+    )
+
+    def colorCurvature(_pos, oldProp):
+        a = max(0, min(1, oldProp[0] / 3 + 0.5))
+        b = a * a * (3 - 2 * a)
+        red = [1, 0, 0]
+        blue = [0, 0, 1]
+        return [(1 - b) * blue[i] + b * red[i] for i in range(3)]
+
+    edges, smoothing = zip(*sharpenedEdges)
+    return (
+        Manifold.smooth(scallop, edges, smoothing)
+        .refine(n)
+        .calculate_curvature(-1, 0)
+        .set_properties(3, colorCurvature)
+    )

bindings/python/examples/split_cube.py

@@ -0,0 +1,23 @@
+"""
+ Copyright 2022 The Manifold Authors.
+
+ Licensed under the Apache License, Version 2.0 (the "License");
+ you may not use this file except in compliance with the License.
+ You may obtain a copy of the License at
+
+      https://www.apache.org/licenses/LICENSE-2.0
+
+ Unless required by applicable law or agreed to in writing, software
+ distributed under the License is distributed on an "AS IS" BASIS,
+ WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ See the License for the specific language governing permissions and
+ limitations under the License.
+ """
+
+from manifold3d import Manifold
+
+
+def run():
+    a = Manifold.cube([1.0, 1.0, 1.0])
+    b = Manifold.cube([1.0, 1.0, 1.0]).rotate([45.0, 45.0, 45.0])
+    return a.split(b)[0]

bindings/python/examples/sponge.py

@@ -4,26 +4,27 @@ import numpy as np
 
 def fractal(holes, hole, w, position, depth, maxDepth):
     w /= 3
-    holes.append(hole.scale([w, w, 1.0]).translate(
-        [position[0], position[1], 0.0]))
+    holes.append(hole.scale([w, w, 1.0]).translate([position[0], position[1], 0.0]))
     if depth == maxDepth:
         return
-    offsets = np.array([
-        [-w, -w],
-        [-w, 0.0],
-        [-w, w],
-        [0.0, w],
-        [w, w],
-        [w, 0.0],
-        [w, -w],
-        [0.0, -w],
-    ])
+    offsets = np.array(
+        [
+            [-w, -w],
+            [-w, 0.0],
+            [-w, w],
+            [0.0, w],
+            [w, w],
+            [w, 0.0],
+            [w, -w],
+            [0.0, -w],
+        ]
+    )
     for offset in offsets:
         fractal(holes, hole, w, position + offset, depth + 1, maxDepth)
 
 
 def posColors(pos, _):
-    return [1 - p / 2 for p in pos] + [1.0]
+    return [-p + 0.5 for p in pos] + [1.0]
 
 
 def run(n=1):
@@ -37,4 +38,6 @@ def run(n=1):
     result -= hole.rotate([90, 0, 0])
     result -= hole.rotate([0, 90, 0])
 
-    return result.trim_by_plane([1, 1, 1], 0).set_properties(4, posColors).scale(100)
+    return (
+        result.trim_by_plane([1, 1, 1], 0).set_properties(4, posColors).scale([100] * 3)
+    )

bindings/python/examples/test_torus_knot.py

@@ -0,0 +1,98 @@
+from manifold3d import *
+import numpy as np
+import pytest
+
+# Creates a classic torus knot, defined as a string wrapping periodically
+# around the surface of an imaginary donut. If p and q have a common
+# factor then you will get multiple separate, interwoven knots. This is
+# an example of using the warp() method, thus avoiding any direct
+# handling of triangles.
+
+
+def run(warp_single=False):
+    # The number of times the thread passes through the donut hole.
+    p = 1
+    # The number of times the thread circles the donut.
+    q = 3
+    # Radius of the interior of the imaginary donut.
+    majorRadius = 25
+    # Radius of the small cross-section of the imaginary donut.
+    minorRadius = 10
+    # Radius of the small cross-section of the actual object.
+    threadRadius = 3.75
+    # Number of linear segments making up the threadRadius circle. Default is
+    # getCircularSegments(threadRadius).
+    circularSegments = -1
+    # Number of segments along the length of the knot. Default makes roughly
+    # square facets.
+    linearSegments = -1
+
+    # These default values recreate Matlab Knot by Emmett Lalish:
+    # https://www.thingiverse.com/thing:7080
+
+    kLoops = np.gcd(p, q)
+    pk = p / kLoops
+    qk = q / kLoops
+    n = (
+        circularSegments
+        if circularSegments > 2
+        else get_circular_segments(threadRadius)
+    )
+    m = linearSegments if linearSegments > 2 else n * qk * majorRadius / threadRadius
+
+    offset = 2
+    circle = CrossSection.circle(1, n).translate([offset, 0])
+
+    def ax_rotate(x, theta):
+        a, b = (x + 1) % 3, (x + 2) % 3
+        s, c = np.sin(theta), np.cos(theta)
+        m = np.zeros((len(theta), 4, 4), dtype=np.float32)
+        m[:, a, a], m[:, a, b] = c, s
+        m[:, b, a], m[:, b, b] = -s, c
+        m[:, x, x], m[:, 3, 3] = 1, 1
+        return m
+
+    def func(pts):
+        npts = pts.shape[0]
+        x, y, z = pts[:, 0], pts[:, 1], pts[:, 2]
+        psi = qk * np.arctan2(x, y)
+        theta = psi * pk / qk
+        x1 = np.sqrt(x * x + y * y)
+        phi = np.arctan2(x1 - offset, z)
+
+        v = np.zeros((npts, 4), dtype=np.float32)
+        v[:, 0] = threadRadius * np.cos(phi)
+        v[:, 2] = threadRadius * np.sin(phi)
+        v[:, 3] = 1
+        r = majorRadius + minorRadius * np.cos(theta)
+
+        m1 = ax_rotate(0, -np.arctan2(pk * minorRadius, qk * r))
+        m1[:, 3, 0] = minorRadius
+        m2 = ax_rotate(1, theta)
+        m2[:, 3, 0] = majorRadius
+        m3 = ax_rotate(2, psi)
+
+        v = v[:, None, :] @ m1 @ m2 @ m3
+        pts[:] = v[:, 0, :3]
+
+    def func_single(v):
+        pts = np.array([v])
+        func(pts)
+        return pts[0]
+
+    if warp_single:
+        return circle.revolve(int(m)).warp(func_single)
+    else:
+        return circle.revolve(int(m)).warp_batch(func)
+
+
+@pytest.mark.parametrize("warp_single", [True, False])
+def test_warp(warp_single):
+    m = run(warp_single=warp_single)
+    assert m.volume() == pytest.approx(20785.76)
+    assert m.surface_area() == pytest.approx(11176.8)
+    assert m.genus() == 1
+
+
+if __name__ == "__main__":
+    run()

bindings/python/examples/union_failure.py

@@ -1,7 +1,8 @@
 from manifold3d import Manifold
 
+
 def run():
     # for some reason this causes collider error
     obj = Manifold.cube()
-    obj += Manifold([Manifold.cube()]).rotate(0, 0, 45+180)
+    obj += Manifold.cube().rotate([0, 0, 45 + 180])
     return obj

bindings/python/gen_docs.py

@@ -0,0 +1,109 @@
+from os.path import dirname
+from hashlib import md5
+import re
+
+base = dirname(dirname(dirname(__file__)))
+
+def snake_case(name):
+    return re.sub("([a-z0-9])([A-Z])", r"\1_\2", name).lower()
+
+
+def python_param_modifier(comment):
+    # p = f":{snake_case(m[0][1:])}:"
+    comment = re.sub(r"@(param \w+)", lambda m: f':{snake_case(m[1])}:', comment)
+    # python API renames `MeshGL` to `Mesh`
+    comment = re.sub('mesh_gl', 'mesh', comment)
+    comment = re.sub('MeshGL', 'Mesh', comment)
+    return comment
+
+
+def method_key(name):
+    name = re.sub("\+", "_plus", name)
+    name = re.sub("\-", "_minus", name)
+    name = re.sub("\^", "_xor", name)
+    name = re.sub("\=", "_eq", name)
+    name = re.sub("\:", "_", name)
+    name = re.sub("\~", "destroy_", name)
+    return name
+
+
+parens_re = re.compile(r"[^(]+\(([^(]*(\(.*\))*[^(]*\))", flags=re.DOTALL)
+args_re = re.compile(
+    r"^[^,^\(^\)]*(\(.*\))*[^,^\(^\)]*[\s\&\*]([0-9\w]+)\s*[,\)]", flags=re.DOTALL
+)
+
+
+def parse_args(s):
+    par = parens_re.match(s)
+    if not par:
+        return None
+    out = []
+    arg_str = par[1]
+    while m := re.search(args_re, arg_str):
+        out += [snake_case(m[2])]
+        arg_str = arg_str[m.span()[1] :]
+    return out
+
+
+def collect(fname, matcher, param_modifier=python_param_modifier):
+    comment = ""
+    with open(fname) as f:
+        for line in f:
+            line = line.lstrip()
+            if line.startswith("/**"):
+                comment = ""
+            elif line.startswith("*/"):
+                pass
+            elif line.startswith("*") and comment is not None:
+                comment += line[1:].lstrip()
+            elif comment and (m := matcher(line)):
+                while (args := parse_args(line)) is None:
+                    line += next(f)
+                    if len(line) > 500:
+                        break
+
+                method = method_key(snake_case(m[1]))
+                # comment = re.sub(param_re, param_modifier, comment)
+                comment = param_modifier(comment)
+                method = "__".join([method, *args])
+                assert method not in comments
+                comments[method] = comment
+                comment = ""
+
+
+comments = {}
+
+method_re = re.compile(r"(\w+::[\w\-\+\^\=\:]+)\(")
+function_re = re.compile(r"([\w\-\+\^\=\:]+)\(")
+
+
+# we don't handle inline functions in classes properly
+# so instead just white-list functions we want
+def select_functions(s):
+    m = function_re.search(s)
+    if m and "Triangulate" in m[0]:
+        return m
+    if m and "Circular" in m[0]:
+        return m
+    return None
+
+
+collect(f"{base}/src/manifold/src/manifold.cpp", lambda s: method_re.search(s))
+collect(f"{base}/src/manifold/src/constructors.cpp", lambda s: method_re.search(s))
+collect(
+    f"{base}/src/cross_section/src/cross_section.cpp", lambda s: method_re.search(s)
+)
+collect(f"{base}/src/polygon/src/polygon.cpp", select_functions)
+collect(f"{base}/src/utilities/include/public.h", select_functions)
+
+comments = dict(sorted(comments.items()))
+
+gen_h = f"autogen_docstrings.inl"
+with open(gen_h, "w") as f:
+    f.write("#pragma once\n\n")
+    f.write("// --- AUTO GENERATED ---\n")
+    f.write("// gen_docs.py is run by cmake build\n\n")
+    f.write("namespace manifold_docstrings {\n")
+    for key, doc in comments.items():
+        f.write(f'const char* {key} = R"___({doc.strip()})___";\n')
+    f.write("} // namespace manifold_docs")

bindings/wasm/CMakeLists.txt

@@ -32,13 +32,21 @@ add_custom_target(js_deps ALL
         DEPENDS manifoldjs
                ${CMAKE_CURRENT_SOURCE_DIR}/manifold*.d.ts)
 
+# copy WASM build back here for publishing to npm
 add_custom_command(
         TARGET js_deps POST_BUILD
         COMMAND ${CMAKE_COMMAND} -E copy
                ${CMAKE_CURRENT_BINARY_DIR}/manifold.*
-               ${CMAKE_CURRENT_SOURCE_DIR}/manifold*.d.ts
+               ${CMAKE_CURRENT_SOURCE_DIR})
+
+# copy WASM build and TS declarations for Vite to package into ManifoldCAD.org
+add_custom_command(
+        TARGET js_deps POST_BUILD
+        COMMAND ${CMAKE_COMMAND} -E copy
+               ${CMAKE_CURRENT_SOURCE_DIR}/manifold*
                ${CMAKE_CURRENT_SOURCE_DIR}/examples/built/)
 
+# copy TS declarations to public so they can be accessed by our editor
 add_custom_command(
         TARGET js_deps POST_BUILD
         COMMAND ${CMAKE_COMMAND} -E copy

bindings/wasm/README.md

@@ -24,7 +24,7 @@ Since Manifold is a WASM module, it does not automatically garbage-collect like
 
 ## Examples
 
-Please see our usage [examples](https://github.com/elalish/manifold/tree/master/bindings/wasm/examples) to see how to interface this library with `three.js`, `<model-viewer>`, and `glTF`. Of particular note are the included libraries for lossless roundtrip of manifold meshes through glTF files, via a new extension: [EXT_manifold](https://github.com/KhronosGroup/glTF/pull/2286). 
+Please see our usage [examples](https://github.com/elalish/manifold/tree/master/bindings/wasm/examples) to see how to interface this library with `three.js`, `<model-viewer>`, and `glTF`. Of particular note are the included libraries for lossless roundtrip of manifold meshes through glTF files, via a new extension: [EXT_mesh_manifold](https://github.com/KhronosGroup/glTF/pull/2286). 
 
 ## About the author
 

bindings/wasm/bindings.cpp

@@ -138,6 +138,8 @@ EMSCRIPTEN_BINDINGS(whatever) {
       .function("_Split", &man_js::Split)
       .function("_SplitByPlane", &man_js::SplitByPlane)
       .function("_TrimByPlane", &Manifold::TrimByPlane)
+      .function("slice", &Manifold::Slice)
+      .function("project", &Manifold::Project)
       .function("hull", select_overload<Manifold() const>(&Manifold::Hull))
       .function("_GetMeshJS", &js::GetMeshJS)
       .function("refine", &Manifold::Refine)

bindings/wasm/examples/editor.css

@@ -65,7 +65,7 @@ button:hover {
   filter: brightness(0.9);
 }
 
-input {
+.name {
   border: none;
   height: 36px;
   background-color: transparent;
@@ -77,7 +77,7 @@ input {
   -webkit-font-smoothing: antialiased;
 }
 
-input:focus {
+.name:focus {
   outline: none;
 }
 
@@ -196,6 +196,14 @@ model-viewer {
   background-image: url(/icons/trash.png);
 }
 
+.play {
+  background-image: url(/icons/play.png);
+}
+
+.pause {
+  background-image: url(/icons/pause.png);
+}
+
 .icon:hover {
   background-color: rgba(0, 0, 0, 0.2);
 }
@@ -221,6 +229,51 @@ model-viewer {
   width: 125px;
 }
 
+#animation {
+  display: none;
+  align-items: center;
+  position: absolute;
+  width: 100%;
+  bottom: 0px;
+  margin: 2px;
+}
+
+#play {
+  position: relative;
+  background-color: #ccc;
+  margin: 8px;
+}
+
+.slider {
+  position: relative;
+  -webkit-appearance: none;
+  appearance: none;
+  width: 100%;
+  height: 10px;
+  border-radius: 5px;
+  margin-right: 16px;
+  background: #ffffff55;
+  outline: none;
+}
+
+.slider::-webkit-slider-thumb {
+  -webkit-appearance: none;
+  appearance: none;
+  width: 25px;
+  height: 25px;
+  border-radius: 50%;
+  background: #ccc;
+  cursor: pointer;
+}
+
+.slider::-moz-range-thumb {
+  width: 25px;
+  height: 25px;
+  border-radius: 50%;
+  background: #ccc;
+  cursor: pointer;
+}
+
 @media only screen and (max-width: 600px) {
   .container {
     flex-direction: column;
@@ -256,4 +309,4 @@ model-viewer {
     width: calc(50% - 16px);
     margin: 8px;
   }
-}
+}

bindings/wasm/examples/editor.js

@@ -90,9 +90,9 @@ function switchTo(scriptName) {
     switching = true;
     currentFileElement.textContent = scriptName;
     setScript('currentName', scriptName);
-    const code =
-        exampleFunctions.get(scriptName) ?? getScript(scriptName) ?? '';
     isExample = exampleFunctions.get(scriptName) != null;
+    const code = isExample ? exampleFunctions.get(scriptName).substring(1) :
+                             getScript(scriptName) ?? '';
     editor.setValue(code);
   }
 }
@@ -110,6 +110,8 @@ function createDropdownItem(name) {
 
   button.onclick = function() {
     saveCurrent();
+    window.location.hash = `#${label.textContent}`;
+    window.location.search = '';
     switchTo(label.textContent);
   };
   // Stop text input spaces from triggering the button
@@ -146,6 +148,7 @@ function addEdit(button) {
     const code = getScript(oldName);
     const form = document.createElement('form');
     const inputElement = document.createElement('input');
+    inputElement.classList.add('name');
     inputElement.value = oldName;
     label.textContent = '';
     button.appendChild(form);
@@ -195,6 +198,7 @@ function addEdit(button) {
       removeScript(label.textContent);
       if (currentFileElement.textContent == label.textContent) {
         switchTo('Intro');
+        window.location.hash = '#Intro';
       }
       const container = button.parentElement;
       container.parentElement.removeChild(container);
@@ -203,16 +207,16 @@ function addEdit(button) {
 }
 
 const newButton = document.querySelector('#new');
-function newItem(code) {
-  const name = uniqueName('New Script');
+function newItem(code, scriptName = undefined) {
+  const name = uniqueName(scriptName ?? 'New Script');
   setScript(name, code);
   const nextButton = createDropdownItem(name);
   newButton.insertAdjacentElement('afterend', nextButton.parentElement);
   addEdit(nextButton);
-  nextButton.click();
+  return {button: nextButton, name};
 };
 newButton.onclick = function() {
-  newItem('');
+  newItem('').button.click();
 };
 
 const runButton = document.querySelector('#compile');
@@ -298,7 +302,20 @@ require(['vs/editor/editor.main'], async function() {
       addEdit(button);
     }
   }
-  switchTo(currentName);
+  const params = new URLSearchParams(window.location.search);
+  if (window.location.hash.length > 0) {
+    const name = unescape(window.location.hash.substring(1));
+    switchTo(name);
+  } else if (params.get('script') != null) {
+    console.log(`Fetching ${params.get('script')}`);
+    autoExecute = false;
+    const response = await fetch(params.get('script'));
+    const code = await response.text();
+    switchTo(newItem(code, params.get('name')).name);
+  } else {
+    switchTo(currentName);
+    window.location.hash = `#${currentName}`;
+  }
 
   if (manifoldInitialized) {
     initializeRun();
@@ -312,7 +329,7 @@ require(['vs/editor/editor.main'], async function() {
     }
     if (isExample) {
       const cursor = editor.getPosition();
-      newItem(editor.getValue());
+      newItem(editor.getValue()).button.click();
       editor.setPosition(cursor);
     }
   });
@@ -322,6 +339,51 @@ require(['vs/editor/editor.main'], async function() {
   };
 });
 
+// Animation ------------------------------------------------------------
+const mv = document.querySelector('model-viewer');
+const animationContainer = document.querySelector('#animation');
+const playButton = document.querySelector('#play');
+const scrubber = document.querySelector('#scrubber');
+let paused = false;
+
+mv.addEventListener('load', () => {
+  const hasAnimation = mv.availableAnimations.length > 0;
+  animationContainer.style.display = hasAnimation ? 'flex' : 'none';
+  if (hasAnimation) {
+    play();
+  }
+});
+
+function play() {
+  mv.play();
+  playButton.classList.remove('play');
+  playButton.classList.add('pause');
+  paused = false;
+  scrubber.classList.add('hide');
+}
+
+function pause() {
+  mv.pause();
+  playButton.classList.remove('pause');
+  playButton.classList.add('play');
+  paused = true;
+  scrubber.max = mv.duration;
+  scrubber.value = mv.currentTime;
+  scrubber.classList.remove('hide');
+}
+
+playButton.onclick = function() {
+  if (paused) {
+    play();
+  } else {
+    pause();
+  }
+};
+
+scrubber.oninput = function() {
+  mv.currentTime = scrubber.value;
+};
+
 // Execution ------------------------------------------------------------
 const consoleElement = document.querySelector('#console');
 const oldLog = console.log;
@@ -357,7 +419,6 @@ function finishRun() {
       `Took ${(Math.round((t1 - t0) / 10) / 100).toLocaleString()} seconds`);
 }
 
-const mv = document.querySelector('model-viewer');
 const output = {
   glbURL: null,
   threeMFURL: null

bindings/wasm/examples/gltf-io.js

@@ -1,285 +0,0 @@
-// Copyright 2023 The Manifold Authors.
-//
-// Licensed under the Apache License, Version 2.0 (the "License");
-// you may not use this file except in compliance with the License.
-// You may obtain a copy of the License at
-//
-//      http://www.apache.org/licenses/LICENSE-2.0
-//
-// Unless required by applicable law or agreed to in writing, software
-// distributed under the License is distributed on an "AS IS" BASIS,
-// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
-// See the License for the specific language governing permissions and
-// limitations under the License.
-
-import {Accessor} from '@gltf-transform/core';
-
-import {EXTManifold} from './manifold-gltf.js';
-
-export const attributeDefs = {
-  'POSITION': {type: Accessor.Type.VEC3, components: 3},
-  'NORMAL': {type: Accessor.Type.VEC3, components: 3},
-  'TANGENT': {type: Accessor.Type.VEC4, components: 4},
-  'TEXCOORD_0': {type: Accessor.Type.VEC2, components: 2},
-  'TEXCOORD_1': {type: Accessor.Type.VEC2, components: 2},
-  'COLOR_0': {type: Accessor.Type.VEC3, components: 3},
-  'JOINTS_0': {type: Accessor.Type.VEC4, components: 4},
-  'WEIGHTS_0': {type: Accessor.Type.VEC4, components: 4},
-  'SKIP_1': {type: null, components: 1},
-  'SKIP_2': {type: null, components: 2},
-  'SKIP_3': {type: null, components: 3},
-  'SKIP_4': {type: null, components: 4},
-};
-
-export function setupIO(io) {
-  return io.registerExtensions([EXTManifold]);
-}
-
-function readPrimitive(primitive, numProp, attributes) {
-  const position = primitive.getAttribute('POSITION');
-  const numVert = position.getCount();
-  const vertProperties = [];
-  let offset = 0;
-  for (const attribute of attributes) {
-    const size = attributeDefs[attribute].components;
-    if (attribute === 'SKIP') {
-      offset += size;
-      continue;
-    }
-    const accessor = primitive.getAttribute(attribute);
-    if (accessor) {
-      const array = accessor.getArray();
-      for (let i = 0; i < numVert; ++i) {
-        for (let j = 0; j < size; ++j) {
-          vertProperties[numProp * i + offset + j] = array[i * size + j];
-        }
-      }
-    } else {
-      for (let i = 0; i < numVert; ++i) {
-        for (let j = 0; j < size; ++j) {
-          vertProperties[numProp * i + offset + j] = 0;
-        }
-      }
-    }
-    offset += size;
-  }
-  return vertProperties;
-}
-
-export function readMesh(mesh, attributes, materials) {
-  const primitives = mesh.listPrimitives();
-  if (primitives.length === 0) {
-    return {};
-  }
-
-  if (attributes.length === 0) {
-    const attributeSet = new Set();
-    for (const primitive of primitives) {
-      const semantics = primitive.listSemantics();
-      for (const semantic of semantics) {
-        attributeSet.add(semantic);
-      }
-    }
-    for (const semantic in attributeDefs) {
-      if (attributeSet.has(semantic)) {
-        attributes.push(semantic);
-        attributeSet.delete(semantic);
-      }
-    }
-    for (const semantic of attributeSet.keys()) {
-      attributes.push(semantic);
-    }
-  }
-
-  if (attributes.length < 1 || attributes[0] !== 'POSITION')
-    throw new Error('First attribute must be "POSITION".');
-
-  const numProp = attributes.map((def) => attributeDefs[def].components)
-                      .reduce((a, b) => a + b);
-
-  const manifoldPrimitive = mesh.getExtension('EXT_manifold');
-
-  let vertPropArray = [];
-  let triVertArray = [];
-  const runIndexArray = [0];
-  const mergeFromVert = [];
-  const mergeToVert = [];
-  if (manifoldPrimitive != null) {
-    vertPropArray = readPrimitive(primitives[0], numProp, attributes);
-    for (const primitive of primitives) {
-      triVertArray = [...triVertArray, ...primitive.getIndices().getArray()];
-      runIndexArray.push(triVertArray.length);
-      materials.push(primitive.getMaterial());
-    }
-    const mergeTriVert = manifoldPrimitive.getMergeIndices()?.getArray() ?? [];
-    const mergeTo = manifoldPrimitive.getMergeValues()?.getArray() ?? [];
-    const vert2merge = new Map();
-    for (const [i, idx] of mergeTriVert.entries()) {
-      vert2merge.set(triVertArray[idx], mergeTo[i]);
-    }
-    for (const [from, to] of vert2merge.entries()) {
-      mergeFromVert.push(from);
-      mergeToVert.push(to);
-    }
-  } else {
-    for (const primitive of primitives) {
-      const numVert = vertPropArray.length / numProp;
-      vertPropArray =
-          [...vertPropArray, ...readPrimitive(primitive, numProp, attributes)];
-      triVertArray = [
-        ...triVertArray,
-        ...primitive.getIndices().getArray().map((i) => i + numVert)
-      ];
-      runIndexArray.push(triVertArray.length);
-      materials.push(primitive.getMaterial());
-    }
-  }
-  const vertProperties = new Float32Array(vertPropArray);
-  const triVerts = new Uint32Array(triVertArray);
-  const runIndex = new Uint32Array(runIndexArray);
-
-  return {
-    numProp,
-    triVerts,
-    vertProperties,
-    runIndex,
-    mergeFromVert,
-    mergeToVert
-  };
-}
-
-export function writeMesh(doc, manifoldMesh, attributes, materials) {
-  if (doc.getRoot().listBuffers().length === 0) {
-    doc.createBuffer();
-  }
-  const buffer = doc.getRoot().listBuffers()[0];
-  const manifoldExtension = doc.createExtension(EXTManifold);
-
-  const attributeUnion = [];
-  for (const matAttributes of attributes) {
-    matAttributes.forEach((attribute, i) => {
-      if (i >= attributeUnion.length) {
-        attributeUnion.push(attribute);
-      } else {
-        const size = attributeDefs[attribute].components;
-        const unionSize = attributeDefs[attributeUnion[i]].components;
-        if (size != unionSize) {
-          throw new Error(
-              'Attribute sizes do not correspond: ', attribute, ' and ',
-              attributeUnion[i]);
-        }
-        if (attributeDefs[attributeUnion[i]].type == null) {
-          attributeUnion[i] = attribute;
-        }
-      }
-    });
-  }
-  if (attributeUnion.length < 1 || attributeUnion[0] !== 'POSITION')
-    throw new Error('First attribute must be "POSITION".');
-
-  const mesh = doc.createMesh();
-  const numPrimitive = manifoldMesh.runIndex.length - 1;
-  for (let run = 0; run < numPrimitive; ++run) {
-    const indices = doc.createAccessor('indices')
-                        .setBuffer(buffer)
-                        .setType(Accessor.Type.SCALAR)
-                        .setArray(new Uint32Array());
-    const primitive = doc.createPrimitive().setIndices(indices);
-    const material = materials[run];
-    if (material) {
-      primitive.setMaterial(material);
-    }
-    mesh.addPrimitive(primitive);
-  }
-
-  const numVert = manifoldMesh.numVert;
-  const numProp = manifoldMesh.numProp;
-  let offset = 0;
-  attributeUnion.forEach((attribute, aIdx) => {
-    if (attribute === 'SKIP') {
-      ++offset;
-      return;
-    }
-    const def = attributeDefs[attribute];
-    if (def == null)
-      throw new Error(attribute + ' is not a recognized attribute.');
-
-    const n = def.components;
-    if (offset + n > numProp) throw new Error('Too many attribute channels.');
-
-    const array = new Float32Array(n * numVert);
-    for (let v = 0; v < numVert; ++v) {
-      for (let i = 0; i < n; ++i) {
-        array[n * v + i] =
-            manifoldMesh.vertProperties[numProp * v + offset + i];
-      }
-    }
-
-    const accessor =
-        doc.createAccessor().setBuffer(buffer).setType(def.type).setArray(
-            array);
-
-    mesh.listPrimitives().forEach((primitive, pIdx) => {
-      if (attributes[pIdx].length > aIdx &&
-          attributeDefs[attributes[pIdx][aIdx]].type != null) {
-        primitive.setAttribute(attribute, accessor);
-      }
-    });
-    offset += n;
-  });
-
-  const manifoldPrimitive = manifoldExtension.createManifoldPrimitive();
-  mesh.setExtension('EXT_manifold', manifoldPrimitive);
-
-  const indices = doc.createAccessor('indices')
-                      .setBuffer(buffer)
-                      .setType(Accessor.Type.SCALAR)
-                      .setArray(manifoldMesh.triVerts);
-  manifoldPrimitive.setIndices(indices);
-  manifoldPrimitive.setRunIndex(manifoldMesh.runIndex);
-
-  const vert2merge = [...Array(manifoldMesh.numVert).keys()];
-  for (const [i, from] of manifoldMesh.mergeFromVert.entries()) {
-    vert2merge[from] = manifoldMesh.mergeToVert[i];
-  }
-  const ind = [];
-  const val = [];
-  for (const [i, vert] of manifoldMesh.triVerts.entries()) {
-    const newVert = vert2merge[vert];
-    if (vert !== newVert) {
-      ind.push(i);
-      val.push(newVert);
-    }
-  }
-  const indicesAccessor = doc.createAccessor()
-                              .setBuffer(buffer)
-                              .setType(Accessor.Type.SCALAR)
-                              .setArray(new Uint32Array(ind));
-  const valuesAccessor = doc.createAccessor()
-                             .setBuffer(buffer)
-                             .setType(Accessor.Type.SCALAR)
-                             .setArray(new Uint32Array(val));
-  manifoldPrimitive.setMerge(indicesAccessor, valuesAccessor);
-
-  return mesh;
-}
-
-export function disposeMesh(mesh) {
-  if (!mesh) return;
-  const primitives = mesh.listPrimitives();
-  for (const primitive of primitives) {
-    primitive.getIndices()?.dispose();
-    for (const accessor of primitive.listAttributes()) {
-      accessor.dispose();
-    }
-  }
-
-  mesh.dispose();
-}
-
-export async function loadTexture(texture, uri) {
-  const response = await fetch(uri);
-  const blob = await response.blob();
-  texture.setMimeType(blob.type);
-  texture.setImage(new Uint8Array(await blob.arrayBuffer()));
-}

bindings/wasm/examples/gltf-io.ts

@@ -0,0 +1,424 @@
+// Copyright 2023 The Manifold Authors.
+//
+// Licensed under the Apache License, Version 2.0 (the "License");
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+//
+//      http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+
+import {Accessor, Document, Material, Mesh, Primitive, Texture, WebIO} from '@gltf-transform/core';
+
+import {EXTManifold, ManifoldPrimitive} from './manifold-gltf';
+import {Mesh as ManifoldMesh, MeshOptions} from './public/manifold-encapsulated-types';
+
+export const attributeDefs = {
+  'POSITION': {type: Accessor.Type.VEC3, components: 3},
+  'NORMAL': {type: Accessor.Type.VEC3, components: 3},
+  'TANGENT': {type: Accessor.Type.VEC4, components: 4},
+  'TEXCOORD_0': {type: Accessor.Type.VEC2, components: 2},
+  'TEXCOORD_1': {type: Accessor.Type.VEC2, components: 2},
+  'COLOR_0': {type: Accessor.Type.VEC3, components: 3},
+  'JOINTS_0': {type: Accessor.Type.VEC4, components: 4},
+  'WEIGHTS_0': {type: Accessor.Type.VEC4, components: 4},
+  'SKIP_1': {type: null, components: 1},
+  'SKIP_2': {type: null, components: 2},
+  'SKIP_3': {type: null, components: 3},
+  'SKIP_4': {type: null, components: 4},
+};
+
+export type Attribute = keyof(typeof attributeDefs);
+
+export interface Properties {
+  material: Material;
+  attributes: Attribute[];
+}
+
+/**
+ * Call this first to register the manifold extension so that readMesh and
+ * writeMesh will work.
+ */
+export function setupIO(io: WebIO) {
+  return io.registerExtensions([EXTManifold]);
+}
+
+/**
+ * Read an input mesh into Manifold-compatible data structures, whether it
+ * contains the EXT_mesh_manifold extension or not.
+ *
+ * @param mesh The Mesh to read.
+ * @param attributes An array of attributes representing the order of desired
+ *     properties returned in the vertProperties array of the output mesh. If
+ *     omitted, this will be populated with the union of all attributes defined
+ *     in the primitives of the input mesh. If present, the first entry must be
+ *     'POSITION', and any attributes in the primitives that are not included in
+ *     this list will be ignored, while those in the list but not defined in a
+ *     primitive will be populated with zeros.
+ * @returns The returned mesh is suitable for initializing a Manifold or Mesh of
+ *     the Manifold library if desired. See Manifold documentation if you prefer
+ *     to use these GL arrays in a different library. The runProperties array
+ *     gives the Material and attributes list associated with each triangle run,
+ *     which in turn corresponds to a primitive of the input mesh. These
+ *     attributes are the intersection of the attributes present on the
+ *     primitive and those requested in the attributes input.
+ */
+export function readMesh(mesh: Mesh, attributes: Attribute[] = []):
+    {mesh: MeshOptions, runProperties: Properties[]}|null {
+  const primitives = mesh.listPrimitives();
+  if (primitives.length === 0) {
+    return null;
+  }
+
+  if (attributes.length === 0) {
+    const attributeSet = new Set<Attribute>();
+    for (const primitive of primitives) {
+      const semantics = primitive.listSemantics() as Attribute[];
+      for (const semantic of semantics) {
+        attributeSet.add(semantic);
+      }
+    }
+    let semantic: Attribute;
+    for (semantic in attributeDefs) {
+      if (attributeSet.has(semantic)) {
+        attributes.push(semantic);
+        attributeSet.delete(semantic);
+      }
+    }
+    for (const semantic of attributeSet.keys()) {
+      attributes.push(semantic);
+    }
+  }
+
+  if (attributes.length < 1 || attributes[0] !== 'POSITION')
+    throw new Error('First attribute must be "POSITION".');
+
+  let numProp = 0;
+  const attributeOffsets = attributes.map((numProp = 0, def => {
+    const last = numProp;
+    numProp += attributeDefs[def].components;
+    return last;
+  }));
+
+  const manifoldPrimitive =
+      mesh.getExtension('EXT_mesh_manifold') as ManifoldPrimitive;
+
+  let vertPropArray: number[] = [];
+  let triVertArray: number[] = [];
+  const runIndexArray = [0];
+  const mergeFromVertArray = [];
+  const mergeToVertArray = [];
+  const runProperties: Properties[] = [];
+  if (manifoldPrimitive != null) {
+    const numVert = primitives[0].getAttribute('POSITION')!.getCount();
+    const foundAttribute = attributes.map((a) => attributeDefs[a].type == null);
+    vertPropArray = new Array<number>(numProp * numVert);
+
+    for (const primitive of primitives) {
+      const indices = primitive.getIndices();
+      if (!indices) {
+        console.log('Skipping non-indexed primitive ', primitive.getName());
+        continue;
+      }
+
+      const attributesIn = primitive.listSemantics() as Attribute[];
+
+      attributes.forEach((attributeOut, idx) => {
+        if (foundAttribute[idx]) {
+          return;
+        }
+        for (const attributeIn of attributesIn) {
+          if (attributeIn === attributeOut) {
+            foundAttribute[idx] = true;
+            const accessor = primitive.getAttribute(attributeIn)!;
+            writeProperties(
+                vertPropArray, accessor, numProp, attributeOffsets[idx]);
+          }
+        }
+      });
+
+      triVertArray = [...triVertArray, ...indices.getArray()!];
+      runIndexArray.push(triVertArray.length);
+      runProperties.push({
+        material: primitive.getMaterial()!,
+        attributes: attributesIn.filter(b => attributes.some(a => a == b))
+      });
+    }
+    const mergeTriVert = manifoldPrimitive.getMergeIndices()?.getArray() ?? [];
+    const mergeTo = manifoldPrimitive.getMergeValues()?.getArray() ?? [];
+    const vert2merge = new Map();
+    for (const [i, idx] of mergeTriVert.entries()) {
+      vert2merge.set(triVertArray[idx], mergeTo[i]);
+    }
+    for (const [from, to] of vert2merge.entries()) {
+      mergeFromVertArray.push(from);
+      mergeToVertArray.push(to);
+    }
+  } else {
+    for (const primitive of primitives) {
+      const indices = primitive.getIndices();
+      if (!indices) {
+        console.log('Skipping non-indexed primitive ', primitive.getName());
+        continue;
+      }
+      const numVert = vertPropArray.length / numProp;
+      vertPropArray =
+          [...vertPropArray, ...readPrimitive(primitive, numProp, attributes)];
+      triVertArray =
+          [...triVertArray, ...indices.getArray()!.map((i) => i + numVert)];
+      runIndexArray.push(triVertArray.length);
+
+      const attributesIn = primitive.listSemantics() as Attribute[];
+      runProperties.push({
+        material: primitive.getMaterial()!,
+        attributes: attributesIn.filter(b => attributes.some(a => a == b))
+      });
+    }
+  }
+  const vertProperties = new Float32Array(vertPropArray);
+  const triVerts = new Uint32Array(triVertArray);
+  const runIndex = new Uint32Array(runIndexArray);
+  const mergeFromVert = new Uint32Array(mergeFromVertArray);
+  const mergeToVert = new Uint32Array(mergeToVertArray);
+
+  const meshOut =
+      {numProp, triVerts, vertProperties, runIndex, mergeFromVert, mergeToVert};
+
+  return {mesh: meshOut, runProperties};
+}
+
+/**
+ * Write a Manifold Mesh into a glTF Mesh object, using the EXT_mesh_manifold
+ * extension to allow for lossless roundtrip of the manifold mesh through the
+ * glTF file.
+ *
+ * @param doc The glTF Document to which this Mesh will be added.
+ * @param manifoldMesh The Manifold Mesh to convert to glTF.
+ * @param id2properties A map from originalID to Properties that include the
+ *     glTF Material and the set of attributes to output. All triangle runs with
+ *     the same originalID will be combined into a single output primitive. Any
+ *     originalIDs not found in the map will have the glTF default material and
+ *     no attributes beyond 'POSITION'. Each attributes array must correspond to
+ *     the manifoldMesh vertProperties, thus the first attribute must always be
+ *     'POSITION'. Any properties that should not be output for a given
+ *     primitive must use the 'SKIP_*' attributes.
+ * @returns The glTF Mesh to add to the Document.
+ */
+export function writeMesh(
+    doc: Document, manifoldMesh: ManifoldMesh,
+    id2properties: Map<number, Properties>): Mesh {
+  if (doc.getRoot().listBuffers().length === 0) {
+    doc.createBuffer();
+  }
+  const buffer = doc.getRoot().listBuffers()[0];
+  const manifoldExtension = doc.createExtension(EXTManifold);
+
+  const mesh = doc.createMesh();
+  const runIndex = [];
+  const attributeUnion: Attribute[] = [];
+  const primitive2attributes = new Map<Primitive, Attribute[]>();
+  const numRun = manifoldMesh.runIndex.length - 1;
+  let lastID = -1;
+  for (let run = 0; run < numRun; ++run) {
+    const id = manifoldMesh.runOriginalID[run];
+    if (id == lastID) {
+      continue;
+    }
+    lastID = id;
+    runIndex.push(manifoldMesh.runIndex[run]);
+
+    const indices = doc.createAccessor('primitive indices of ID ' + id)
+                        .setBuffer(buffer)
+                        .setType(Accessor.Type.SCALAR)
+                        .setArray(new Uint32Array(1));
+    const primitive = doc.createPrimitive().setIndices(indices);
+
+    const properties = id2properties.get(id);
+    if (properties) {
+      const {material, attributes} = properties;
+      if (attributes.length < 1 || attributes[0] !== 'POSITION')
+        throw new Error('First attribute must be "POSITION".');
+
+      primitive.setMaterial(material);
+      primitive2attributes.set(primitive, attributes);
+
+      properties.attributes.forEach((attribute, i) => {
+        if (i >= attributeUnion.length) {
+          attributeUnion.push(attribute);
+        } else {
+          const size = attributeDefs[attribute].components;
+          const unionSize = attributeDefs[attributeUnion[i]].components;
+          if (size != unionSize) {
+            throw new Error(
+                'Attribute sizes do not correspond: ' + attribute + ' and ' +
+                attributeUnion[i]);
+          }
+          if (attributeDefs[attributeUnion[i]].type == null) {
+            attributeUnion[i] = attribute;
+          }
+        }
+      });
+    } else {
+      primitive2attributes.set(primitive, ['POSITION']);
+    }
+
+    mesh.addPrimitive(primitive);
+  }
+  runIndex.push(manifoldMesh.runIndex[numRun]);
+
+  const numVert = manifoldMesh.numVert;
+  const numProp = manifoldMesh.numProp;
+  let offset = 0;
+  attributeUnion.forEach((attribute, aIdx) => {
+    const def = attributeDefs[attribute];
+    if (def == null)
+      throw new Error(attribute + ' is not a recognized attribute.');
+
+    if (def.type == null) {
+      ++offset;
+      return;
+    }
+
+    const n = def.components;
+    if (offset + n > numProp) throw new Error('Too many attribute channels.');
+
+    const array = new Float32Array(n * numVert);
+    for (let v = 0; v < numVert; ++v) {
+      for (let i = 0; i < n; ++i) {
+        let x = manifoldMesh.vertProperties[numProp * v + offset + i];
+        if (attribute == 'COLOR_0') {
+          x = Math.max(0, Math.min(1, x));
+        }
+        array[n * v + i] = x;
+      }
+    }
+
+    const accessor = doc.createAccessor(attribute)
+                         .setBuffer(buffer)
+                         .setType(def.type)
+                         .setArray(array);
+
+    for (const primitive of mesh.listPrimitives()) {
+      const attributes = primitive2attributes.get(primitive)!;
+      if (attributes.length > aIdx &&
+          attributeDefs[attributes[aIdx]].type != null) {
+        primitive.setAttribute(attribute, accessor);
+      }
+    }
+    offset += n;
+  });
+
+  const manifoldPrimitive = manifoldExtension.createManifoldPrimitive();
+  mesh.setExtension('EXT_mesh_manifold', manifoldPrimitive);
+
+  const indices = doc.createAccessor('manifold indices')
+                      .setBuffer(buffer)
+                      .setType(Accessor.Type.SCALAR)
+                      .setArray(manifoldMesh.triVerts);
+  manifoldPrimitive.setIndices(indices);
+  manifoldPrimitive.setRunIndex(runIndex);
+
+  const vert2merge = [...Array(manifoldMesh.numVert).keys()];
+  const ind = [];
+  const val = [];
+  if (manifoldMesh.mergeFromVert && manifoldMesh.mergeToVert) {
+    for (const [i, from] of manifoldMesh.mergeFromVert.entries()) {
+      vert2merge[from] = manifoldMesh.mergeToVert[i];
+    }
+
+    for (const [i, vert] of manifoldMesh.triVerts.entries()) {
+      const newVert = vert2merge[vert];
+      if (vert !== newVert) {
+        ind.push(i);
+        val.push(newVert);
+      }
+    }
+  }
+  if (ind.length > 0) {
+    const indicesAccessor = doc.createAccessor('merge from')
+                                .setBuffer(buffer)
+                                .setType(Accessor.Type.SCALAR)
+                                .setArray(new Uint32Array(ind));
+    const valuesAccessor = doc.createAccessor('merge to')
+                               .setBuffer(buffer)
+                               .setType(Accessor.Type.SCALAR)
+                               .setArray(new Uint32Array(val));
+    manifoldPrimitive.setMerge(indicesAccessor, valuesAccessor);
+  }
+
+  return mesh;
+}
+
+/**
+ * Helper function to dispose of a Mesh, useful when replacing an existing Mesh
+ * with one from writeMesh.
+ */
+export function disposeMesh(mesh: Mesh) {
+  if (!mesh) return;
+  const primitives = mesh.listPrimitives();
+  for (const primitive of primitives) {
+    primitive.getIndices()?.dispose();
+    for (const accessor of primitive.listAttributes()) {
+      accessor.dispose();
+    }
+  }
+
+  const manifoldPrimitive =
+      mesh.getExtension('EXT_mesh_manifold') as ManifoldPrimitive;
+  if (manifoldPrimitive) {
+    manifoldPrimitive.getIndices()?.dispose();
+    manifoldPrimitive.getMergeIndices()?.dispose();
+    manifoldPrimitive.getMergeValues()?.dispose();
+  }
+
+  mesh.dispose();
+}
+
+/**
+ * Helper function to download an image and apply it to the given texture.
+ *
+ * @param texture The texture to update
+ * @param uri The location of the image to download
+ */
+export async function loadTexture(texture: Texture, uri: string) {
+  const response = await fetch(uri);
+  const blob = await response.blob();
+  texture.setMimeType(blob.type);
+  texture.setImage(new Uint8Array(await blob.arrayBuffer()));
+}
+
+function writeProperties(
+    vertProperties: number[], accessor: Accessor, numProp: number,
+    offset: number) {
+  const array = accessor.getArray()!;
+  const size = accessor.getElementSize();
+  const numVert = accessor.getCount();
+  for (let i = 0; i < numVert; ++i) {
+    for (let j = 0; j < size; ++j) {
+      vertProperties[numProp * i + offset + j] = array[i * size + j];
+    }
+  }
+}
+
+function readPrimitive(
+    primitive: Primitive, numProp: number, attributes: Attribute[]) {
+  const vertProperties: number[] = [];
+  let offset = 0;
+  for (const attribute of attributes) {
+    const size = attributeDefs[attribute].components;
+    if (attributeDefs[attribute].type == null) {
+      offset += size;
+      continue;
+    }
+    const accessor = primitive.getAttribute(attribute);
+    if (accessor) {
+      writeProperties(vertProperties, accessor, numProp, offset);
+    }
+    offset += size;
+  }
+  return vertProperties;
+}

bindings/wasm/examples/index.html

@@ -40,7 +40,7 @@
   <div class="page">
     <div class="container header" style="background-color: #A288B3;">
       <p style="flex: 2;"><a target="_blank" href="https://github.com/elalish/manifold">Manifold</a>CAD </p>
-      <div class="margin" id="file-margin" >
+      <div class="margin" id="file-margin">
         <button type="button" id="file" class="blue item"><b class="uparrow"></b>&nbsp;&nbsp;<span
             id="current">Intro</span></button>
         <div class="dropdown" id="fileDropdown">
@@ -70,6 +70,10 @@
           <div class="center" slot="poster">
             <p id="poster" style="text-align: center; line-height: 36px;"> Loading... </p>
           </div>
+          <div id="animation">
+            <button type="button" id="play" class="icon pause"></button>
+            <input type="range" min="0" max="1" step="0.01" value="0" class="slider" id="scrubber">
+          </div>
         </model-viewer>
         <div id="console"></div>
       </div>
@@ -84,6 +88,6 @@
   integrity="sha512-iQEIc0rsSDujsfjtD+lfyJ1W23Bh/lbgriubKDAym6VlEIDRj9rrbSIyJRyshOrl8s0yRcQ0+gyrZfSLyjJGWQ=="
   crossorigin="anonymous" referrerpolicy="no-referrer">
 <script type="module" src="editor.js"></script>
-<script type="module" src="https://ajax.googleapis.com/ajax/libs/model-viewer/3.1.1/model-viewer.min.js"></script>
+<script type="module" src="https://ajax.googleapis.com/ajax/libs/model-viewer/3.3.0/model-viewer.min.js"></script>
 
-</html>
+</html>

bindings/wasm/examples/make-manifold.html

@@ -46,7 +46,7 @@
 <body>
   <p>Load a glTF/GLB model and our <a href="https://github.com/elalish/manifold">Manifold</a> library will attempt to
     merge it into a set of manifold objects. If the model is water-tight, you can download the new GLB which will have
-    our <a href="https://github.com/KhronosGroup/glTF/pull/2286">EXT_manifold</a> extension, thus preserving the
+    our <a href="https://github.com/KhronosGroup/glTF/pull/2286">EXT_mesh_manifold</a> extension, thus preserving the
     manifold data without losing any mesh properties.</p>
   <p> If the View Manifold GLB checkbox is enabled, then some meshes in the model have open edges and don't represent a
     solid object. Check this box to see only the manifold parts, which will also enable them to be downloaded. If the
@@ -59,13 +59,13 @@
     <span id="poster" slot="poster">Drop a GLB here</span>
   </model-viewer>
 </body>
-<script type="module" src="https://ajax.googleapis.com/ajax/libs/model-viewer/2.0.1/model-viewer.min.js"></script>
+<script type="module" src="https://ajax.googleapis.com/ajax/libs/model-viewer/3.3.0/model-viewer.min.js"></script>
 <script type="module">
   import { WebIO } from '@gltf-transform/core';
   import { prune } from '@gltf-transform/functions';
   import { KHRONOS_EXTENSIONS } from '@gltf-transform/extensions';
   import { SimpleDropzone } from 'simple-dropzone';
-  import { readMesh, writeMesh, setupIO } from './gltf-io.js';
+  import { readMesh, writeMesh, setupIO, disposeMesh } from './gltf-io';
   import Module from './built/manifold.js';
 
   const io = setupIO(new WebIO());
@@ -147,31 +147,30 @@
     const docIn = await io.read(url);
     const nodes = docIn.getRoot().listNodes();
     for (const node of nodes) {
-      const attributes = [];
       const mesh = node.getMesh();
       if (!mesh) {
         continue;
       }
 
-      const materials = [];
-      const tmpMesh = readMesh(mesh, attributes, materials);
-      tmpMesh.runOriginalID = [];
-      const firstID = Manifold.reserveIDs(materials.length);
-      for (let i = 0; i < materials.length; ++i) {
-        tmpMesh.runOriginalID.push(firstID + i);
-      }
-      const manifoldMesh = new Mesh(tmpMesh);
-      mesh.dispose();
+      const tmp = readMesh(mesh);
 
-      const matAttributes = [];
-      materials.forEach(() => matAttributes.push(attributes));
+      tmp.mesh.runOriginalID = [];
+      const id2properties = new Map();
+      const numID = tmp.runProperties.length;
+      const firstID = Manifold.reserveIDs(numID);
+      for (let i = 0; i < numID; ++i) {
+        tmp.mesh.runOriginalID.push(firstID + i);
+        id2properties.set(firstID + i, tmp.runProperties[i]);
+      }
+      const manifoldMesh = new Mesh(tmp.mesh);
+      disposeMesh(mesh);
 
       manifoldMesh.merge();
 
       try {
         // Test manifoldness - will throw if not.
         const manifold = Manifold(manifoldMesh);
-        node.setMesh(writeMesh(docIn, manifold.getMesh(), matAttributes, materials));
+        node.setMesh(writeMesh(docIn, manifold.getMesh(), id2properties));
         manifold.delete();
         anyManifold = true;
       } catch (e) {

bindings/wasm/examples/manifold-gltf.ts

@@ -12,15 +12,32 @@
 // See the License for the specific language governing permissions and
 // limitations under the License.
 
-import {Extension, ExtensionProperty, PropertyType} from '@gltf-transform/core';
+import {Accessor, Extension, ExtensionProperty, GLTF, IProperty, PropertyType, ReaderContext, WriterContext} from '@gltf-transform/core';
 
-const NAME = 'EXT_manifold';
+const NAME = 'EXT_mesh_manifold';
+const MERGE = 'MERGE';
 
-export class ManifoldPrimitive extends ExtensionProperty {
+interface IManifoldPrimitive extends IProperty {
+  mergeIndices: Accessor;
+  mergeValues: Accessor;
+  indices: Accessor;
+  runIndex: number[]|Uint32Array;
+}
+
+interface ManifoldDef {
+  manifoldPrimitive: GLTF.IMeshPrimitive;
+  mergeIndices: number;
+  mergeValues: number;
+}
+
+export class ManifoldPrimitive extends ExtensionProperty<IManifoldPrimitive> {
   static EXTENSION_NAME = NAME;
+  public declare extensionName: typeof NAME;
+  public declare propertyType: 'ManifoldPrimitive';
+  public declare parentTypes: [PropertyType.MESH];
 
   init() {
-    this.EXTENSION_NAME = NAME;
+    ManifoldPrimitive.EXTENSION_NAME = NAME;
     this.propertyType = 'ManifoldPrimitive';
     this.parentTypes = [PropertyType.MESH];
   }
@@ -39,7 +56,7 @@ export class ManifoldPrimitive extends ExtensionProperty {
     return this.getRef('mergeValues');
   }
 
-  setMerge(indicesAccessor, valuesAccessor) {
+  setMerge(indicesAccessor: Accessor, valuesAccessor: Accessor) {
     if (indicesAccessor.getCount() !== valuesAccessor.getCount())
       throw new Error('merge vectors must be the same length.');
     this.setRef('mergeIndices', indicesAccessor);
@@ -50,28 +67,29 @@ export class ManifoldPrimitive extends ExtensionProperty {
     return this.get('runIndex');
   }
 
-  setRunIndex(runIndex) {
+  setRunIndex(runIndex: number[]|Uint32Array) {
     return this.set('runIndex', runIndex);
   }
 
-  setIndices(indices) {
+  setIndices(indices: Accessor) {
     return this.setRef('indices', indices);
   }
 
   getIndices() {
-    return this.getRef('indices');
+    return this.getRef('indices')!;
   }
 }
 
 export class EXTManifold extends Extension {
   extensionName = NAME;
+  prewriteTypes = [PropertyType.ACCESSOR];
   static EXTENSION_NAME = NAME;
 
   createManifoldPrimitive() {
     return new ManifoldPrimitive(this.document.getGraph());
   }
 
-  read(context) {
+  read(context: ReaderContext) {
     const {json} = context.jsonDoc;
     const meshDefs = json.meshes || [];
 
@@ -82,20 +100,27 @@ export class EXTManifold extends Extension {
       const manifoldPrimitive = this.createManifoldPrimitive();
       mesh.setExtension(NAME, manifoldPrimitive);
 
-      const manifoldDef = meshDef.extensions[NAME];
+      const manifoldDef = meshDef.extensions[NAME] as ManifoldDef;
 
       if (manifoldDef.manifoldPrimitive) {
         let count = 0;
         const runIndex = [];
         runIndex.push(count);
         for (const primitive of mesh.listPrimitives()) {
-          count += primitive.getIndices().getCount();
+          const indices = primitive.getIndices();
+          if (!indices) {
+            console.log('Skipping non-indexed primitive ', primitive.getName());
+            continue;
+          }
+          count += indices.getCount();
           runIndex.push(count);
         }
         manifoldPrimitive.setRunIndex(runIndex);
+        manifoldPrimitive.setIndices(
+            context.accessors[manifoldDef.manifoldPrimitive.indices!]);
       }
 
-      if (manifoldDef.mergeIndices && manifoldDef.mergeValues) {
+      if (manifoldDef.mergeIndices != null && manifoldDef.mergeValues != null) {
         manifoldPrimitive.setMerge(
             context.accessors[manifoldDef.mergeIndices],
             context.accessors[manifoldDef.mergeValues]);
@@ -105,15 +130,34 @@ export class EXTManifold extends Extension {
     return this;
   }
 
-  write(context) {
+  prewrite(context: WriterContext) {
+    this.document.getRoot().listMeshes().forEach((mesh) => {
+      const manifoldPrimitive = mesh.getExtension(NAME) as ManifoldPrimitive;
+      if (!manifoldPrimitive) return;
+
+      const indices = manifoldPrimitive.getIndices();
+      context.addAccessorToUsageGroup(
+          indices, WriterContext.BufferViewUsage.ELEMENT_ARRAY_BUFFER);
+
+      const mergeFrom = manifoldPrimitive.getMergeIndices();
+      const mergeTo = manifoldPrimitive.getMergeValues();
+      if (!mergeFrom || !mergeTo) return;
+
+      context.addAccessorToUsageGroup(mergeFrom, MERGE);
+      context.addAccessorToUsageGroup(mergeTo, MERGE);
+    });
+    return this;
+  }
+
+  write(context: WriterContext) {
     const {json} = context.jsonDoc;
 
     this.document.getRoot().listMeshes().forEach((mesh) => {
-      const manifoldPrimitive = mesh.getExtension(NAME);
+      const manifoldPrimitive = mesh.getExtension(NAME) as ManifoldPrimitive;
       if (!manifoldPrimitive) return;
 
-      const meshIndex = context.meshIndexMap.get(mesh);
-      const meshDef = json.meshes[meshIndex];
+      const meshIndex = context.meshIndexMap.get(mesh)!;
+      const meshDef = json.meshes![meshIndex];
 
       const runIndex = manifoldPrimitive.getRunIndex();
       const numPrimitive = runIndex.length - 1;
@@ -124,20 +168,20 @@ export class EXTManifold extends Extension {
       }
 
       const mergeIndicesIndex =
-          context.accessorIndexMap.get(manifoldPrimitive.getMergeIndices());
+          context.accessorIndexMap.get(manifoldPrimitive.getMergeIndices()!)!;
       const mergeValuesIndex =
-          context.accessorIndexMap.get(manifoldPrimitive.getMergeValues());
-      const mergeIndices = json.accessors[mergeIndicesIndex];
-      const mergeValues = json.accessors[mergeValuesIndex];
+          context.accessorIndexMap.get(manifoldPrimitive.getMergeValues()!)!;
+      const mergeIndices = json.accessors![mergeIndicesIndex];
+      const mergeValues = json.accessors![mergeValuesIndex];
 
       const existingPrimitive = meshDef.primitives[0];
       const primitive = {
-        indices: context.accessorIndexMap.get(manifoldPrimitive.getIndices()),
+        indices: context.accessorIndexMap.get(manifoldPrimitive.getIndices())!,
         mode: existingPrimitive.mode,
         attributes: {'POSITION': existingPrimitive.attributes['POSITION']}
       };
 
-      const indices = json.accessors[primitive.indices];
+      const indices = json.accessors![primitive.indices];
       if (!indices) {
         return;
       }
@@ -146,21 +190,21 @@ export class EXTManifold extends Extension {
         indices.sparse = {
           count: mergeIndices.count,
           indices: {
-            bufferView: mergeIndices.bufferView,
+            bufferView: mergeIndices.bufferView!,
             byteOffset: mergeIndices.byteOffset,
             componentType: mergeIndices.componentType
           },
           values: {
-            bufferView: mergeValues.bufferView,
+            bufferView: mergeValues.bufferView!,
             byteOffset: mergeValues.byteOffset,
           }
         };
       }
 
       for (let i = 0; i < numPrimitive; ++i) {
-        const accessor = json.accessors[meshDef.primitives[i].indices];
+        const accessor = json.accessors![meshDef.primitives[i].indices!];
         accessor.bufferView = indices.bufferView;
-        accessor.byteOffset = indices.byteOffset + 4 * runIndex[i];
+        accessor.byteOffset = indices.byteOffset! + 4 * runIndex[i];
         accessor.count = runIndex[i + 1] - runIndex[i];
       }
 

bindings/wasm/examples/model-viewer.html

@@ -48,11 +48,11 @@
   <button id="download">Download manifold GLB</button>
   <model-viewer camera-controls shadow-intensity="1" alt="Output of mesh Boolean operation"></model-viewer>
 </body>
-<script type="module" src="https://ajax.googleapis.com/ajax/libs/model-viewer/2.0.1/model-viewer.min.js"></script>
+<script type="module" src="https://ajax.googleapis.com/ajax/libs/model-viewer/3.3.0/model-viewer.min.js"></script>
 <script type="module">
   import { Document, WebIO } from '@gltf-transform/core';
   import { clearNodeTransform, flatten, prune } from '@gltf-transform/functions';
-  import { loadTexture, readMesh, writeMesh, disposeMesh, setupIO } from './gltf-io.js';
+  import { loadTexture, readMesh, writeMesh, disposeMesh, setupIO } from './gltf-io';
   import Module from './built/manifold.js';
 
   const io = setupIO(new WebIO());
@@ -63,10 +63,9 @@
 
   const { Manifold, Mesh } = wasm;
 
-  const attributes = ['POSITION', 'TEXCOORD_0'];
-  const materialMap = new Map();
+  const id2properties = new Map();
 
-  async function readGLB(url, attributes, materialMap) {
+  async function readGLB(url) {
     const manifolds = [];
     const docIn = await io.read(url);
     await docIn.transform(flatten());
@@ -75,30 +74,33 @@
     for (const node of nodes) {
       clearNodeTransform(node);
       const materials = [];
-      const manifoldMesh = readMesh(node.getMesh(), attributes, materials);
-      manifoldMesh.runOriginalID = [];
-      const firstID = Manifold.reserveIDs(materials.length);
-      for (let i = 0; i < materials.length; ++i) {
-        manifoldMesh.runOriginalID.push(firstID + i);
+      const tmp = readMesh(node.getMesh());
+      tmp.mesh.runOriginalID = [];
+      const numID = tmp.runProperties.length;
+      const firstID = Manifold.reserveIDs(numID);
+      for (let i = 0; i < numID; ++i) {
+        tmp.mesh.runOriginalID.push(firstID + i);
         ids.push(firstID + i);
+        id2properties.set(firstID + i, tmp.runProperties[i]);
       }
-      manifolds.push(Manifold(new Mesh(manifoldMesh)));
+      manifolds.push(Manifold(new Mesh(tmp.mesh)));
     }
     // pull in materials, TODO: replace with transfer() when available
     const startIdx = doc.getRoot().listMaterials().length;
     doc.merge(docIn);
     doc.getRoot().listScenes().forEach((s) => s.dispose());
     doc.getRoot().listBuffers().forEach((s) => s.dispose());
+    doc.getRoot().listAccessors().forEach((s) => s.dispose());
     for (const [i, id] of ids.entries()) {
       const material = doc.getRoot().listMaterials()[startIdx + i];
-      materialMap.set(id, material);
+      id2properties.get(id).material = material;
     }
 
     return Manifold.union(manifolds);
   }
 
-  const space = await readGLB('/models/space.glb', attributes, materialMap);
-  const moon = await readGLB('/models/moon.glb', attributes, materialMap);
+  const space = await readGLB('/models/space.glb');
+  const moon = await readGLB('/models/moon.glb');
 
   const node = doc.createNode();
   doc.createScene().addChild(node);
@@ -129,13 +131,7 @@
     const manifoldMesh = manifold.rotate([-90, 0, 0]).getMesh();
 
     disposeMesh(node.getMesh());
-    const materials = [];
-    const matAttributes = [];
-    for (const id of manifoldMesh.runOriginalID) {
-      materials.push(materialMap.get(id));
-      matAttributes.push(attributes);
-    }
-    const mesh = writeMesh(doc, manifoldMesh, matAttributes, materials);
+    const mesh = writeMesh(doc, manifoldMesh, id2properties);
     node.setMesh(mesh);
     await doc.transform(prune());
 

bindings/wasm/examples/package.json

@@ -8,25 +8,25 @@
     "build": "tsc --noEmit && vite build",
     "build:watch": "vite build --watch",
     "preview": "vite preview",
-    "test": "vitest",
+    "test": "vitest run",
     "test:ui": "vitest --ui",
-    "test:run": "vitest run"
+    "test:dev": "vitest"
   },
   "dependencies": {
-    "@gltf-transform/core": "^3.2.1",
-    "@gltf-transform/extensions": "^3.2.1",
-    "@gltf-transform/functions": "^3.2.1",
-    "@jscadui/3mf-export": "^0.4.2",
+    "@gltf-transform/core": "^3.8.0",
+    "@gltf-transform/extensions": "^3.8.0",
+    "@gltf-transform/functions": "^3.8.0",
+    "@jscadui/3mf-export": "^0.5.0",
     "fflate": "^0.8.0",
     "gl-matrix": "^3.4.3",
     "simple-dropzone": "0.8.3",
-    "three": "0.151.2"
+    "three": "0.158.0"
   },
   "devDependencies": {
     "@vitest/ui": "^0.31.1",
     "@vitest/web-worker": "^0.31.1",
-    "typescript": "5.0.4",
-    "vite": "^4.3.2",
+    "typescript": "5.2.2",
+    "vite": "^4.5.0",
     "vitest": "^0.31.1"
   }
-}
+}

bindings/wasm/examples/public/editor.d.ts

@@ -25,8 +25,11 @@ declare class GLTFNode {
   clone(parent?: GLTFNode): GLTFNode;
 }
 
+type Attribute = 'POSITION'|'NORMAL'|'TANGENT'|'TEXCOORD_0'|'TEXCOORD_1'|
+    'COLOR_0'|'JOINTS_0'|'WEIGHTS_0'|'SKIP_1'|'SKIP_2'|'SKIP_3'|'SKIP_4';
+
 declare class GLTFMaterial {
-  attributes?: string[];
+  attributes?: Attribute[];
   roughness?: number;
   metallic?: number;
   baseColorFactor?: [number, number, number];
@@ -35,15 +38,48 @@ declare class GLTFMaterial {
   name?: string;
 }
 
+declare const globalDefaults: {
+  roughness: number,
+  metallic: number,
+  baseColorFactor: [number, number, number],
+  alpha: number,
+  unlit: boolean,
+  animationLength: number,
+  animationMode: 'loop'|'ping-pong';
+}
+
 /**
- * Set material properties on the input manifold. They will be carried along
- * through operations.
+ * Returns a shallow copy of the input manifold with the given material
+ * properties applied. They will be carried along through operations.
  *
- * @param manifold The object to add properties to - returned for chaining.
+ * @param manifold The input object.
  * @param material A set of material properties to apply to this manifold.
  */
-declare function setMaterial(
-    manifold: Manifold, material: GLTFMaterial): Manifold;
+declare function setMaterial(manifold: Manifold, material: GLTFMaterial):
+    Manifold;
+
+/**
+ * Apply a morphing animation to the input manifold. Specify the start
+ * function which will be applied to the vertex positions of the first frame and
+ * linearly interpolated across the length of the overall animation. This
+ * animation will only be shown if this manifold is used directly on a GLTFNode.
+ *
+ * @param manifold The object to add morphing animation to.
+ * @param func A warping function to apply to the first animation frame.
+ */
+declare function setMorphStart(
+    manifold: Manifold, func: (v: Vec3) => void): void;
+
+/**
+ * Apply a morphing animation to the input manifold. Specify the end
+ * function which will be applied to the vertex positions of the last frame and
+ * linearly interpolated across the length of the overall animation. This
+ * animation will only be shown if this manifold is used directly on a GLTFNode.
+ *
+ * @param manifold The object to add morphing animation to.
+ * @param func A warping function to apply to the last animation frame.
+ */
+declare function setMorphEnd(manifold: Manifold, func: (v: Vec3) => void): void;
 
 /**
  * Wrap any shape object with this method to display it and any copies in

bindings/wasm/examples/public/examples.js

@@ -50,9 +50,9 @@ export const examples = {
 
     TetrahedronPuzzle: function() {
       // A tetrahedron cut into two identical halves that can screw together as
-      // a puzzle. This only outputs one of the halves. This demonstrates how
-      // redundant points along a polygon can be used to make twisted extrusions
-      // smoother. Based on the screw puzzle by George Hart:
+      // a puzzle. This demonstrates how redundant points along a polygon can be
+      // used to make twisted extrusions smoother. It also showcases animated
+      // assemblies. Based on the screw puzzle by George Hart:
       // https://www.thingiverse.com/thing:186372
 
       const edgeLength = 50;  // Length of each edge of the overall tetrahedron.
@@ -61,19 +61,21 @@ export const examples = {
 
       const scale = edgeLength / (2 * Math.sqrt(2));
 
-      const tet = Manifold.tetrahedron().scale(scale);
+      const tet = Manifold.tetrahedron().intersect(
+          Manifold.tetrahedron().rotate([0, 0, 90]).scale(2.5));
 
       const box = [];
-      box.push([1, -1], [1, 1]);
+      box.push([2, -2], [2, 2]);
       for (let i = 0; i <= nDivisions; ++i) {
-        box.push([gap / (4 * scale), 1 - i * 2 / nDivisions]);
+        box.push([gap / (2 * scale), 2 - i * 4 / nDivisions]);
       }
 
       const cyan = [0, 1, 1];
       const magenta = [1, 0, 1];
       const fade = (color, pos) => {
         for (let i = 0; i < 3; ++i) {
-          color[i] = cyan[i] * pos[2] + magenta[i] * (1 - pos[2]);
+          const x = pos[2] / 2;
+          color[i] = cyan[i] * x + magenta[i] * (1 - x);
         }
       };
 
@@ -82,25 +84,43 @@ export const examples = {
       // channels as colors, and sets the factor to white, since our default is
       // yellow.
       const screw = setMaterial(
-          Manifold.extrude(box, 1, nDivisions, 270).setProperties(3, fade),
+          Manifold.extrude(box, 2, nDivisions, 270).setProperties(3, fade),
           {baseColorFactor: [1, 1, 1], attributes: ['COLOR_0']});
 
-      const result = tet.intersect(
-          screw.rotate([0, 0, -45]).translate([0, 0, -0.5]).scale(2 * scale));
+      const result =
+          tet.intersect(screw.rotate([0, 0, -45]).translate([0, 0, -1]))
+              .scale(scale);
 
       // Assigned materials are only applied to a GLTFNode. Note that material
       // definitions cascade, applying recursively to all child surfaces, but
-      // overridden by any materials defined lower down. Our default material:
-      // {
-      //   roughness = 0.2,
-      //   metallic = 1,
-      //   baseColorFactor = [1, 1, 0],
-      //   alpha = 1,
-      //   unlit = false,
-      //   name = ''
-      // }
-      const node = new GLTFNode();
-      node.manifold = result;
+      // overridden by any materials defined lower down. The default material
+      // properties, as well as animation parameters can be set via
+      // globalDefaults.
+
+      const layFlat = new GLTFNode();
+      layFlat.rotation =
+          [45, -Math.atan(1 / Math.sqrt(2)) * 180 / Math.PI, 120];
+      layFlat.translation = [0, 0, scale * Math.sqrt(3) / 3];
+
+      const fixed = new GLTFNode(layFlat);
+      fixed.manifold = result;
+      fixed.rotation = [0, 0, 180];
+
+      // For 3MF export, only top-level objects are independently arrangeable.
+      const layFlat2 = layFlat.clone();
+
+      const moving = new GLTFNode(layFlat2);
+      moving.manifold = result;
+      // Use functions to create animation, which runs from t=0 to t=1.
+      moving.translation = (t) => {
+        const a = 1 - t;
+        const x = a > 0.5 ? scale * 2 * (0.5 - a) : 0;
+        return [x, x, -2 * x + scale * 4 * (0.5 - Math.abs(0.5 - a))];
+      };
+      moving.rotation = (t) => [0, 0, 270 * 2 * (1 - t)];
+
+      globalDefaults.animationLength = 10;  // seconds
+      globalDefaults.animationMode = 'ping-pong';
       return result;
     },
 
@@ -139,6 +159,7 @@ export const examples = {
       // Smooth, complex manifolds can be created using the warp() function.
       // This example recreates the Exploitable Heart by Emmett Lalish:
       // https://www.thingiverse.com/thing:6190
+      // It also demonstrates the use of setMorph to animate a warping function.
 
       const func = (v) => {
         const x2 = v[0] * v[0];
@@ -174,9 +195,20 @@ export const examples = {
       };
 
       const ball = Manifold.sphere(1, 200);
+
       const heart = ball.warp(func);
       const box = heart.boundingBox();
-      const result = heart.scale(100 / (box.max[0] - box.min[0]));
+      const scale = 100 / (box.max[0] - box.min[0]);
+
+      setMorphStart(ball, func);
+      const node = new GLTFNode();
+      node.manifold = ball;
+      node.scale = [scale, scale, scale];
+
+      globalDefaults.animationLength = 5;  // seconds
+      globalDefaults.animationMode = 'ping-pong';
+
+      const result = heart.scale(scale);
       return result;
     },
 
@@ -253,71 +285,68 @@ export const examples = {
       // an example of using the warp() method, thus avoiding any direct
       // handling of triangles.
 
-      // @param p The number of times the thread passes through the donut hole.
-      // @param q The number of times the thread circles the donut.
-      // @param majorRadius Radius of the interior of the imaginary donut.
-      // @param minorRadius Radius of the small cross-section of the imaginary
-      //   donut.
-      // @param threadRadius Radius of the small cross-section of the actual
-      //   object.
-      // @param circularSegments Number of linear segments making up the
-      //   threadRadius circle. Default is getCircularSegments(threadRadius).
-      // @param linearSegments Number of segments along the length of the knot.
-      //   Default makes roughly square facets.
+      // The number of times the thread passes through the donut hole.
+      const p = 1;
+      // The number of times the thread circles the donut.
+      const q = 3;
+      // Radius of the interior of the imaginary donut.
+      const majorRadius = 25;
+      // Radius of the small cross-section of the imaginary donut.
+      const minorRadius = 10;
+      // Radius of the small cross-section of the actual object.
+      const threadRadius = 3.75;
+      // Number of linear segments making up the threadRadius circle. Default is
+      // getCircularSegments(threadRadius).
+      const circularSegments = -1;
+      // Number of segments along the length of the knot. Default makes roughly
+      // square facets.
+      const linearSegments = -1;
 
-      function torusKnot(
-          p, q, majorRadius, minorRadius, threadRadius, circularSegments = 0,
-          linearSegments = 0) {
-        const {vec3} = glMatrix;
-
-        function gcd(a, b) {
-          return b == 0 ? a : gcd(b, a % b);
-        }
-
-        const kLoops = gcd(p, q);
-        p /= kLoops;
-        q /= kLoops;
-        const n = circularSegments > 2 ? circularSegments :
-                                         getCircularSegments(threadRadius);
-        const m = linearSegments > 2 ? linearSegments :
-                                       n * q * majorRadius / threadRadius;
-
-        const offset = 2
-        const circle = CrossSection.circle(1, n).translate([offset, 0]);
-
-        const func = (v) => {
-          const psi = q * Math.atan2(v[0], v[1]);
-          const theta = psi * p / q;
-          const x1 = Math.sqrt(v[0] * v[0] + v[1] * v[1]);
-          const phi = Math.atan2(x1 - offset, v[2]);
-          vec3.set(
-              v, threadRadius * Math.cos(phi), 0, threadRadius * Math.sin(phi));
-          const center = vec3.fromValues(0, 0, 0);
-          const r = majorRadius + minorRadius * Math.cos(theta);
-          vec3.rotateX(v, v, center, -Math.atan2(p * minorRadius, q * r));
-          v[0] += minorRadius;
-          vec3.rotateY(v, v, center, theta);
-          v[0] += majorRadius;
-          vec3.rotateZ(v, v, center, psi);
-        };
-
-        let knot = Manifold.revolve(circle, m).warp(func);
-
-        if (kLoops > 1) {
-          const knots = [];
-          for (let k = 0; k < kLoops; ++k) {
-            knots.push(knot.rotate([0, 0, 360 * (k / kLoops) * (q / p)]));
-          }
-          knot = Manifold.compose(knots);
-        }
-
-        return knot;
-      }
-
-      // This recreates Matlab Knot by Emmett Lalish:
+      // These default values recreate Matlab Knot by Emmett Lalish:
       // https://www.thingiverse.com/thing:7080
 
-      const result = torusKnot(1, 3, 25, 10, 3.75);
+      const {vec3} = glMatrix;
+
+      function gcd(a, b) {
+        return b == 0 ? a : gcd(b, a % b);
+      }
+
+      const kLoops = gcd(p, q);
+      const pk = p / kLoops;
+      const qk = q / kLoops;
+      const n = circularSegments > 2 ? circularSegments :
+                                       getCircularSegments(threadRadius);
+      const m = linearSegments > 2 ? linearSegments :
+                                     n * qk * majorRadius / threadRadius;
+
+      const offset = 2
+      const circle = CrossSection.circle(1, n).translate([offset, 0]);
+
+      const func = (v) => {
+        const psi = qk * Math.atan2(v[0], v[1]);
+        const theta = psi * pk / qk;
+        const x1 = Math.sqrt(v[0] * v[0] + v[1] * v[1]);
+        const phi = Math.atan2(x1 - offset, v[2]);
+        vec3.set(
+            v, threadRadius * Math.cos(phi), 0, threadRadius * Math.sin(phi));
+        const center = vec3.fromValues(0, 0, 0);
+        const r = majorRadius + minorRadius * Math.cos(theta);
+        vec3.rotateX(v, v, center, -Math.atan2(pk * minorRadius, qk * r));
+        v[0] += minorRadius;
+        vec3.rotateY(v, v, center, theta);
+        v[0] += majorRadius;
+        vec3.rotateZ(v, v, center, psi);
+      };
+
+      const result = Manifold.revolve(circle, m).warp(func);
+
+      for (let k = 0; k < kLoops; ++k) {
+        const node = new GLTFNode();
+        node.manifold = result;
+        node.rotation = [0, 0, 360 * (k / kLoops) * (qk / pk)];
+        node.translation = [0, 0, result.boundingBox().min[2]];
+      }
+
       return result;
     },
 
@@ -350,12 +379,12 @@ export const examples = {
 
         const hole = Manifold.compose(holes);
 
-        result = Manifold.difference(
-            result,
-            hole,
-            hole.rotate([90, 0, 0]),
-            hole.rotate([0, 90, 0]),
-        );
+        result = Manifold.difference([
+          result,
+          hole,
+          hole.rotate([90, 0, 0]),
+          hole.rotate([0, 90, 0]),
+        ]);
         return result;
       }
 

bindings/wasm/examples/worker.test.js

@@ -18,7 +18,7 @@ import {WebIO} from '@gltf-transform/core';
 import {expect, suite, test} from 'vitest';
 
 import Module from './built/manifold.js';
-import {readMesh, setupIO} from './gltf-io.js';
+import {readMesh, setupIO} from './gltf-io';
 import {examples} from './public/examples.js';
 import ManifoldWorker from './worker?worker';
 
@@ -63,18 +63,16 @@ async function runExample(name) {
         const docIn = await io.read(glbURL);
         const nodes = docIn.getRoot().listNodes();
         for (const node of nodes) {
-          const mesh = node.getMesh();
-          if (!mesh) {
+          const docMesh = node.getMesh();
+          if (!docMesh) {
             continue;
           }
-          const attributes = [];
-          const materials = [];
-          const manifoldMesh = readMesh(mesh, attributes, materials);
-          const manifold = wasm.Manifold(manifoldMesh);
+          const {mesh} = readMesh(docMesh);
+          const manifold = wasm.Manifold(mesh);
           const prop = manifold.getProperties();
           const genus = manifold.genus();
-          console.log(genus);
           manifold.delete();
+          // Return properties of first mesh encountered.
           resolve({...prop, genus});
         }
       } catch (e) {
@@ -97,8 +95,8 @@ suite('Examples', () => {
   test('Tetrahedron Puzzle', async () => {
     const result = await runExample('Tetrahedron Puzzle');
     expect(result.genus).to.equal(0, 'Genus');
-    expect(result.volume).to.be.closeTo(7297, 1, 'Volume');
-    expect(result.surfaceArea).to.be.closeTo(3303, 1, 'Surface Area');
+    expect(result.volume).to.be.closeTo(7240, 1, 'Volume');
+    expect(result.surfaceArea).to.be.closeTo(3235, 1, 'Surface Area');
   });
 
   test('Rounded Frame', async () => {
@@ -111,8 +109,8 @@ suite('Examples', () => {
   test('Heart', async () => {
     const result = await runExample('Heart');
     expect(result.genus).to.equal(0, 'Genus');
-    expect(result.volume).to.be.closeTo(282743, 10, 'Volume');
-    expect(result.surfaceArea).to.be.closeTo(22187, 1, 'Surface Area');
+    expect(result.volume).to.be.closeTo(3.342, 0.001, 'Volume');
+    expect(result.surfaceArea).to.be.closeTo(11.51, 0.01, 'Surface Area');
   });
 
   test('Scallop', async () => {

bindings/wasm/examples/worker.ts

@@ -12,27 +12,38 @@
 // See the License for the specific language governing permissions and
 // limitations under the License.
 
-import {Document, mat4, Material, Node, WebIO} from '@gltf-transform/core';
+import {Accessor, Animation, AnimationSampler, Document, mat4, Material, Mesh as GLTFMesh, Node, WebIO} from '@gltf-transform/core';
 import {KHRMaterialsUnlit, KHRONOS_EXTENSIONS} from '@gltf-transform/extensions';
-import {fileForContentTypes, to3dmodel} from '@jscadui/3mf-export';
+import {fileForContentTypes, FileForRelThumbnail, to3dmodel} from '@jscadui/3mf-export';
 import {strToU8, Zippable, zipSync} from 'fflate'
 import * as glMatrix from 'gl-matrix';
 
 import Module from './built/manifold';
-//@ts-ignore
-import {setupIO, writeMesh} from './gltf-io';
-import type {GLTFMaterial, Quat} from './public/editor';
+import {Properties, setupIO, writeMesh} from './gltf-io';
+import {GLTFMaterial, Quat} from './public/editor';
 import type {CrossSection, Manifold, ManifoldToplevel, Mesh, Vec3} from './public/manifold';
 
+interface GlobalDefaults {
+  roughness: number;
+  metallic: number;
+  baseColorFactor: [number, number, number];
+  alpha: number;
+  unlit: boolean;
+  animationLength: number;
+  animationMode: 'loop'|'ping-pong';
+}
+
 interface WorkerStatic extends ManifoldToplevel {
   GLTFNode: typeof GLTFNode;
   show(manifold: Manifold): Manifold;
   only(manifold: Manifold): Manifold;
-  setMaterial(manifold: Manifold, material: GLTFMaterial): void;
+  setMaterial(manifold: Manifold, material: GLTFMaterial): Manifold;
+  setMorphStart(manifold: Manifold, func: (v: Vec3) => void): void;
+  setMorphEnd(manifold: Manifold, func: (v: Vec3) => void): void;
   cleanup(): void;
 }
 
-const module = await Module() as unknown as WorkerStatic;
+const module = await Module() as WorkerStatic;
 module.setup();
 
 // Faster on modern browsers than Float32Array
@@ -50,7 +61,7 @@ const manifoldStaticFunctions = [
 const manifoldMemberFunctions = [
   'add', 'subtract', 'intersect', 'decompose', 'warp', 'transform', 'translate',
   'rotate', 'scale', 'mirror', 'refine', 'setProperties', 'asOriginal',
-  'trimByPlane', 'split', 'splitByPlane', 'hull'
+  'trimByPlane', 'split', 'splitByPlane', 'slice', 'project', 'hull'
 ];
 // CrossSection static methods (that return a new cross-section)
 const crossSectionStaticFunctions = [
@@ -66,7 +77,8 @@ const crossSectionMemberFunctions = [
 const toplevelConstructors = ['show', 'only', 'setMaterial'];
 const toplevel = [
   'setMinCircularAngle', 'setMinCircularEdgeLength', 'setCircularSegments',
-  'getCircularSegments', 'Mesh', 'GLTFNode', 'Manifold', 'CrossSection'
+  'getCircularSegments', 'Mesh', 'GLTFNode', 'Manifold', 'CrossSection',
+  'setMorphStart', 'setMorphEnd'
 ];
 const exposedFunctions = toplevelConstructors.concat(toplevel);
 
@@ -83,15 +95,13 @@ function addMembers(
   const cls = module[className];
   const obj = areStatic ? cls : cls.prototype;
   for (const name of methodNames) {
-    if (name != 'cylinder') {
-      const originalFn = obj[name];
-      obj[name] = function(...args: any) {
-        //@ts-ignore
-        const result = originalFn(...args);
-        memoryRegistry.push(result);
-        return result;
-      };
-    }
+    const originalFn = obj[name];
+    obj[name] = function(...args: any) {
+      //@ts-ignore
+      const result = originalFn(...args);
+      memoryRegistry.push(result);
+      return result;
+    };
   }
 }
 
@@ -127,6 +137,18 @@ let ghost = false;
 const shown = new Map<number, Mesh>();
 const singles = new Map<number, Mesh>();
 
+const FPS = 30;
+
+const GLOBAL_DEFAULTS = {
+  roughness: 0.2,
+  metallic: 1,
+  baseColorFactor: [1, 1, 0] as [number, number, number],
+  alpha: 1,
+  unlit: false,
+  animationLength: 1,
+  animationMode: 'loop'
+};
+
 const SHOW = {
   baseColorFactor: [1, 0, 0],
   alpha: 0.25,
@@ -144,8 +166,15 @@ const GHOST = {
 const nodes = new Array<GLTFNode>();
 const id2material = new Map<number, GLTFMaterial>();
 const materialCache = new Map<GLTFMaterial, Material>();
-let nextGlobalID = 0;
 const object2globalID = new Map<GLTFNode|Manifold, number>();
+const manifold2morph = new Map<Manifold, Morph>();
+// lib3mf doesn't like objectid=0
+let nextGlobalID = 1;
+let animation: Animation;
+let timesAccessor: Accessor;
+let weightsAccessor: Accessor;
+let weightsSampler: AnimationSampler;
+let hasAnimation: boolean;
 
 function cleanup() {
   ghost = false;
@@ -155,7 +184,13 @@ function cleanup() {
   id2material.clear();
   materialCache.clear();
   object2globalID.clear();
-  nextGlobalID = 0;
+  manifold2morph.clear();
+  nextGlobalID = 1;
+}
+
+interface Morph {
+  start?: (v: Vec3) => void;
+  end?: (v: Vec3) => void;
 }
 
 interface Mesh3MF {
@@ -198,9 +233,9 @@ interface To3MF {
 class GLTFNode {
   private _parent?: GLTFNode;
   manifold?: Manifold;
-  translation?: Vec3;
-  rotation?: Vec3;
-  scale?: Vec3;
+  translation?: Vec3|((t: number) => Vec3);
+  rotation?: Vec3|((t: number) => Vec3);
+  scale?: Vec3|((t: number) => Vec3);
   material?: GLTFMaterial;
   name?: string;
 
@@ -221,12 +256,32 @@ class GLTFNode {
 
 module.GLTFNode = GLTFNode;
 
+const globalDefaults = {...GLOBAL_DEFAULTS};
+
 module.setMaterial = (manifold: Manifold, material: GLTFMaterial): Manifold => {
   const out = manifold.asOriginal();
   id2material.set(out.originalID(), material);
   return out;
 };
 
+module.setMorphStart = (manifold: Manifold, func: (v: Vec3) => void): void => {
+  const morph = manifold2morph.get(manifold);
+  if (morph != null) {
+    morph.start = func;
+  } else {
+    manifold2morph.set(manifold, {start: func});
+  }
+};
+
+module.setMorphEnd = (manifold: Manifold, func: (v: Vec3) => void): void => {
+  const morph = manifold2morph.get(manifold);
+  if (morph != null) {
+    morph.end = func;
+  } else {
+    manifold2morph.set(manifold, {end: func});
+  }
+};
+
 function debug(manifold: Manifold, map: Map<number, Mesh>) {
   let result = manifold.asOriginal();
   map.set(result.originalID(), result.getMesh());
@@ -271,8 +326,8 @@ function log(...args: any[]) {
 }
 
 self.onmessage = async (e) => {
-  const content = e.data +
-      '\nreturn exportModels(typeof result === "undefined" ? undefined : result);\n';
+  const content = 'const globalDefaults = {};\n' + e.data +
+      '\nreturn exportModels(globalDefaults, typeof result === "undefined" ? undefined : result);\n';
   try {
     const f = new Function(
         'exportModels', 'glMatrix', 'module', ...exposedFunctions, content);
@@ -288,23 +343,146 @@ self.onmessage = async (e) => {
   }
 };
 
-function createGLTFnode(doc: Document, node: GLTFNode) {
+function euler2quat(rotation: Vec3): Quat {
+  const {quat} = glMatrix;
+  const deg2rad = Math.PI / 180;
+  const q = quat.create() as Quat;
+  quat.rotateZ(q, q, deg2rad * rotation[2]);
+  quat.rotateY(q, q, deg2rad * rotation[1]);
+  quat.rotateX(q, q, deg2rad * rotation[0]);
+  return q;
+}
+
+function addMotion(
+    doc: Document, type: 'translation'|'rotation'|'scale', node: GLTFNode,
+    out: Node): Vec3|null {
+  const motion = node[type];
+  if (motion == null) {
+    return null;
+  }
+  if (typeof motion !== 'function') {
+    return motion;
+  }
+
+  const nFrames = timesAccessor.getCount();
+  const nEl = type == 'rotation' ? 4 : 3;
+  const frames = new Float32Array(nEl * nFrames);
+  for (let i = 0; i < nFrames; ++i) {
+    const x = i / (nFrames - 1);
+    const m = motion(
+        globalDefaults.animationMode !== 'ping-pong' ?
+            x :
+            (1 - Math.cos(x * 2 * Math.PI)) / 2);
+    frames.set(nEl === 4 ? euler2quat(m) : m, nEl * i);
+  }
+
+  const framesAccessor =
+      doc.createAccessor(node.name + ' ' + type + ' frames')
+          .setBuffer(doc.getRoot().listBuffers()[0])
+          .setArray(frames)
+          .setType(nEl === 4 ? Accessor.Type.VEC4 : Accessor.Type.VEC3);
+  const sampler = doc.createAnimationSampler()
+                      .setInput(timesAccessor)
+                      .setOutput(framesAccessor)
+                      .setInterpolation('LINEAR');
+  const channel = doc.createAnimationChannel()
+                      .setTargetPath(type)
+                      .setTargetNode(out)
+                      .setSampler(sampler);
+  animation.addSampler(sampler);
+  animation.addChannel(channel);
+  hasAnimation = true;
+  return motion(0);
+}
+
+function setMorph(doc: Document, node: Node, manifold: Manifold) {
+  if (manifold2morph.has(manifold)) {
+    const channel = doc.createAnimationChannel()
+                        .setTargetPath('weights')
+                        .setTargetNode(node)
+                        .setSampler(weightsSampler);
+    animation.addChannel(channel);
+    hasAnimation = true;
+  }
+}
+
+function morphStart(manifoldMesh: Mesh, morph?: Morph): number[] {
+  const inputPositions: number[] = [];
+  if (morph == null) {
+    return inputPositions;
+  }
+
+  for (let i = 0; i < manifoldMesh.numVert; ++i) {
+    for (let j = 0; j < 3; ++j)
+      inputPositions[i * 3 + j] =
+          manifoldMesh.vertProperties[i * manifoldMesh.numProp + j];
+  }
+  if (morph.start) {
+    for (let i = 0; i < manifoldMesh.numVert; ++i) {
+      const vertProp = manifoldMesh.vertProperties;
+      const offset = i * manifoldMesh.numProp;
+      const pos = inputPositions.slice(offset, offset + 3) as Vec3;
+      morph.start(pos);
+      for (let j = 0; j < 3; ++j) vertProp[offset + j] = pos[j];
+    }
+  }
+  return inputPositions;
+}
+
+function morphEnd(
+    doc: Document, manifoldMesh: Mesh, mesh: GLTFMesh, inputPositions: number[],
+    morph?: Morph) {
+  if (morph == null) {
+    return;
+  }
+
+  mesh.setWeights([0]);
+
+  mesh.listPrimitives().forEach((primitive, i) => {
+    if (morph.end) {
+      for (let i = 0; i < manifoldMesh.numVert; ++i) {
+        const pos = inputPositions.slice(3 * i, 3 * (i + 1)) as Vec3;
+        morph.end(pos);
+        inputPositions.splice(3 * i, 3, ...pos);
+      }
+    }
+
+    const startPosition = primitive.getAttribute('POSITION')!.getArray()!;
+    const array = new Float32Array(startPosition.length);
+
+    const offset = manifoldMesh.runIndex[i];
+    for (let j = 0; j < array.length; ++j) {
+      array[j] = inputPositions[offset + j] - startPosition[j];
+    }
+
+    const morphAccessor = doc.createAccessor(mesh.getName() + ' morph target')
+                              .setBuffer(doc.getRoot().listBuffers()[0])
+                              .setArray(array)
+                              .setType(Accessor.Type.VEC3);
+    const morphTarget =
+        doc.createPrimitiveTarget().setAttribute('POSITION', morphAccessor);
+    primitive.addTarget(morphTarget);
+  });
+}
+
+function createGLTFnode(doc: Document, node: GLTFNode): Node {
   const out = doc.createNode(node.name);
-  if (node.translation) {
-    out.setTranslation(node.translation);
+
+  const pos = addMotion(doc, 'translation', node, out);
+  if (pos != null) {
+    out.setTranslation(pos);
   }
-  if (node.rotation) {
-    const {quat} = glMatrix;
-    const deg2rad = Math.PI / 180;
-    const q = quat.create() as Quat;
-    quat.rotateX(q, q, deg2rad * node.rotation[0]);
-    quat.rotateY(q, q, deg2rad * node.rotation[1]);
-    quat.rotateZ(q, q, deg2rad * node.rotation[2]);
-    out.setRotation(q);
+
+  const rot = addMotion(doc, 'rotation', node, out);
+  if (rot != null) {
+    out.setRotation(euler2quat(rot));
   }
-  if (node.scale) {
-    out.setScale(node.scale);
+
+  const scale = addMotion(doc, 'scale', node, out);
+  if (scale != null) {
+    out.setScale(scale);
   }
+
   return out;
 }
 
@@ -318,15 +496,13 @@ function getBackupMaterial(node?: GLTFNode): GLTFMaterial {
   return node.material;
 }
 
-function makeDefaultedMaterial(doc: Document, {
-  roughness = 0.2,
-  metallic = 1,
-  baseColorFactor = [1, 1, 0],
-  alpha = 1,
-  unlit = false,
-  name = ''
-}: GLTFMaterial = {}): Material {
-  const material = doc.createMaterial(name);
+function makeDefaultedMaterial(
+    doc: Document, matIn: GLTFMaterial = {}): Material {
+  const defaults = {...globalDefaults};
+  Object.assign(defaults, matIn);
+  const {roughness, metallic, baseColorFactor, alpha, unlit} = defaults;
+
+  const material = doc.createMaterial(matIn.name ?? '');
 
   if (unlit) {
     const unlit = doc.createExtension(KHRMaterialsUnlit).createUnlit();
@@ -373,19 +549,23 @@ function addMesh(
   // From Z-up to Y-up (glTF)
   const manifoldMesh = manifold.getMesh();
 
-  const materials = new Array<GLTFMaterial>();
-  const attributes = new Array<Array<string>>();
+  const id2properties = new Map<number, Properties>();
   for (const id of manifoldMesh.runOriginalID!) {
     const material = id2material.get(id) || backupMaterial;
-    materials.push(material);
-    attributes.push(['POSITION', ...material.attributes ?? []]);
+    id2properties.set(id, {
+      material: getCachedMaterial(doc, ghost ? GHOST : material),
+      attributes: ['POSITION', ...material.attributes ?? []]
+    });
   }
 
-  const gltfMaterials = materials.map((matDef) => {
-    return getCachedMaterial(doc, ghost ? GHOST : matDef);
-  });
+  const morph = manifold2morph.get(manifold);
+  const inputPositions = morphStart(manifoldMesh, morph);
 
-  node.setMesh(writeMesh(doc, manifoldMesh, attributes, gltfMaterials));
+  const mesh = writeMesh(doc, manifoldMesh, id2properties);
+  mesh.setName(node.getName());
+  node.setMesh(mesh);
+
+  morphEnd(doc, manifoldMesh, mesh, inputPositions, morph);
 
   const vertices = manifoldMesh.numProp === 3 ?
       manifoldMesh.vertProperties :
@@ -402,21 +582,18 @@ function addMesh(
       {vertices, indices: manifoldMesh.triVerts, id: `${nextGlobalID++}`});
 
   for (const [run, id] of manifoldMesh.runOriginalID!.entries()) {
-    let inMesh = shown.get(id);
-    let single = false;
-    if (inMesh == null) {
-      single = true;
-      inMesh = singles.get(id);
-    }
+    const show = shown.has(id);
+    const inMesh = show ? shown.get(id) : singles.get(id);
     if (inMesh == null) {
       continue;
     }
-    const mat = single ? materials[run] : SHOW;
-    const debugNode =
-        doc.createNode('debug')
-            .setMesh(writeMesh(
-                doc, inMesh, attributes, [getCachedMaterial(doc, mat)]))
-            .setMatrix(manifoldMesh.transform(run));
+
+    id2properties.get(id)!.material = getCachedMaterial(
+        doc, show ? SHOW : (id2material.get(id) || backupMaterial));
+
+    const debugNode = doc.createNode('debug')
+                          .setMesh(writeMesh(doc, inMesh, id2properties))
+                          .setMatrix(manifoldMesh.transform(run));
     node.addChild(debugNode);
   }
 }
@@ -450,6 +627,7 @@ function createNodeFromCache(
   const node = createGLTFnode(doc, nodeDef);
   const {manifold} = nodeDef;
   if (manifold != null) {
+    setMorph(doc, node, manifold);
     const backupMaterial = getBackupMaterial(nodeDef);
     const cachedNodes = manifold2node.get(manifold);
     if (cachedNodes == null) {
@@ -481,7 +659,10 @@ function createNodeFromCache(
   return node;
 }
 
-async function exportModels(manifold?: Manifold) {
+async function exportModels(defaults: GlobalDefaults, manifold?: Manifold) {
+  Object.assign(globalDefaults, GLOBAL_DEFAULTS);
+  Object.assign(globalDefaults, defaults);
+
   const doc = new Document();
   const halfRoot2 = Math.sqrt(2) / 2;
   const mm2m = 1 / 1000;
@@ -490,6 +671,32 @@ async function exportModels(manifold?: Manifold) {
                       .setScale([mm2m, mm2m, mm2m]);
   doc.createScene().addChild(wrapper);
 
+  animation = doc.createAnimation('');
+  hasAnimation = false;
+  const nFrames = Math.round(globalDefaults.animationLength * FPS) + 1;
+  const times = new Float32Array(nFrames);
+  const weights = new Float32Array(nFrames);
+  for (let i = 0; i < nFrames; ++i) {
+    const x = i / (nFrames - 1);
+    times[i] = x * globalDefaults.animationLength;
+    weights[i] = globalDefaults.animationMode !== 'ping-pong' ?
+        x :
+        (1 - Math.cos(x * 2 * Math.PI)) / 2;
+  }
+  timesAccessor = doc.createAccessor('animation times')
+                      .setBuffer(doc.createBuffer())
+                      .setArray(times)
+                      .setType(Accessor.Type.SCALAR);
+  weightsAccessor = doc.createAccessor('animation weights')
+                        .setBuffer(doc.getRoot().listBuffers()[0])
+                        .setArray(weights)
+                        .setType(Accessor.Type.SCALAR);
+  weightsSampler = doc.createAnimationSampler()
+                       .setInput(timesAccessor)
+                       .setOutput(weightsAccessor)
+                       .setInterpolation('LINEAR');
+  animation.addSampler(weightsSampler);
+
   const to3mf = {
     meshes: [],
     components: [],
@@ -547,13 +754,24 @@ async function exportModels(manifold?: Manifold) {
     to3mf.items.push({objectID: `${object2globalID.get(manifold)}`});
   }
 
+  if (!hasAnimation) {
+    timesAccessor.dispose();
+    weightsAccessor.dispose();
+    weightsSampler.dispose();
+    animation.dispose();
+  }
+
   const glb = await io.writeBinary(doc);
   const blobGLB = new Blob([glb], {type: 'application/octet-stream'});
 
+  const fileForRelThumbnail = new FileForRelThumbnail();
+  fileForRelThumbnail.add3dModel('3D/3dmodel.model')
+
   const model = to3dmodel(to3mf);
   const files: Zippable = {};
   files['3D/3dmodel.model'] = strToU8(model);
   files[fileForContentTypes.name] = strToU8(fileForContentTypes.content);
+  files[fileForRelThumbnail.name] = strToU8(fileForRelThumbnail.content);
   const zipFile = zipSync(files);
   const blob3MF = new Blob(
       [zipFile],

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

@@ -12,7 +12,7 @@
 // See the License for the specific language governing permissions and
 // limitations under the License.
 
-import {Box, FillRule, JoinType, Mat4, Polygons, Properties, Rect, SealedFloat32Array, SealedUint32Array, SimplePolygon, Smoothness, Vec2, Vec3} from './manifold-global-types';
+import {Box, FillRule, JoinType, Mat3, Mat4, Polygons, Properties, Rect, SealedFloat32Array, SealedUint32Array, SimplePolygon, Smoothness, Vec2, Vec3} from './manifold-global-types';
 
 /**
  * Triangulates a set of /epsilon-valid polygons.
@@ -685,6 +685,22 @@ export class Manifold {
    */
   trimByPlane(normal: Vec3, originOffset: number): Manifold;
 
+  /**
+   * Returns the cross section of this object parallel to the X-Y plane at the
+   * specified height. Using a height equal to the bottom
+   * of the bounding box will return the bottom faces, while using a height
+   * equal to the top of the bounding box will return empty.
+   *
+   * @param height Z-level of slice.
+   */
+  slice(height: number): CrossSection;
+
+  /**
+   * Returns a cross section representing the projected outline of this object
+   * onto the X-Y plane.
+   */
+  project(): CrossSection;
+
   // Convex Hulls
 
   /**
@@ -838,19 +854,7 @@ export class Manifold {
   delete(): void;
 }
 
-export class Mesh {
-  constructor(options: {
-    numProp: number,
-    vertProperties: Float32Array,
-    triVerts: Uint32Array,
-    mergeFromVert?: Uint32Array,
-    mergeToVert?: Uint32Array,
-    runIndex?: Uint32Array,
-    runOriginalID?: Uint32Array,
-    runTransform?: Float32Array,
-    faceID?: Uint32Array,
-    halfedgeTangent?: Float32Array
-  });
+export interface MeshOptions {
   numProp: number;
   vertProperties: Float32Array;
   triVerts: Uint32Array;
@@ -861,6 +865,20 @@ export class Mesh {
   runTransform?: Float32Array;
   faceID?: Uint32Array;
   halfedgeTangent?: Float32Array;
+}
+
+export class Mesh {
+  constructor(options: MeshOptions);
+  numProp: number;
+  vertProperties: Float32Array;
+  triVerts: Uint32Array;
+  mergeFromVert: Uint32Array;
+  mergeToVert: Uint32Array;
+  runIndex: Uint32Array;
+  runOriginalID: Uint32Array;
+  runTransform: Float32Array;
+  faceID: Uint32Array;
+  halfedgeTangent: Float32Array;
   get numTri(): number;
   get numVert(): number;
   get numRun(): number;

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

@@ -22,6 +22,18 @@ export interface SealedFloat32Array<N extends number> extends Float32Array {
 
 export type Vec2 = [number, number];
 export type Vec3 = [number, number, number];
+// 3x3 matrix stored in column-major order
+export type Mat3 = [
+  number,
+  number,
+  number,
+  number,
+  number,
+  number,
+  number,
+  number,
+  number,
+];
 // 4x4 matrix stored in column-major order
 export type Mat4 = [
   number,

bindings/wasm/package.json

@@ -1,6 +1,6 @@
 {
   "name": "manifold-3d",
-  "version": "2.1.0",
+  "version": "2.3.0",
   "description": "Geometry library for topological robustness",
   "main": "manifold.js",
   "files": [

extras/CMakeLists.txt

@@ -41,8 +41,9 @@ endif()
 if(BUILD_TEST_CGAL)
     add_executable(perfTestCGAL perf_test_cgal.cpp)
     find_package(CGAL REQUIRED COMPONENTS Core)
+    find_package(Boost REQUIRED COMPONENTS thread)
     target_compile_definitions(perfTestCGAL PRIVATE CGAL_USE_GMPXX)
-    target_link_libraries(perfTestCGAL manifold CGAL::CGAL CGAL::CGAL_Core)
+    target_link_libraries(perfTestCGAL manifold CGAL::CGAL CGAL::CGAL_Core Boost::thread)
     target_compile_options(perfTestCGAL PRIVATE ${MANIFOLD_FLAGS})
     target_compile_features(perfTestCGAL PUBLIC cxx_std_17)
 endif()

flake.nix

@@ -52,7 +52,7 @@
             }: pkgs.stdenv.mkDerivation {
               inherit doCheck;
               pname = "manifold-${parallel-backend}";
-              version = "2.2.0";
+              version = "2.3.0";
               src = self;
               nativeBuildInputs = (with pkgs; [
                 cmake
@@ -60,7 +60,7 @@
                 glm
                 ninja
                 (python39.withPackages
-                  (ps: with ps; [ trimesh ]))
+                  (ps: with ps; [ trimesh pytest ]))
                 gtest
               ]) ++ build-tools;
               cmakeFlags = [
@@ -74,6 +74,7 @@
                 ./manifold_test
                 cd ../../
                 PYTHONPATH=$PYTHONPATH:$(pwd)/build/bindings/python python3 bindings/python/examples/run_all.py
+                PYTHONPATH=$PYTHONPATH:$(pwd)/build/bindings/python python3 -m pytest
                 cd build
               '';
             };
@@ -103,7 +104,7 @@
               parallelBackends)) // {
             manifold-js = pkgs.buildEmscriptenPackage {
               name = "manifold-js";
-              version = "beta";
+              version = "2.3.0";
               src = self;
               nativeBuildInputs = (with pkgs; [ cmake python39 ]);
               buildInputs = [ pkgs.nodejs ];
@@ -134,7 +135,7 @@
             # but how should we make it work with other python versions?
             manifold3d = with pkgs.python3Packages; buildPythonPackage {
               pname = "manifold3d";
-              version = "2.2.0";
+              version = "2.3.0";
               src = self;
               propagatedBuildInputs = [
                 numpy
@@ -155,12 +156,14 @@
               ];
               checkInputs = [
                 trimesh
+                pytest
               ];
               format = "pyproject";
               dontUseCmakeConfigure = true;
               doCheck = true;
               checkPhase = ''
                 python3 bindings/python/examples/run_all.py
+                python3 -m pytest
               '';
             };
           };

format.sh

@@ -0,0 +1,16 @@
+#!/usr/bin/env bash
+
+shopt -s extglob
+if [ -z "$CLANG_FORMAT" ]; then
+CLANG_FORMAT=$(dirname $(which clang-11))/clang-format
+fi
+
+$CLANG_FORMAT -i extras/*.cpp
+$CLANG_FORMAT -i meshIO/**/*.{h,cpp}
+$CLANG_FORMAT -i samples/**/*.{h,cpp}
+$CLANG_FORMAT -i test/*.{h,cpp}
+$CLANG_FORMAT -i bindings/*/*.cpp
+$CLANG_FORMAT -i bindings/c/include/*.h
+$CLANG_FORMAT -i bindings/wasm/**/*.{js,ts,html}
+$CLANG_FORMAT -i src/!(third_party)/*/*.{h,cpp}
+black bindings/python/examples/*.py

meshIO/src/meshIO.cpp

@@ -234,6 +234,7 @@ void ExportMesh(const std::string& filename, const MeshGL& mesh,
                    ? 1
                    : mesh.vertProperties[i * mesh.numProp +
                                          options.mat.colorChannels[j]];
+      c = glm::saturate(c);
       mesh_out->mColors[0][i] = aiColor4D(c.r, c.g, c.b, c.a);
     }
   }

minimizer.sh

@@ -0,0 +1,59 @@
+#!/usr/bin/env bash
+
+# Modified from https://github.com/google/fuzztest/blob/main/tools/minimizer.sh
+# If you want to use a specific reproducer, pass it in FUZZTEST_MINIMIZE_REPRODUCER
+# environment variable.
+# Otherwise, just run it as
+# $ ./minimizer.sh ./build/test/polygon_fuzz --fuzz=PolygonFuzz.TriangulationNoCrash
+# If you want to kill it, try SIGKILL, SIGTERM is not very useful
+# If you want to end the loop early, rename $(pwd)/reproducers into something
+# else...
+
+if [ -z "${FUZZTEST_MINIMIZE_REPRODUCER}"]; then
+  mkdir -p $(pwd)/reproducers
+  FUZZTEST_REPRODUCERS_OUT_DIR=$(pwd)/reproducers "$@"
+  wildcard="$(pwd)/reproducers/*"
+  reproducers=($wildcard)
+  FUZZTEST_MINIMIZE_REPRODUCER="${reproducers[0]}"
+fi
+readonly ORIGINAL_REPRODUCER="${FUZZTEST_MINIMIZE_REPRODUCER}"
+
+for i in {0001..9999}; do
+  echo
+  echo "╔════════════════════════════════════════════════╗"
+  echo "║ Minimization round: ${i}                       ║"
+  echo "╚════════════════════════════════════════════════╝"
+  echo "Note that to terminate early, simply move $ORIGINAL_REPRODUCER to somewhere else..."
+  echo
+
+  if [ ! -f $ORIGINAL_REPRODUCER ]; then
+    echo "Terminated by the user"
+    break
+  fi
+
+  TEMP_DIR=$(mktemp -d)
+  FUZZTEST_REPRODUCERS_OUT_DIR="${TEMP_DIR}" \
+  FUZZTEST_MINIMIZE_REPRODUCER="${FUZZTEST_MINIMIZE_REPRODUCER}" \
+  "$@"
+
+  if [ $? -eq 130 ]; then
+    echo
+    echo "╔═══════════════════════════════════════════════╗"
+    echo "║ Minimization terminated.                      ║"
+    echo "╚═══════════════════════════════════════════════╝"
+    echo
+    echo "Find the smallest reproducer at:"
+    echo
+    echo "${FUZZTEST_MINIMIZE_REPRODUCER}"
+
+    rm -rf "${TEMP_DIR}"
+    break
+  fi
+
+  SMALLER_REPRODUCER=$(find "${TEMP_DIR}" -type f)
+  NEW_NAME="${ORIGINAL_REPRODUCER}-min-${i}"
+  mv "${SMALLER_REPRODUCER}" "${NEW_NAME}"
+  FUZZTEST_MINIMIZE_REPRODUCER="${NEW_NAME}"
+
+  rm -rf "${TEMP_DIR}"
+done

pyproject.toml

@@ -1,13 +1,13 @@
 [project]
 name = "manifold3d"
-version = "2.2.0"
+version = "2.3.0"
 authors = [
   { name="Emmett Lalish", email="elalish@gmail.com" },
 ]
 description = "Library for geometric robustness"
 readme = "README.md"
 classifiers = [
-    "Development Status :: 4 - Beta",
+    "Development Status :: 5 - Production/Stable",
     "License :: OSI Approved :: Apache Software License",
     "Operating System :: OS Independent",
     "Programming Language :: C++",

pytest.ini

@@ -0,0 +1,3 @@
+[pytest]
+testpaths =
+    bindings/python/examples

src/collider/include/collider.h

@@ -38,7 +38,9 @@ class Collider {
   Vec<thrust::pair<int, int>> internalChildren_;
 
   int NumInternal() const { return internalChildren_.size(); };
-  int NumLeaves() const { return NumInternal() + 1; };
+  int NumLeaves() const {
+    return internalChildren_.empty() ? 0 : (NumInternal() + 1);
+  };
 };
 
 }  // namespace manifold

src/collider/src/collider.cpp

@@ -276,6 +276,7 @@ namespace manifold {
  */
 Collider::Collider(const VecView<const Box>& leafBB,
                    const VecView<const uint32_t>& leafMorton) {
+  ZoneScoped;
   ASSERT(leafBB.size() == leafMorton.size(), userErr,
          "vectors must be the same length");
   int num_nodes = 2 * leafBB.size() - 1;
@@ -299,6 +300,7 @@ Collider::Collider(const VecView<const Box>& leafBB,
  */
 template <const bool selfCollision, const bool inverted, typename T>
 SparseIndices Collider::Collisions(const VecView<const T>& queriesIn) const {
+  ZoneScoped;
   // note that the length is 1 larger than the number of queries so the last
   // element can store the sum when using exclusive scan
   if (queriesIn.size() < kSequentialThreshold) {
@@ -336,6 +338,7 @@ SparseIndices Collider::Collisions(const VecView<const T>& queriesIn) const {
  * hierarchy.
  */
 void Collider::UpdateBoxes(const VecView<const Box>& leafBB) {
+  ZoneScoped;
   ASSERT(leafBB.size() == NumLeaves(), userErr,
          "must have the same number of updated boxes as original");
   // copy in leaf node Boxes
@@ -355,6 +358,7 @@ void Collider::UpdateBoxes(const VecView<const Box>& leafBB) {
  * axis-aligned, abort and return false to indicate recalculation is necessary.
  */
 bool Collider::Transform(glm::mat4x3 transform) {
+  ZoneScoped;
   bool axisAligned = true;
   for (int row : {0, 1, 2}) {
     int count = 0;

src/cross_section/include/cross_section.h

@@ -21,6 +21,7 @@
 #include "glm/ext/matrix_float3x2.hpp"
 #include "glm/ext/vector_float2.hpp"
 #include "public.h"
+#include "vec_view.h"
 
 namespace manifold {
 
@@ -97,6 +98,8 @@ class CrossSection {
   CrossSection Mirror(const glm::vec2 ax) const;
   CrossSection Transform(const glm::mat3x2& m) const;
   CrossSection Warp(std::function<void(glm::vec2&)> warpFunc) const;
+  CrossSection WarpBatch(
+      std::function<void(VecView<glm::vec2>)> warpFunc) const;
   CrossSection Simplify(double epsilon = 1e-6) const;
 
   // Adapted from Clipper2 docs:

src/cross_section/src/cross_section.cpp

@@ -146,6 +146,18 @@ void decompose_hole(const C2::PolyTreeD* outline,
   }
 }
 
+void flatten(const C2::PolyTreeD* tree, C2::PathsD& polys, int i) {
+  auto n_outlines = tree->Count();
+  if (i < n_outlines) {
+    auto outline = tree->Child(i);
+    flatten(outline, polys, 0);
+    polys.push_back(outline->Polygon());
+    if (i < n_outlines - 1) {
+      flatten(tree, polys, i + 1);
+    }
+  }
+}
+
 bool V2Lesser(glm::vec2 a, glm::vec2 b) {
   if (a.x == b.x) return a.y < b.y;
   return a.x < b.x;
@@ -555,21 +567,42 @@ CrossSection CrossSection::Transform(const glm::mat3x2& m) const {
  */
 CrossSection CrossSection::Warp(
     std::function<void(glm::vec2&)> warpFunc) const {
-  auto paths = GetPaths();
-  auto warped = C2::PathsD();
-  warped.reserve(paths->paths_.size());
-  for (auto path : paths->paths_) {
-    auto sz = path.size();
-    auto s = C2::PathD(sz);
-    for (int i = 0; i < sz; ++i) {
-      auto v = v2_of_pd(path[i]);
-      warpFunc(v);
-      s[i] = v2_to_pd(v);
+  return WarpBatch([&warpFunc](VecView<glm::vec2> vecs) {
+    for (glm::vec2& p : vecs) {
+      warpFunc(p);
+    }
+  });
+}
+
+/**
+ * Same as CrossSection::Warp but calls warpFunc with
+ * a VecView which is roughly equivalent to std::span
+ * pointing to all vec2 elements to be modified in-place
+ *
+ * @param warpFunc A function that modifies multiple vertex positions.
+ */
+CrossSection CrossSection::WarpBatch(
+    std::function<void(VecView<glm::vec2>)> warpFunc) const {
+  std::vector<glm::vec2> tmp_verts;
+  C2::PathsD paths = GetPaths()->paths_;  // deep copy
+  for (C2::PathD const& path : paths) {
+    for (C2::PointD const& p : path) {
+      tmp_verts.push_back(v2_of_pd(p));
     }
-    warped.push_back(s);
   }
+
+  warpFunc(VecView<glm::vec2>(tmp_verts.data(), tmp_verts.size()));
+
+  auto cursor = tmp_verts.begin();
+  for (C2::PathD& path : paths) {
+    for (C2::PointD& p : path) {
+      p = v2_to_pd(*cursor);
+      ++cursor;
+    }
+  }
+
   return CrossSection(
-      shared_paths(C2::Union(warped, C2::FillRule::Positive, precision_)));
+      shared_paths(C2::Union(paths, C2::FillRule::Positive, precision_)));
 }
 
 /**
@@ -585,7 +618,28 @@ CrossSection CrossSection::Warp(
  * offseting operations are to be performed, which would compound the issue.
  */
 CrossSection CrossSection::Simplify(double epsilon) const {
-  auto ps = SimplifyPaths(GetPaths()->paths_, epsilon, false);
+  C2::PolyTreeD tree;
+  C2::BooleanOp(C2::ClipType::Union, C2::FillRule::Positive, GetPaths()->paths_,
+                C2::PathsD(), tree, precision_);
+
+  C2::PathsD polys;
+  flatten(&tree, polys, 0);
+
+  // Filter out contours less than epsilon wide.
+  C2::PathsD filtered;
+  for (C2::PathD poly : polys) {
+    auto area = C2::Area(poly);
+    Rect box;
+    for (auto vert : poly) {
+      box.Union(glm::vec2(vert.x, vert.y));
+    }
+    glm::vec2 size = box.Size();
+    if (glm::abs(area) > glm::max(size.x, size.y) * epsilon) {
+      filtered.push_back(poly);
+    }
+  }
+
+  auto ps = SimplifyPaths(filtered, epsilon, true);
   return CrossSection(shared_paths(ps));
 }
 

src/manifold/include/manifold.h

@@ -18,6 +18,7 @@
 
 #include "cross_section.h"
 #include "public.h"
+#include "vec_view.h"
 
 namespace manifold {
 
@@ -206,6 +207,7 @@ class Manifold {
   Manifold Transform(const glm::mat4x3&) const;
   Manifold Mirror(glm::vec3) const;
   Manifold Warp(std::function<void(glm::vec3&)>) const;
+  Manifold WarpBatch(std::function<void(VecView<glm::vec3>)>) const;
   Manifold SetProperties(
       int, std::function<void(float*, glm::vec3, const float*)>) const;
   Manifold CalculateCurvature(int gaussianIdx, int meanIdx) const;
@@ -234,9 +236,20 @@ class Manifold {
   Manifold TrimByPlane(glm::vec3 normal, float originOffset) const;
   ///@}
 
+  /** @name 2D from 3D
+   */
+  ///@{
+  CrossSection Slice(float height = 0) const;
+  CrossSection Project() const;
+  ///@}
+
+  /** @name Convex hull
+   */
+  ///@{
   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.

src/manifold/src/boolean3.cpp

@@ -92,6 +92,7 @@ struct CopyFaceEdges {
 
 SparseIndices Filter11(const Manifold::Impl &inP, const Manifold::Impl &inQ,
                        const SparseIndices &p1q2, const SparseIndices &p2q1) {
+  ZoneScoped;
   SparseIndices p1q1(3 * p1q2.size() + 3 * p2q1.size());
   for_each_n(autoPolicy(p1q2.size()), zip(countAt(0), countAt(0)), p1q2.size(),
              CopyFaceEdges<false>({p1q2, p1q1, inQ.halfedge_}));
@@ -273,6 +274,7 @@ std::tuple<Vec<int>, Vec<glm::vec4>> Shadow11(SparseIndices &p1q1,
                                               const Manifold::Impl &inP,
                                               const Manifold::Impl &inQ,
                                               float expandP) {
+  ZoneScoped;
   Vec<int> s11(p1q1.size());
   Vec<glm::vec4> xyzz11(p1q1.size());
 
@@ -367,6 +369,7 @@ std::tuple<Vec<int>, Vec<float>> Shadow02(const Manifold::Impl &inP,
                                           const Manifold::Impl &inQ,
                                           SparseIndices &p0q2, bool forward,
                                           float expandP) {
+  ZoneScoped;
   Vec<int> s02(p0q2.size());
   Vec<float> z02(p0q2.size());
 
@@ -476,6 +479,7 @@ std::tuple<Vec<int>, Vec<glm::vec3>> Intersect12(
     const SparseIndices &p0q2, const Vec<int> &s11, const SparseIndices &p1q1,
     const Vec<float> &z02, const Vec<glm::vec4> &xyzz11, SparseIndices &p1q2,
     bool forward) {
+  ZoneScoped;
   Vec<int> x12(p1q2.size());
   Vec<glm::vec3> v12(p1q2.size());
 
@@ -492,6 +496,7 @@ std::tuple<Vec<int>, Vec<glm::vec3>> Intersect12(
 
 Vec<int> Winding03(const Manifold::Impl &inP, Vec<int> &vertices, Vec<int> &s02,
                    bool reverse) {
+  ZoneScoped;
   // verts that are not shadowed (not in p0q2) have winding number zero.
   Vec<int> w03(inP.NumVert(), 0);
   // checking is slow, so just sort and reduce

src/manifold/src/boolean_result.cpp

@@ -101,6 +101,7 @@ std::tuple<Vec<int>, Vec<int>> SizeOutput(
     const Vec<int> &i03, const Vec<int> &i30, const Vec<int> &i12,
     const Vec<int> &i21, const SparseIndices &p1q2, const SparseIndices &p2q1,
     bool invertQ) {
+  ZoneScoped;
   Vec<int> sidesPerFacePQ(inP.NumTri() + inQ.NumTri(), 0);
   auto sidesPerFaceP = sidesPerFacePQ.view(0, inP.NumTri());
   auto sidesPerFaceQ = sidesPerFacePQ.view(inP.NumTri(), inQ.NumTri());
@@ -167,6 +168,7 @@ void AddNewEdgeVerts(
     concurrent_map<std::pair<int, int>, std::vector<EdgePos>> &edgesNew,
     const SparseIndices &p1q2, const Vec<int> &i12, const Vec<int> &v12R,
     const Vec<Halfedge> &halfedgeP, bool forward) {
+  ZoneScoped;
   // For each edge of P that intersects a face of Q (p1q2), add this vertex to
   // P's corresponding edge vector and to the two new edges, which are
   // intersections between the face of Q and the two faces of P attached to the
@@ -256,6 +258,7 @@ void AppendPartialEdges(Manifold::Impl &outR, Vec<char> &wholeHalfedgeP,
                         Vec<TriRef> &halfedgeRef, const Manifold::Impl &inP,
                         const Vec<int> &i03, const Vec<int> &vP2R,
                         const Vec<int>::IterC faceP2R, bool forward) {
+  ZoneScoped;
   // Each edge in the map is partially retained; for each of these, look up
   // their original verts and include them based on their winding number (i03),
   // while remapping them to the output using vP2R. Use the verts position
@@ -282,7 +285,6 @@ void AppendPartialEdges(Manifold::Impl &outR, Vec<char> &wholeHalfedgeP,
     }
 
     int inclusion = i03[vStart];
-    bool reversed = inclusion < 0;
     EdgePos edgePos = {vP2R[vStart],
                        glm::dot(outR.vertPos_[vP2R[vStart]], edgeVec),
                        inclusion > 0};
@@ -292,7 +294,6 @@ void AppendPartialEdges(Manifold::Impl &outR, Vec<char> &wholeHalfedgeP,
     }
 
     inclusion = i03[vEnd];
-    reversed |= inclusion < 0;
     edgePos = {vP2R[vEnd], glm::dot(outR.vertPos_[vP2R[vEnd]], edgeVec),
                inclusion < 0};
     for (int j = 0; j < glm::abs(inclusion); ++j) {
@@ -337,6 +338,7 @@ void AppendNewEdges(
     Manifold::Impl &outR, Vec<int> &facePtrR,
     concurrent_map<std::pair<int, int>, std::vector<EdgePos>> &edgesNew,
     Vec<TriRef> &halfedgeRef, const Vec<int> &facePQ2R, const int numFaceP) {
+  ZoneScoped;
   // Pair up each edge's verts and distribute to faces based on indices in key.
   Vec<Halfedge> &halfedgeR = outR.halfedge_;
   Vec<glm::vec3> &vertPosR = outR.vertPos_;
@@ -399,7 +401,6 @@ struct DuplicateHalfedges {
     if (!thrust::get<0>(in)) return;
     Halfedge halfedge = thrust::get<1>(in);
     if (!halfedge.IsForward()) return;
-    const int edgeP = thrust::get<2>(in);
 
     const int inclusion = i03[halfedge.startVert];
     if (inclusion == 0) return;
@@ -442,6 +443,7 @@ void AppendWholeEdges(Manifold::Impl &outR, Vec<int> &facePtrR,
                       const Vec<char> wholeHalfedgeP, const Vec<int> &i03,
                       const Vec<int> &vP2R, VecView<const int> faceP2R,
                       bool forward) {
+  ZoneScoped;
   for_each_n(ManifoldParams().deterministic ? ExecutionPolicy::Seq
                                             : autoPolicy(inP.halfedge_.size()),
              zip(wholeHalfedgeP.begin(), inP.halfedge_.begin(), countAt(0)),
@@ -463,9 +465,8 @@ struct MapTriRef {
   }
 };
 
-Vec<TriRef> UpdateReference(Manifold::Impl &outR, const Manifold::Impl &inP,
-                            const Manifold::Impl &inQ, bool invertQ) {
-  Vec<TriRef> refPQ = outR.meshRelation_.triRef;
+void UpdateReference(Manifold::Impl &outR, const Manifold::Impl &inP,
+                     const Manifold::Impl &inQ, bool invertQ) {
   const int offsetQ = Manifold::Impl::meshIDCounter_;
   for_each_n(
       autoPolicy(outR.NumTri()), outR.meshRelation_.triRef.begin(),
@@ -480,7 +481,6 @@ Vec<TriRef> UpdateReference(Manifold::Impl &outR, const Manifold::Impl &inP,
     outR.meshRelation_.meshIDtransform[pair.first + offsetQ].backSide ^=
         invertQ;
   }
-  return refPQ;
 }
 
 struct Barycentric {
@@ -514,8 +514,9 @@ struct Barycentric {
   }
 };
 
-void CreateProperties(Manifold::Impl &outR, const Vec<TriRef> &refPQ,
-                      const Manifold::Impl &inP, const Manifold::Impl &inQ) {
+void CreateProperties(Manifold::Impl &outR, const Manifold::Impl &inP,
+                      const Manifold::Impl &inQ) {
+  ZoneScoped;
   const int numPropP = inP.NumProp();
   const int numPropQ = inQ.NumProp();
   const int numProp = glm::max(numPropP, numPropQ);
@@ -526,7 +527,8 @@ void CreateProperties(Manifold::Impl &outR, const Vec<TriRef> &refPQ,
   outR.meshRelation_.triProperties.resize(numTri);
 
   Vec<glm::vec3> bary(outR.halfedge_.size());
-  for_each_n(autoPolicy(numTri), zip(countAt(0), refPQ.cbegin()), numTri,
+  for_each_n(autoPolicy(numTri),
+             zip(countAt(0), outR.meshRelation_.triRef.cbegin()), numTri,
              Barycentric({bary, inP.vertPos_, inQ.vertPos_, outR.vertPos_,
                           inP.halfedge_, inQ.halfedge_, outR.halfedge_,
                           outR.precision_}));
@@ -534,6 +536,9 @@ void CreateProperties(Manifold::Impl &outR, const Vec<TriRef> &refPQ,
   using Entry = std::pair<glm::ivec3, int>;
   int idMissProp = outR.NumVert();
   std::vector<std::vector<Entry>> propIdx(outR.NumVert() + 1);
+  std::vector<int> propMissIdx[2];
+  propMissIdx[0].resize(inQ.NumPropVert(), -1);
+  propMissIdx[1].resize(inP.NumPropVert(), -1);
   outR.meshRelation_.properties.reserve(outR.NumVert() * numProp);
   int idx = 0;
 
@@ -541,14 +546,14 @@ void CreateProperties(Manifold::Impl &outR, const Vec<TriRef> &refPQ,
     // Skip collapsed triangles
     if (outR.halfedge_[3 * tri].startVert < 0) continue;
 
-    const int triPQ = refPQ[tri].tri;
-    const bool PQ = refPQ[tri].meshID == 0;
+    const TriRef ref = outR.meshRelation_.triRef[tri];
+    const bool PQ = ref.meshID == 0;
     const int oldNumProp = PQ ? numPropP : numPropQ;
     const auto &properties =
         PQ ? inP.meshRelation_.properties : inQ.meshRelation_.properties;
     const glm::ivec3 &triProp = oldNumProp == 0 ? glm::ivec3(-1)
-                                : PQ ? inP.meshRelation_.triProperties[triPQ]
-                                     : inQ.meshRelation_.triProperties[triPQ];
+                                : PQ ? inP.meshRelation_.triProperties[ref.tri]
+                                     : inQ.meshRelation_.triProperties[ref.tri];
 
     for (const int i : {0, 1, 2}) {
       const int vert = outR.halfedge_[3 * tri + i].startVert;
@@ -556,8 +561,7 @@ void CreateProperties(Manifold::Impl &outR, const Vec<TriRef> &refPQ,
 
       glm::ivec4 key(PQ, idMissProp, -1, -1);
       if (oldNumProp > 0) {
-        key[1] = vert;
-        int edge = -1;
+        int edge = -2;
         for (const int j : {0, 1, 2}) {
           if (uvw[j] == 1) {
             // On a retained vert, the propVert must also match
@@ -571,24 +575,37 @@ void CreateProperties(Manifold::Impl &outR, const Vec<TriRef> &refPQ,
           // On an edge, both propVerts must match
           const int p0 = triProp[Next3(edge)];
           const int p1 = triProp[Prev3(edge)];
+          key[1] = vert;
           key[2] = glm::min(p0, p1);
           key[3] = glm::max(p0, p1);
+        } else if (edge == -2) {
+          key[1] = vert;
         }
       }
 
-      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 (key.y == idMissProp && key.z >= 0) {
+        // only key.x/key.z matters
+        auto &entry = propMissIdx[key.x][key.z];
+        if (entry >= 0) {
+          outR.meshRelation_.triProperties[tri][i] = entry;
+          continue;
         }
+        entry = idx;
+      } else {
+        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 (bFound) continue;
+        bin.push_back(std::make_pair(glm::ivec3(key.x, key.z, key.w), 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) {
         if (p < oldNumProp) {
           glm::vec3 oldProps;
@@ -724,7 +741,6 @@ Manifold::Impl Boolean3::Result(OpType op) const {
   std::tie(faceEdge, facePQ2R) =
       SizeOutput(outR, inP_, inQ_, i03, i30, i12, i21, p1q2_, p2q1_, invertQ);
 
-  const int numFaceR = faceEdge.size() - 1;
   // This gets incremented for each halfedge that's added to a face so that the
   // next one knows where to slot in.
   Vec<int> facePtrR = faceEdge;
@@ -770,12 +786,12 @@ Manifold::Impl Boolean3::Result(OpType op) const {
   if (ManifoldParams().intermediateChecks)
     ASSERT(outR.IsManifold(), logicErr, "triangulated mesh is not manifold!");
 
-  Vec<TriRef> refPQ = UpdateReference(outR, inP_, inQ_, invertQ);
+  CreateProperties(outR, inP_, inQ_);
+
+  UpdateReference(outR, inP_, inQ_, invertQ);
 
   outR.SimplifyTopology();
 
-  CreateProperties(outR, refPQ, inP_, inQ_);
-
   if (ManifoldParams().intermediateChecks)
     ASSERT(outR.Is2Manifold(), logicErr, "simplified mesh is not 2-manifold!");
 

src/manifold/src/constructors.cpp

@@ -220,6 +220,9 @@ Manifold Manifold::Sphere(float radius, int circularSegments) {
 /**
  * Constructs a manifold from a set of polygons by extruding them along the
  * Z-axis.
+ * Note that high twistDegrees with small nDivisions may cause
+ * self-intersection. This is not checked here and it is up to the user to
+ * choose the correct parameters.
  *
  * @param crossSection A set of non-overlapping polygons to extrude.
  * @param height Z-extent of extrusion.
@@ -230,11 +233,13 @@ Manifold Manifold::Sphere(float radius, int circularSegments) {
  * bottom, interpolated linearly for the divisions in between.
  * @param scaleTop Amount to scale the top (independently in X and Y). If the
  * scale is {0, 0}, a pure cone is formed with only a single vertex at the top.
+ * Note that scale is applied after twist.
  * Default {1, 1}.
  */
 Manifold Manifold::Extrude(const CrossSection& crossSection, float height,
                            int nDivisions, float twistDegrees,
                            glm::vec2 scaleTop) {
+  ZoneScoped;
   auto polygons = crossSection.ToPolygons();
   if (polygons.size() == 0 || height <= 0.0f) {
     return Invalid();
@@ -264,9 +269,9 @@ Manifold Manifold::Extrude(const CrossSection& crossSection, float height,
   for (int i = 1; i < nDivisions + 1; ++i) {
     float alpha = i / float(nDivisions);
     float phi = alpha * twistDegrees;
-    glm::mat2 transform(cosd(phi), sind(phi), -sind(phi), cosd(phi));
     glm::vec2 scale = glm::mix(glm::vec2(1.0f), scaleTop, alpha);
-    transform = transform * glm::mat2(scale.x, 0.0f, 0.0f, scale.y);
+    glm::mat2 rotation(cosd(phi), sind(phi), -sind(phi), cosd(phi));
+    glm::mat2 transform = glm::mat2(scale.x, 0.0f, 0.0f, scale.y) * rotation;
     int j = 0;
     int idx = 0;
     for (const auto& poly : polygons) {
@@ -320,6 +325,7 @@ Manifold Manifold::Extrude(const CrossSection& crossSection, float height,
  */
 Manifold Manifold::Revolve(const CrossSection& crossSection,
                            int circularSegments, float revolveDegrees) {
+  ZoneScoped;
   Polygons polygons = crossSection.ToPolygons();
 
   if (polygons.size() == 0) {
@@ -327,8 +333,9 @@ Manifold Manifold::Revolve(const CrossSection& crossSection,
   }
 
   const Rect bounds = crossSection.Bounds();
+  const float radius = bounds.max.x;
 
-  if (bounds.max.x <= 0) {
+  if (radius <= 0) {
     return Invalid();
   } else if (bounds.min.x < 0) {
     // Take the x>=0 slice.
@@ -340,21 +347,15 @@ Manifold Manifold::Revolve(const CrossSection& crossSection,
     polygons = (crossSection ^ posBoundingBox).ToPolygons();
   }
 
-  float radius = 0.0f;
-  for (const auto& poly : polygons) {
-    for (const auto& vert : poly) {
-      radius = fmax(radius, vert.x);
-    }
-  }
-
   if (revolveDegrees > 360.0f) {
     revolveDegrees = 360.0f;
   }
   const bool isFullRevolution = revolveDegrees == 360.0f;
 
-  const int nDivisions = circularSegments > 2
-                             ? circularSegments
-                             : Quality::GetCircularSegments(radius);
+  const int nDivisions =
+      circularSegments > 2
+          ? circularSegments
+          : Quality::GetCircularSegments(radius) * revolveDegrees / 360;
 
   auto pImpl_ = std::make_shared<Impl>();
   auto& vertPos = pImpl_->vertPos_;
@@ -384,14 +385,14 @@ Manifold Manifold::Revolve(const CrossSection& crossSection,
 
       if (!isFullRevolution) startPoses.push_back(startPosIndex);
 
+      const glm::vec2 currPolyVertex = poly[polyVert];
+      const glm::vec2 prevPolyVertex =
+          poly[polyVert == 0 ? poly.size() - 1 : polyVert - 1];
+
       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];
+          (polyVert == 0 ? nRevolveAxisVerts + (nSlices * nPosVerts) : 0) +
+          (prevPolyVertex.x == 0.0 ? -1 : -nSlices);
 
       for (int slice = 0; slice < nSlices; ++slice) {
         const float phi = slice * dPhi;
@@ -464,6 +465,7 @@ Manifold Manifold::Compose(const std::vector<Manifold>& manifolds) {
  * containing a copy of the original. It is the inverse operation of Compose().
  */
 std::vector<Manifold> Manifold::Decompose() const {
+  ZoneScoped;
   UnionFind<> uf(NumVert());
   // Graph graph;
   auto pImpl_ = GetCsgLeafNode().GetImpl();

src/manifold/src/csg_tree.cpp

@@ -172,6 +172,7 @@ CsgNodeType CsgLeafNode::GetNodeType() const { return CsgNodeType::Leaf; }
  */
 Manifold::Impl CsgLeafNode::Compose(
     const std::vector<std::shared_ptr<CsgLeafNode>> &nodes) {
+  ZoneScoped;
   float precision = -1;
   int numVert = 0;
   int numEdge = 0;
@@ -454,6 +455,7 @@ std::shared_ptr<CsgLeafNode> CsgOpNode::ToLeafNode() const {
 std::shared_ptr<Manifold::Impl> CsgOpNode::BatchBoolean(
     OpType operation,
     std::vector<std::shared_ptr<const Manifold::Impl>> &results) {
+  ZoneScoped;
   auto getImplPtr = GetImplPtr();
   ASSERT(operation != OpType::Subtract, logicErr,
          "BatchBoolean doesn't support Difference.");
@@ -529,6 +531,7 @@ std::shared_ptr<Manifold::Impl> CsgOpNode::BatchBoolean(
  * `BatchBoolean` instead of using `Compose` for non-intersecting manifolds.
  */
 void CsgOpNode::BatchUnion() const {
+  ZoneScoped;
   // INVARIANT: children_ is a vector of leaf nodes
   // this kMaxUnionSize is a heuristic to avoid the pairwise disjoint check
   // with O(n^2) complexity to take too long.

src/manifold/src/edge_op.cpp

@@ -154,20 +154,23 @@ void Manifold::Impl::SimplifyTopology() {
   // verts. They must be removed before edge collapse.
   SplitPinchedVerts();
 
-  Vec<SortEntry> entries;
-  entries.reserve(nbEdges / 2);
-  for (int i = 0; i < nbEdges; ++i) {
-    if (halfedge_[i].IsForward()) {
-      entries.push_back({halfedge_[i].startVert, halfedge_[i].endVert, i});
+  {
+    ZoneScopedN("DedupeEdge");
+    Vec<SortEntry> entries;
+    entries.reserve(nbEdges / 2);
+    for (int i = 0; i < nbEdges; ++i) {
+      if (halfedge_[i].IsForward()) {
+        entries.push_back({halfedge_[i].startVert, halfedge_[i].endVert, i});
+      }
     }
-  }
 
-  stable_sort(policy, entries.begin(), entries.end());
-  for (int i = 0; i < entries.size() - 1; ++i) {
-    if (entries[i].start == entries[i + 1].start &&
-        entries[i].end == entries[i + 1].end) {
-      DedupeEdge(entries[i].index);
-      numFlagged++;
+    stable_sort(policy, entries.begin(), entries.end());
+    for (int i = 0; i < entries.size() - 1; ++i) {
+      if (entries[i].start == entries[i + 1].start &&
+          entries[i].end == entries[i + 1].end) {
+        DedupeEdge(entries[i].index);
+        numFlagged++;
+      }
     }
   }
 
@@ -178,9 +181,14 @@ void Manifold::Impl::SimplifyTopology() {
   }
 #endif
 
+  if (!ManifoldParams().cleanupTriangles) {
+    return;
+  }
+
   std::vector<int> scratchBuffer;
   scratchBuffer.reserve(10);
   {
+    ZoneScopedN("CollapseShortEdge");
     numFlagged = 0;
     ShortEdge se{halfedge_, vertPos_, precision_};
     for_each_n(policy, countAt(0), nbEdges, [&](int i) { bflags[i] = se(i); });
@@ -201,6 +209,7 @@ void Manifold::Impl::SimplifyTopology() {
 #endif
 
   {
+    ZoneScopedN("CollapseFlaggedEdge");
     numFlagged = 0;
     FlagEdge se{halfedge_, meshRelation_.triRef};
     for_each_n(policy, countAt(0), nbEdges, [&](int i) { bflags[i] = se(i); });
@@ -221,6 +230,7 @@ void Manifold::Impl::SimplifyTopology() {
 #endif
 
   {
+    ZoneScopedN("RecursiveEdgeSwap");
     numFlagged = 0;
     SwappableEdge se{halfedge_, vertPos_, faceNormal_, precision_};
     for_each_n(policy, countAt(0), nbEdges, [&](int i) { bflags[i] = se(i); });
@@ -501,25 +511,20 @@ void Manifold::Impl::CollapseEdge(const int edge, std::vector<int>& edges) {
   const int tri1 = toRemove.pairedHalfedge / 3;
   const int triVert0 = (edge + 1) % 3;
   const int triVert1 = toRemove.pairedHalfedge % 3;
-  const int prop0 = triProp.size() > 0 ? triProp[tri0][edge % 3] : -1;
-  const int prop1 = triProp.size() > 0
-                        ? triProp[tri1][(toRemove.pairedHalfedge + 1) % 3]
-                        : -1;
   current = start;
   while (current != tri0edge[2]) {
     current = NextHalfedge(current);
 
-    // Update the shifted triangles to the vertBary of endVert
-    const int tri = current / 3;
-    const int vIdx = current - 3 * tri;
-
-    if (!shortEdge) {
-      if (triProp.size() > 0) {
-        if (triProp[tri][vIdx] == prop0) {
-          triProp[tri][vIdx] = triProp[tri0][triVert0];
-        } else if (triProp[tri][vIdx] == prop1) {
-          triProp[tri][vIdx] = triProp[tri1][triVert1];
-        }
+    if (triProp.size() > 0) {
+      // Update the shifted triangles to the vertBary of endVert
+      const int tri = current / 3;
+      const int vIdx = current - 3 * tri;
+      if (triRef[tri].meshID == triRef[tri0].meshID &&
+          triRef[tri].tri == triRef[tri0].tri) {
+        triProp[tri][vIdx] = triProp[tri0][triVert0];
+      } else if (triRef[tri].meshID == triRef[tri1].meshID &&
+                 triRef[tri].tri == triRef[tri1].tri) {
+        triProp[tri][vIdx] = triProp[tri1][triVert1];
       }
     }
 
@@ -656,6 +661,7 @@ void Manifold::Impl::RecursiveEdgeSwap(const int edge, int& tag,
 }
 
 void Manifold::Impl::SplitPinchedVerts() {
+  ZoneScoped;
   std::vector<bool> vertProcessed(NumVert(), false);
   std::vector<bool> halfedgeProcessed(halfedge_.size(), false);
   for (int i = 0; i < halfedge_.size(); ++i) {

src/manifold/src/face_op.cpp

@@ -17,7 +17,7 @@
 #define TBB_PREVIEW_CONCURRENT_ORDERED_CONTAINERS 1
 #include <tbb/concurrent_map.h>
 #endif
-#include <map>
+#include <unordered_set>
 
 #include "impl.h"
 #include "polygon.h"
@@ -39,6 +39,7 @@ using AddTriangle = std::function<void(int, glm::ivec3, glm::vec3, TriRef)>;
  */
 void Manifold::Impl::Face2Tri(const Vec<int>& faceEdge,
                               const Vec<TriRef>& halfedgeRef) {
+  ZoneScoped;
   Vec<glm::ivec3> triVerts;
   Vec<glm::vec3> triNormal;
   Vec<TriRef>& triRef = meshRelation_.triRef;
@@ -126,7 +127,9 @@ void Manifold::Impl::Face2Tri(const Vec<int>& faceEdge,
   auto generalTriangulation = [&](int face) {
     const glm::vec3 normal = faceNormal_[face];
     const glm::mat3x2 projection = GetAxisAlignedProjection(normal);
-    const PolygonsIdx polys = Face2Polygons(face, projection, faceEdge);
+    const PolygonsIdx polys =
+        Face2Polygons(halfedge_.cbegin() + faceEdge[face],
+                      halfedge_.cbegin() + faceEdge[face + 1], projection);
     return TriangulateIdx(polys, precision_);
   };
 #if MANIFOLD_PAR == 'T' && __has_include(<tbb/tbb.h>)
@@ -191,20 +194,16 @@ void Manifold::Impl::Face2Tri(const Vec<int>& faceEdge,
 }
 
 /**
- * For the input face index, return a set of 2D polygons formed by the input
- * projection of the vertices.
+ * Returns a set of 2D polygons formed by the input projection of the vertices
+ * of the list of Halfedges, which must be an even-manifold, meaning each vert
+ * must be referenced the same number of times as a startVert and endVert.
  */
-PolygonsIdx Manifold::Impl::Face2Polygons(int face, glm::mat3x2 projection,
-                                          const Vec<int>& faceEdge) const {
-  const int firstEdge = faceEdge[face];
-  const int lastEdge = faceEdge[face + 1];
-
-  std::map<int, int> vert_edge;
-  for (int edge = firstEdge; edge < lastEdge; ++edge) {
-    const bool inserted =
-        vert_edge.emplace(std::make_pair(halfedge_[edge].startVert, edge))
-            .second;
-    ASSERT(inserted, topologyErr, "face has duplicate vertices.");
+PolygonsIdx Manifold::Impl::Face2Polygons(VecView<Halfedge>::IterC start,
+                                          VecView<Halfedge>::IterC end,
+                                          glm::mat3x2 projection) const {
+  std::multimap<int, int> vert_edge;
+  for (auto edge = start; edge != end; ++edge) {
+    vert_edge.emplace(std::make_pair(edge->startVert, edge - start));
   }
 
   PolygonsIdx polys;
@@ -217,13 +216,104 @@ PolygonsIdx Manifold::Impl::Face2Polygons(int face, glm::mat3x2 projection,
       thisEdge = startEdge;
       polys.push_back({});
     }
-    int vert = halfedge_[thisEdge].startVert;
+    int vert = (start + thisEdge)->startVert;
     polys.back().push_back({projection * vertPos_[vert], vert});
-    const auto result = vert_edge.find(halfedge_[thisEdge].endVert);
+    const auto result = vert_edge.find((start + thisEdge)->endVert);
     ASSERT(result != vert_edge.end(), topologyErr, "non-manifold edge");
     thisEdge = result->second;
     vert_edge.erase(result);
   }
   return polys;
 }
+
+CrossSection Manifold::Impl::Slice(float height) const {
+  Box plane = bBox_;
+  plane.min.z = plane.max.z = height;
+  Vec<Box> query;
+  query.push_back(plane);
+  const SparseIndices collisions =
+      collider_.Collisions<false, false>(query.cview());
+
+  std::unordered_set<int> tris;
+  for (int i = 0; i < collisions.size(); ++i) {
+    const int tri = collisions.Get(i, 1);
+    float min = std::numeric_limits<float>::infinity();
+    float max = -std::numeric_limits<float>::infinity();
+    for (const int j : {0, 1, 2}) {
+      const float z = vertPos_[halfedge_[3 * tri + j].startVert].z;
+      min = glm::min(min, z);
+      max = glm::max(max, z);
+    }
+
+    if (min <= height && max > height) {
+      tris.insert(tri);
+    }
+  }
+
+  Polygons polys;
+  while (!tris.empty()) {
+    const int startTri = *tris.begin();
+    SimplePolygon poly;
+
+    int k = 0;
+    for (const int j : {0, 1, 2}) {
+      if (vertPos_[halfedge_[3 * startTri + j].startVert].z > height &&
+          vertPos_[halfedge_[3 * startTri + Next3(j)].startVert].z <= height) {
+        k = Next3(j);
+        break;
+      }
+    }
+
+    int tri = startTri;
+    do {
+      tris.erase(tris.find(tri));
+      if (vertPos_[halfedge_[3 * tri + k].endVert].z <= height) {
+        k = Next3(k);
+      }
+
+      Halfedge up = halfedge_[3 * tri + k];
+      const glm::vec3 below = vertPos_[up.startVert];
+      const glm::vec3 above = vertPos_[up.endVert];
+      const float a = (height - below.z) / (above.z - below.z);
+      poly.push_back(glm::vec2(glm::mix(below, above, a)));
+
+      const int pair = up.pairedHalfedge;
+      tri = pair / 3;
+      k = Next3(pair % 3);
+    } while (tri != startTri);
+
+    polys.push_back(poly);
+  }
+
+  return CrossSection(polys);
+}
+
+CrossSection Manifold::Impl::Project() const {
+  const glm::mat3x2 projection = GetAxisAlignedProjection({0, 0, 1});
+  auto policy = autoPolicy(halfedge_.size());
+
+  Vec<Halfedge> cusps(NumEdge());
+  cusps.resize(copy_if<decltype(cusps.begin())>(
+                   policy, halfedge_.cbegin(), halfedge_.cend(), cusps.begin(),
+                   [&](Halfedge edge) {
+                     return faceNormal_[edge.face].z >= 0 &&
+                            faceNormal_[halfedge_[edge.pairedHalfedge].face].z <
+                                0;
+                   }) -
+               cusps.begin());
+
+  PolygonsIdx polysIndexed =
+      Face2Polygons(cusps.cbegin(), cusps.cend(), projection);
+
+  Polygons polys;
+  for (const auto& poly : polysIndexed) {
+    SimplePolygon simple;
+    for (const PolyVert& polyVert : poly) {
+      simple.push_back(polyVert.pos);
+    }
+    polys.push_back(simple);
+  }
+
+  return CrossSection(polys).Simplify(precision_);
+}
 }  // namespace manifold

src/manifold/src/impl.cpp

@@ -22,6 +22,7 @@
 #include "hashtable.h"
 #include "mesh_fixes.h"
 #include "par.h"
+#include "svd.h"
 
 namespace {
 using namespace manifold;
@@ -320,6 +321,7 @@ int GetLabels(std::vector<int>& components,
 
 void DedupePropVerts(manifold::Vec<glm::ivec3>& triProp,
                      const Vec<thrust::pair<int, int>>& vert2vert) {
+  ZoneScoped;
   std::vector<int> vertLabels;
   const int numLabels = GetLabels(vertLabels, vert2vert, vert2vert.size());
 
@@ -581,6 +583,7 @@ Manifold::Impl::Impl(Shape shape) {
 }
 
 void Manifold::Impl::RemoveUnreferencedVerts(Vec<glm::ivec3>& triVerts) {
+  ZoneScoped;
   Vec<int> vertOld2New(NumVert() + 1, 0);
   auto policy = autoPolicy(NumVert());
   for_each(policy, triVerts.cbegin(), triVerts.cend(),
@@ -613,6 +616,7 @@ void Manifold::Impl::InitializeOriginal() {
 }
 
 void Manifold::Impl::CreateFaces(const std::vector<float>& propertyTolerance) {
+  ZoneScoped;
   Vec<float> propertyToleranceD =
       propertyTolerance.empty() ? Vec<float>(meshRelation_.numProp, kTolerance)
                                 : propertyTolerance;
@@ -661,6 +665,7 @@ void Manifold::Impl::CreateFaces(const std::vector<float>& propertyTolerance) {
  * Create the halfedge_ data structure from an input triVerts array like Mesh.
  */
 void Manifold::Impl::CreateHalfedges(const Vec<glm::ivec3>& triVerts) {
+  ZoneScoped;
   const int numTri = triVerts.size();
   const int numHalfedge = 3 * numTri;
   // drop the old value first to avoid copy
@@ -714,7 +719,14 @@ void Manifold::Impl::MarkFailure(Error status) {
 }
 
 void Manifold::Impl::Warp(std::function<void(glm::vec3&)> warpFunc) {
-  thrust::for_each_n(thrust::host, vertPos_.begin(), NumVert(), warpFunc);
+  WarpBatch([&warpFunc](VecView<glm::vec3> vecs) {
+    thrust::for_each(thrust::host, vecs.begin(), vecs.end(), warpFunc);
+  });
+}
+
+void Manifold::Impl::WarpBatch(
+    std::function<void(VecView<glm::vec3>)> warpFunc) {
+  warpFunc(vertPos_.view());
   CalculateBBox();
   if (!IsFinite()) {
     MarkFailure(Error::NonFiniteVertex);
@@ -725,11 +737,11 @@ void Manifold::Impl::Warp(std::function<void(glm::vec3&)> warpFunc) {
   CalculateNormals();
   SetPrecision();
   CreateFaces();
-  SimplifyTopology();
   Finish();
 }
 
 Manifold::Impl Manifold::Impl::Transform(const glm::mat4x3& transform_) const {
+  ZoneScoped;
   if (transform_ == glm::mat4x3(1.0f)) return *this;
   auto policy = autoPolicy(NumVert());
   Impl result;
@@ -776,11 +788,8 @@ Manifold::Impl Manifold::Impl::Transform(const glm::mat4x3& transform_) const {
   if (!result.collider_.Transform(transform_)) result.Update();
 
   result.CalculateBBox();
-  float scale = 0;
-  for (int i : {0, 1, 2})
-    scale =
-        glm::max(scale, transform_[0][i] + transform_[1][i] + transform_[2][i]);
-  result.precision_ *= scale;
+  // Scale the precision by the norm of the 3x3 portion of the transform.
+  result.precision_ *= SpectralNorm(glm::mat3(transform_));
   // Maximum of inherited precision loss and translational precision loss.
   result.SetPrecision(result.precision_);
   return result;
@@ -807,6 +816,7 @@ void Manifold::Impl::SetPrecision(float minPrecision) {
  * recalculation.
  */
 void Manifold::Impl::CalculateNormals() {
+  ZoneScoped;
   vertNormal_.resize(NumVert());
   auto policy = autoPolicy(NumTri());
   fill(policy, vertNormal_.begin(), vertNormal_.end(), glm::vec3(0));
@@ -849,6 +859,7 @@ void Manifold::Impl::IncrementMeshIDs() {
  */
 SparseIndices Manifold::Impl::EdgeCollisions(const Impl& Q,
                                              bool inverted) const {
+  ZoneScoped;
   Vec<TmpEdge> edges = CreateTmpEdges(Q.halfedge_);
   const int numEdge = edges.size();
   Vec<Box> QedgeBB(numEdge);
@@ -877,6 +888,7 @@ SparseIndices Manifold::Impl::EdgeCollisions(const Impl& Q,
  */
 SparseIndices Manifold::Impl::VertexCollisionsZ(
     VecView<const glm::vec3> vertsIn, bool inverted) const {
+  ZoneScoped;
   if (inverted)
     return collider_.Collisions<false, true>(vertsIn);
   else

src/manifold/src/impl.h

@@ -76,6 +76,7 @@ struct Manifold::Impl {
   void Update();
   void MarkFailure(Error status);
   void Warp(std::function<void(glm::vec3&)> warpFunc);
+  void WarpBatch(std::function<void(VecView<glm::vec3>)> warpFunc);
   Impl Transform(const glm::mat4x3& transform) const;
   SparseIndices EdgeCollisions(const Impl& B, bool inverted = false) const;
   SparseIndices VertexCollisionsZ(VecView<const glm::vec3> vertsIn,
@@ -115,8 +116,11 @@ struct Manifold::Impl {
 
   // face_op.cu
   void Face2Tri(const Vec<int>& faceEdge, const Vec<TriRef>& halfedgeRef);
-  PolygonsIdx Face2Polygons(int face, glm::mat3x2 projection,
-                            const Vec<int>& faceEdge) const;
+  PolygonsIdx Face2Polygons(VecView<Halfedge>::IterC start,
+                            VecView<Halfedge>::IterC end,
+                            glm::mat3x2 projection) const;
+  CrossSection Slice(float height) const;
+  CrossSection Project() const;
 
   // edge_op.cu
   void SimplifyTopology();

src/manifold/src/manifold.cpp

@@ -155,6 +155,7 @@ Manifold::Manifold(const Mesh& mesh) {
  * saving or other operations outside of the context of this library.
  */
 Mesh Manifold::GetMesh() const {
+  ZoneScoped;
   const Impl& impl = *GetCsgLeafNode().GetImpl();
 
   Mesh result;
@@ -189,6 +190,7 @@ Mesh Manifold::GetMesh() const {
  * MeshGLs must use the same channels for their normals.
  */
 MeshGL Manifold::GetMeshGL(glm::ivec3 normalIdx) const {
+  ZoneScoped;
   const Impl& impl = *GetCsgLeafNode().GetImpl();
 
   const int numProp = NumProp();
@@ -289,7 +291,7 @@ MeshGL Manifold::GetMeshGL(glm::ivec3 normalIdx) const {
   // Duplicate verts with different props
   std::vector<int> vert2idx(impl.NumVert(), -1);
   std::vector<std::vector<glm::ivec2>> vertPropPair(impl.NumVert());
-  out.vertProperties.reserve(numVert * out.numProp);
+  out.vertProperties.reserve(numVert * static_cast<size_t>(out.numProp));
 
   for (int run = 0; run < out.runOriginalID.size(); ++run) {
     for (int tri = out.runIndex[run] / 3; tri < out.runIndex[run + 1] / 3;
@@ -484,7 +486,7 @@ int Manifold::NumOverlaps(const Manifold& other) const {
 
   overlaps = other.GetCsgLeafNode().GetImpl()->EdgeCollisions(
       *GetCsgLeafNode().GetImpl());
-  return num_overlaps += overlaps.size();
+  return num_overlaps + overlaps.size();
 }
 
 /**
@@ -567,6 +569,20 @@ Manifold Manifold::Warp(std::function<void(glm::vec3&)> warpFunc) const {
   return Manifold(std::make_shared<CsgLeafNode>(pImpl));
 }
 
+/**
+ * Same as Manifold::Warp but calls warpFunc with with
+ * a VecView which is roughly equivalent to std::span
+ * pointing to all vec3 elements to be modified in-place
+ *
+ * @param warpFunc A function that modifies multiple vertex positions.
+ */
+Manifold Manifold::WarpBatch(
+    std::function<void(VecView<glm::vec3>)> warpFunc) const {
+  auto pImpl = std::make_shared<Impl>(*GetCsgLeafNode().GetImpl());
+  pImpl->WarpBatch(warpFunc);
+  return Manifold(std::make_shared<CsgLeafNode>(pImpl));
+}
+
 /**
  * Create a new copy of this manifold with updated vertex properties by
  * supplying a function that takes the existing position and properties as
@@ -613,7 +629,6 @@ Manifold Manifold::SetProperties(
 
   pImpl->meshRelation_.numProp = numProp;
   pImpl->CreateFaces();
-  pImpl->SimplifyTopology();
   pImpl->Finish();
   return Manifold(std::make_shared<CsgLeafNode>(pImpl));
 }
@@ -777,6 +792,24 @@ Manifold Manifold::TrimByPlane(glm::vec3 normal, float originOffset) const {
   return *this ^ Halfspace(BoundingBox(), normal, originOffset);
 }
 
+/**
+ * Returns the cross section of this object parallel to the X-Y plane at the
+ * specified Z height, defaulting to zero. Using a height equal to the bottom of
+ * the bounding box will return the bottom faces, while using a height equal to
+ * the top of the bounding box will return empty.
+ */
+CrossSection Manifold::Slice(float height) const {
+  return GetCsgLeafNode().GetImpl()->Slice(height);
+}
+
+/**
+ * Returns a cross section representing the projected outline of this object
+ * onto the X-Y plane.
+ */
+CrossSection Manifold::Project() const {
+  return GetCsgLeafNode().GetImpl()->Project();
+}
+
 ExecutionParams& ManifoldParams() { return manifoldParams; }
 
 /**
@@ -787,6 +820,7 @@ ExecutionParams& ManifoldParams() { return manifoldParams; }
  * hull.
  */
 Manifold Manifold::Hull(const std::vector<glm::vec3>& pts) {
+  ZoneScoped;
   const int numVert = pts.size();
   if (numVert < 4) return Manifold();
 

src/manifold/src/properties.cpp

@@ -316,6 +316,7 @@ int Manifold::Impl::NumDegenerateTris() const {
 }
 
 Properties Manifold::Impl::GetProperties() const {
+  ZoneScoped;
   if (IsEmpty()) return {0, 0};
   // Kahan summation
   float area = 0;
@@ -339,6 +340,7 @@ Properties Manifold::Impl::GetProperties() const {
 }
 
 void Manifold::Impl::CalculateCurvature(int gaussianIdx, int meanIdx) {
+  ZoneScoped;
   if (IsEmpty()) return;
   if (gaussianIdx < 0 && meanIdx < 0) return;
   Vec<float> vertMeanCurvature(NumVert(), 0);
@@ -370,7 +372,6 @@ void Manifold::Impl::CalculateCurvature(int gaussianIdx, int meanIdx) {
                         numProp, gaussianIdx, meanIdx}));
 
   CreateFaces();
-  SimplifyTopology();
   Finish();
 }
 

src/manifold/src/smoothing.cpp

@@ -196,10 +196,13 @@ struct SmoothBezier {
     const glm::vec3 edgeNormal =
         (triNormal[edge.face] + triNormal[halfedge[edge.pairedHalfedge].face]) /
         2.0f;
-    glm::vec3 dir = glm::normalize(glm::cross(glm::cross(edgeNormal, edgeVec),
-                                              vertNormal[edge.startVert]));
+    glm::vec3 dir = SafeNormalize(glm::cross(glm::cross(edgeNormal, edgeVec),
+                                             vertNormal[edge.startVert]));
 
-    const float weight = glm::abs(glm::dot(dir, glm::normalize(edgeVec)));
+    float weight = glm::abs(glm::dot(dir, SafeNormalize(edgeVec)));
+    if (weight == 0) {
+      weight = 1;
+    }
     // Quadratic weighted bezier for circular interpolation
     const glm::vec4 bz2 =
         weight *
@@ -320,6 +323,7 @@ namespace manifold {
  */
 void Manifold::Impl::CreateTangents(
     const std::vector<Smoothness>& sharpenedEdges) {
+  ZoneScoped;
   const int numHalfedge = halfedge_.size();
   halfedgeTangent_.resize(numHalfedge);
 
@@ -423,6 +427,7 @@ void Manifold::Impl::CreateTangents(
  * refinement (smoothing).
  */
 Vec<Barycentric> Manifold::Impl::Subdivide(int n) {
+  ZoneScoped;
   if (n < 2) return Vec<Barycentric>();
   faceNormal_.resize(0);
   vertNormal_.resize(0);

src/manifold/src/sort.cpp

@@ -300,6 +300,7 @@ void Manifold::Impl::Finish() {
  * Sorts the vertices according to their Morton code.
  */
 void Manifold::Impl::SortVerts() {
+  ZoneScoped;
   const int numVert = NumVert();
   Vec<uint32_t> vertMorton(numVert);
   auto policy = autoPolicy(numVert);
@@ -338,6 +339,7 @@ void Manifold::Impl::SortVerts() {
  * also given.
  */
 void Manifold::Impl::ReindexVerts(const Vec<int>& vertNew2Old, int oldNumVert) {
+  ZoneScoped;
   Vec<int> vertOld2New(oldNumVert);
   scatter(autoPolicy(oldNumVert), countAt(0), countAt(NumVert()),
           vertNew2Old.begin(), vertOld2New.begin());
@@ -349,6 +351,7 @@ void Manifold::Impl::ReindexVerts(const Vec<int>& vertNew2Old, int oldNumVert) {
  * Removes unreferenced property verts and reindexes triProperties.
  */
 void Manifold::Impl::CompactProps() {
+  ZoneScoped;
   if (meshRelation_.numProp == 0) return;
 
   const int numVerts = meshRelation_.properties.size() / meshRelation_.numProp;
@@ -377,6 +380,7 @@ void Manifold::Impl::CompactProps() {
  */
 void Manifold::Impl::GetFaceBoxMorton(Vec<Box>& faceBox,
                                       Vec<uint32_t>& faceMorton) const {
+  ZoneScoped;
   faceBox.resize(NumTri());
   faceMorton.resize(NumTri());
   for_each_n(autoPolicy(NumTri()),
@@ -389,6 +393,7 @@ void Manifold::Impl::GetFaceBoxMorton(Vec<Box>& faceBox,
  * bounding box and Morton code arrays are also sorted accordingly.
  */
 void Manifold::Impl::SortFaces(Vec<Box>& faceBox, Vec<uint32_t>& faceMorton) {
+  ZoneScoped;
   Vec<int> faceNew2Old(NumTri());
   auto policy = autoPolicy(faceNew2Old.size());
   sequence(policy, faceNew2Old.begin(), faceNew2Old.end());
@@ -419,6 +424,7 @@ void Manifold::Impl::SortFaces(Vec<Box>& faceBox, Vec<uint32_t>& faceMorton) {
  * faces to gather into this.
  */
 void Manifold::Impl::GatherFaces(const Vec<int>& faceNew2Old) {
+  ZoneScoped;
   const int numTri = faceNew2Old.size();
   if (meshRelation_.triRef.size() == NumTri())
     Permute(meshRelation_.triRef, faceNew2Old);
@@ -441,6 +447,7 @@ void Manifold::Impl::GatherFaces(const Vec<int>& faceNew2Old) {
 }
 
 void Manifold::Impl::GatherFaces(const Impl& old, const Vec<int>& faceNew2Old) {
+  ZoneScoped;
   const int numTri = faceNew2Old.size();
   auto policy = autoPolicy(numTri);
 
@@ -510,6 +517,7 @@ MeshGL::MeshGL(const Mesh& mesh) {
  * will report a Status if it is not manifold.
  */
 bool MeshGL::Merge() {
+  ZoneScoped;
   std::multiset<std::pair<int, int>> openEdges;
 
   std::vector<int> merge(NumVert());

src/polygon/src/polygon.cpp

@@ -13,26 +13,17 @@
 // limitations under the License.
 
 #include "polygon.h"
-#if MANIFOLD_PAR == 'T'
-#include "tbb/tbb.h"
-#endif
 
 #include <algorithm>
-#include <numeric>
-#if MANIFOLD_PAR == 'T' && TBB_INTERFACE_VERSION >= 10000 && \
-    __has_include(<pstl/glue_execution_defs.h>)
-#include <execution>
-#endif
 #include <list>
 #include <map>
-#if __has_include(<memory_resource>)
-#include <memory_resource>
-#endif
+#include <numeric>
 #include <queue>
 #include <set>
 #include <stack>
 
 #include "optional_assert.h"
+#include "utils.h"
 
 namespace {
 using namespace manifold;
@@ -181,6 +172,8 @@ void PrintFailure(const std::exception &e, const PolygonsIdx &polys,
 class EarClip {
  public:
   EarClip(const PolygonsIdx &polys, float precision) : precision_(precision) {
+    ZoneScoped;
+
     int numVert = 0;
     for (const SimplePolygonIdx &poly : polys) {
       numVert += poly.size();
@@ -199,6 +192,8 @@ class EarClip {
   }
 
   std::vector<glm::ivec3> Triangulate() {
+    ZoneScoped;
+
     for (const VertItr start : holes_) {
       CutKeyhole(start);
     }
@@ -286,7 +281,8 @@ class EarClip {
       VertItr center = tail;
       VertItr last = center;
 
-      while (nextL != nextR && tail != nextR) {
+      while (nextL != nextR && tail != nextR &&
+             nextL != (toLeft ? right : left)) {
         const glm::vec2 edgeL = nextL->pos - center->pos;
         const float l2 = glm::dot(edgeL, edgeL);
         if (l2 <= p2) {
@@ -550,6 +546,9 @@ class EarClip {
       auto vert = poly.begin();
       polygon_.push_back({vert->idx, earsQueue_.end(), vert->pos});
       const VertItr first = std::prev(polygon_.end());
+
+      bound = glm::max(
+          bound, glm::max(glm::abs(first->pos.x), glm::abs(first->pos.y)));
       VertItr last = first;
       // This is not the real rightmost start, but just an arbitrary vert for
       // now to identify each polygon.
@@ -578,6 +577,8 @@ class EarClip {
   // Find the actual rightmost starts after degenerate removal. Also calculate
   // the polygon bounding boxes.
   void FindStart(VertItr first) {
+    const glm::vec2 origin = first->pos;
+
     VertItr start = first;
     float maxX = -std::numeric_limits<float>::infinity();
     Rect bBox;
@@ -587,7 +588,8 @@ class EarClip {
 
     auto AddPoint = [&](VertItr v) {
       bBox.Union(v->pos);
-      const double area1 = glm::determinant(glm::dmat2(v->pos, v->right->pos));
+      const double area1 =
+          glm::determinant(glm::dmat2(v->pos - origin, v->right->pos - origin));
       const double t1 = area + area1;
       areaCompensation += (area - t1) + area1;
       area = t1;
@@ -605,9 +607,9 @@ class EarClip {
 
     area += areaCompensation;
     const glm::vec2 size = bBox.Size();
-    const double minArea = precision_ * glm::max(size.x, size.y);
+    const float minArea = precision_ * glm::max(size.x, size.y);
 
-    if (area < -minArea) {
+    if (glm::isfinite(maxX) && area < -minArea) {
       holes_.insert(start);
       hole2BBox_.insert({start, bBox});
     } else {
@@ -623,7 +625,6 @@ class EarClip {
   // incorrect due to precision, we check for polygon edges both ahead and
   // behind to ensure all valid options are found.
   void CutKeyhole(const VertItr start) {
-    const float startX = start->pos.x;
     const Rect bBox = hole2BBox_[start];
     const int onTop = start->pos.y >= bBox.max.y - precision_   ? 1
                       : start->pos.y <= bBox.min.y + precision_ ? -1
@@ -668,7 +669,6 @@ class EarClip {
   // and returns it. It does so by finding any reflex verts inside the triangle
   // containing the best connection and the initial horizontal line.
   VertItr FindCloserBridge(VertItr start, VertItr edge, int onTop) {
-    const float p2 = precision_ * precision_;
     VertItr best = edge->pos.x > edge->right->pos.x ? edge : edge->right;
     const float maxX = best->pos.x;
     const float above = best->pos.y > start->pos.y ? 1 : -1;
@@ -743,6 +743,8 @@ class EarClip {
   // The main ear-clipping loop. This is called once for each simple polygon -
   // all holes have already been key-holed and joined to an outer polygon.
   void TriangulatePoly(VertItr start) {
+    ZoneScoped;
+
     // A simple polygon always creates two fewer triangles than it has verts.
     int numTri = -2;
     earsQueue_.clear();

src/sdf/src/sdf.cpp

@@ -151,6 +151,7 @@ struct ComputeVerts {
   }
 
   inline void operator()(Uint64 mortonCode) {
+    ZoneScoped;
     if (gridVerts.Full()) return;
 
     const glm::ivec4 gridIndex = DecodeMorton(mortonCode);
@@ -205,6 +206,7 @@ struct BuildTris {
   }
 
   void operator()(int idx) {
+    ZoneScoped;
     Uint64 basekey = gridVerts.KeyAt(idx);
     if (basekey == kOpen) return;
 
@@ -305,7 +307,6 @@ Mesh LevelSet(std::function<float(glm::vec3)> sdf, Box bounds, float edgeLength,
   Mesh out;
 
   const glm::vec3 dim = bounds.Size();
-  const float maxDim = std::max(dim[0], std::max(dim[1], dim[2]));
   const glm::ivec3 gridSize(dim / edgeLength);
   const glm::vec3 spacing = dim / (glm::vec3(gridSize));
 

src/utilities/CMakeLists.txt

@@ -32,21 +32,22 @@ FetchContent_Declare(Thrust
     GIT_TAG 2.1.0
     GIT_SHALLOW TRUE
     GIT_PROGRESS TRUE
-    FIND_PACKAGE_ARGS NAMES Thrust
 )
-FetchContent_MakeAvailable(Thrust)
+find_package(Thrust QUIET)
+if(NOT Thrust_FOUND AND NOT DEFINED thrust_SOURCE_DIR)
+    FetchContent_Populate(Thrust)
+endif()
 
 if (TRACY_ENABLE)
     include(FetchContent)
     FetchContent_Declare(tracy
         GIT_REPOSITORY https://github.com/wolfpld/tracy.git
-        GIT_TAG v0.9.1
+        GIT_TAG v0.10
         GIT_SHALLOW TRUE
         GIT_PROGRESS TRUE
     )
     FetchContent_MakeAvailable(tracy)
     target_link_libraries(${PROJECT_NAME} INTERFACE TracyClient)
-    target_compile_options(${PROJECT_NAME} INTERFACE -DTRACY_NO_CODE_TRANSFER=1)
 endif()
 
 if(MANIFOLD_PAR STREQUAL "TBB")
@@ -61,12 +62,12 @@ if(MANIFOLD_PAR STREQUAL "TBB")
         set(TBB_STRICT OFF CACHE INTERNAL "" FORCE)
         FetchContent_Declare(TBB
             GIT_REPOSITORY https://github.com/oneapi-src/oneTBB.git
-            GIT_TAG        v2021.10.0
+            GIT_TAG        v2021.11.0
             GIT_SHALLOW    TRUE
             GIT_PROGRESS   TRUE
         )
         FetchContent_MakeAvailable(TBB)
-        set_property(DIRECTORY ${TBB_SOURCE_DIR} PROPERTY EXCLUDE_FROM_ALL YES)
+        set_property(DIRECTORY ${tbb_SOURCE_DIR} PROPERTY EXCLUDE_FROM_ALL YES)
     endif()
     target_compile_options(${PROJECT_NAME} INTERFACE -DMANIFOLD_PAR='T')
     if(TARGET TBB::tbb)
@@ -84,15 +85,15 @@ endif()
 
 target_include_directories(${PROJECT_NAME} INTERFACE ${PROJECT_SOURCE_DIR}/include)
 target_link_libraries(${PROJECT_NAME} INTERFACE glm::glm)
-if(DEFINED Thrust_SOURCE_DIR)
+if(DEFINED thrust_SOURCE_DIR)
 target_include_directories(${PROJECT_NAME} INTERFACE
-    ${Thrust_SOURCE_DIR}
-    ${Thrust_SOURCE_DIR}/dependencies/libcudacxx/include
+    ${thrust_SOURCE_DIR}
+    ${thrust_SOURCE_DIR}/dependencies/libcudacxx/include
 )
 else()
 target_include_directories(${PROJECT_NAME} INTERFACE
-    ${THRUST_ROOT}
-    ${THRUST_ROOT}/dependencies/libcudacxx/include
+    ${_THRUST_INCLUDE_DIR}
+    ${_THRUST_INCLUDE_DIR}/dependencies/libcudacxx/include
 )
 endif()
 

src/utilities/include/public.h

@@ -613,6 +613,8 @@ struct ExecutionParams {
   bool suppressErrors = false;
   /// Deterministic outputs. Will disable some parallel optimizations.
   bool deterministic = false;
+  /// Perform optional but recommended triangle cleanups in SimplifyTopology()
+  bool cleanupTriangles = true;
 };
 
 #ifdef MANIFOLD_DEBUG
@@ -633,6 +635,13 @@ inline std::ostream& operator<<(std::ostream& stream, const glm::tvec4<T>& v) {
                 << ", w = " << v.w;
 }
 
+inline std::ostream& operator<<(std::ostream& stream, const glm::mat3& mat) {
+  glm::mat3 tam = glm::transpose(mat);
+  return stream << tam[0] << std::endl
+                << tam[1] << std::endl
+                << tam[2] << std::endl;
+}
+
 inline std::ostream& operator<<(std::ostream& stream, const glm::mat4x3& mat) {
   glm::mat3x4 tam = glm::transpose(mat);
   return stream << tam[0] << std::endl

src/utilities/include/sparse.h

@@ -48,7 +48,11 @@ class SparseIndices {
   }
 
   SparseIndices() = default;
-  SparseIndices(size_t size) : data_(size * sizeof(int64_t)) {}
+  SparseIndices(size_t size) {
+    if (size >= std::numeric_limits<int32_t>::max() / sizeof(int64_t))
+      throw std::out_of_range("SparseIndices too large");
+    data_ = Vec<char>(size * sizeof(int64_t));
+  }
 
   int size() const { return data_.size() / sizeof(int64_t); }
 

src/utilities/include/svd.h

@@ -0,0 +1,310 @@
+// MIT License
+
+// Copyright (c) 2019 wi-re
+// Copyright 2023 The Manifold Authors.
+
+// Permission is hereby granted, free of charge, to any person obtaining a copy
+// of this software and associated documentation files (the "Software"), to deal
+// in the Software without restriction, including without limitation the rights
+// to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+// copies of the Software, and to permit persons to whom the Software is
+// furnished to do so, subject to the following conditions:
+
+// The above copyright notice and this permission notice shall be included in
+// all copies or substantial portions of the Software.
+
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+// SOFTWARE.
+
+// Modified from https://github.com/wi-re/tbtSVD, removing CUDA dependence and
+// approximate inverse square roots.
+
+#include <algorithm>
+#include <cmath>
+#include <glm/glm.hpp>
+#include <tuple>
+#include <utility>
+
+namespace {
+// Constants used for calculation of Givens quaternions
+inline constexpr float _gamma = 5.828427124f;   // sqrt(8)+3;
+inline constexpr float _cStar = 0.923879532f;   // cos(pi/8)
+inline constexpr float _sStar = 0.3826834323f;  // sin(pi/8)
+// Threshold value
+inline constexpr float _SVD_EPSILON = 1e-6f;
+// Iteration counts for Jacobi Eigen Analysis, influences precision
+inline constexpr int JACOBI_STEPS = 12;
+
+// Helper function used to swap X with Y and Y with  X if c == true
+inline void CondSwap(bool c, float& X, float& Y) {
+  float Z = X;
+  X = c ? Y : X;
+  Y = c ? Z : Y;
+}
+// Helper function used to swap X with Y and Y with -X if c == true
+inline void CondNegSwap(bool c, float& X, float& Y) {
+  float Z = -X;
+  X = c ? Y : X;
+  Y = c ? Z : Y;
+}
+// A simple symmetric 3x3 Matrix class (contains no storage for (0, 1) (0, 2)
+// and (1, 2)
+struct Symmetric3x3 {
+  float m_00 = 1.f;
+  float m_10 = 0.f, m_11 = 1.f;
+  float m_20 = 0.f, m_21 = 0.f, m_22 = 1.f;
+
+  Symmetric3x3(float a11 = 1.f, float a21 = 0.f, float a22 = 1.f,
+               float a31 = 0.f, float a32 = 0.f, float a33 = 1.f)
+      : m_00(a11), m_10(a21), m_11(a22), m_20(a31), m_21(a32), m_22(a33) {}
+  Symmetric3x3(glm::mat3 o)
+      : m_00(o[0][0]),
+        m_10(o[0][1]),
+        m_11(o[1][1]),
+        m_20(o[0][2]),
+        m_21(o[1][2]),
+        m_22(o[2][2]) {}
+};
+// Helper struct to store 2 floats to avoid OUT parameters on functions
+struct Givens {
+  float ch = _cStar;
+  float sh = _sStar;
+};
+// Helper struct to store 2 Matrices to avoid OUT parameters on functions
+struct QR {
+  glm::mat3 Q, R;
+};
+// Calculates the squared norm of the vector.
+inline float Dist2(glm::vec3 v) { return glm::dot(v, v); }
+// For an explanation of the math see
+// http://pages.cs.wisc.edu/~sifakis/papers/SVD_TR1690.pdf Computing the
+// Singular Value Decomposition of 3 x 3 matrices with minimal branching and
+// elementary floating point operations See Algorithm 2 in reference. Given a
+// matrix A this function returns the Givens quaternion (x and w component, y
+// and z are 0)
+inline Givens ApproximateGivensQuaternion(Symmetric3x3& A) {
+  Givens g{2.f * (A.m_00 - A.m_11), A.m_10};
+  bool b = _gamma * g.sh * g.sh < g.ch * g.ch;
+  float w = 1.f / sqrt(fmaf(g.ch, g.ch, g.sh * g.sh));
+  if (w != w) b = 0;
+  return Givens{b ? w * g.ch : (float)_cStar, b ? w * g.sh : (float)_sStar};
+}
+// Function used to apply a Givens rotation S. Calculates the weights and
+// updates the quaternion to contain the cumulative rotation
+inline void JacobiConjugation(const int32_t x, const int32_t y, const int32_t z,
+                              Symmetric3x3& S, glm::vec4& q) {
+  auto g = ApproximateGivensQuaternion(S);
+  float scale = 1.f / fmaf(g.ch, g.ch, g.sh * g.sh);
+  float a = fmaf(g.ch, g.ch, -g.sh * g.sh) * scale;
+  float b = 2.f * g.sh * g.ch * scale;
+  Symmetric3x3 _S = S;
+  // perform conjugation S = Q'*S*Q
+  S.m_00 = fmaf(a, fmaf(a, _S.m_00, b * _S.m_10),
+                b * (fmaf(a, _S.m_10, b * _S.m_11)));
+  S.m_10 = fmaf(a, fmaf(-b, _S.m_00, a * _S.m_10),
+                b * (fmaf(-b, _S.m_10, a * _S.m_11)));
+  S.m_11 = fmaf(-b, fmaf(-b, _S.m_00, a * _S.m_10),
+                a * (fmaf(-b, _S.m_10, a * _S.m_11)));
+  S.m_20 = fmaf(a, _S.m_20, b * _S.m_21);
+  S.m_21 = fmaf(-b, _S.m_20, a * _S.m_21);
+  S.m_22 = _S.m_22;
+  // update cumulative rotation qV
+  glm::vec3 tmp = g.sh * glm::vec3(q);
+  g.sh *= q[3];
+  // (x,y,z) corresponds to ((0,1,2),(1,2,0),(2,0,1)) for (p,q) =
+  // ((0,1),(1,2),(0,2))
+  q[z] = fmaf(q[z], g.ch, g.sh);
+  q[3] = fmaf(q[3], g.ch, -tmp[z]);  // w
+  q[x] = fmaf(q[x], g.ch, tmp[y]);
+  q[y] = fmaf(q[y], g.ch, -tmp[x]);
+  // re-arrange matrix for next iteration
+  _S.m_00 = S.m_11;
+  _S.m_10 = S.m_21;
+  _S.m_11 = S.m_22;
+  _S.m_20 = S.m_10;
+  _S.m_21 = S.m_20;
+  _S.m_22 = S.m_00;
+  S.m_00 = _S.m_00;
+  S.m_10 = _S.m_10;
+  S.m_11 = _S.m_11;
+  S.m_20 = _S.m_20;
+  S.m_21 = _S.m_21;
+  S.m_22 = _S.m_22;
+}
+// Function used to contain the Givens permutations and the loop of the jacobi
+// steps controlled by JACOBI_STEPS Returns the quaternion q containing the
+// cumulative result used to reconstruct S
+inline glm::mat3 JacobiEigenAnalysis(Symmetric3x3 S) {
+  glm::vec4 q(0, 0, 0, 1);
+  for (int32_t i = 0; i < JACOBI_STEPS; i++) {
+    JacobiConjugation(0, 1, 2, S, q);
+    JacobiConjugation(1, 2, 0, S, q);
+    JacobiConjugation(2, 0, 1, S, q);
+  }
+  return glm::mat3(1.f - 2.f * (fmaf(q.y, q.y, q.z * q.z)),  //
+                   2.f * fmaf(q.x, q.y, +q.w * q.z),         //
+                   2.f * fmaf(q.x, q.z, -q.w * q.y),         //
+                   2 * fmaf(q.x, q.y, -q.w * q.z),           //
+                   1 - 2 * fmaf(q.x, q.x, q.z * q.z),        //
+                   2 * fmaf(q.y, q.z, q.w * q.x),            //
+                   2 * fmaf(q.x, q.z, q.w * q.y),            //
+                   2 * fmaf(q.y, q.z, -q.w * q.x),           //
+                   1 - 2 * fmaf(q.x, q.x, q.y * q.y));
+}
+// Implementation of Algorithm 3
+inline void SortSingularValues(glm::mat3& B, glm::mat3& V) {
+  float rho1 = Dist2(B[0]);
+  float rho2 = Dist2(B[1]);
+  float rho3 = Dist2(B[2]);
+  bool c;
+  c = rho1 < rho2;
+  CondNegSwap(c, B[0][0], B[1][0]);
+  CondNegSwap(c, V[0][0], V[1][0]);
+  CondNegSwap(c, B[0][1], B[1][1]);
+  CondNegSwap(c, V[0][1], V[1][1]);
+  CondNegSwap(c, B[0][2], B[1][2]);
+  CondNegSwap(c, V[0][2], V[1][2]);
+  CondSwap(c, rho1, rho2);
+  c = rho1 < rho3;
+  CondNegSwap(c, B[0][0], B[2][0]);
+  CondNegSwap(c, V[0][0], V[2][0]);
+  CondNegSwap(c, B[0][1], B[2][1]);
+  CondNegSwap(c, V[0][1], V[2][1]);
+  CondNegSwap(c, B[0][2], B[2][2]);
+  CondNegSwap(c, V[0][2], V[2][2]);
+  CondSwap(c, rho1, rho3);
+  c = rho2 < rho3;
+  CondNegSwap(c, B[1][0], B[2][0]);
+  CondNegSwap(c, V[1][0], V[2][0]);
+  CondNegSwap(c, B[1][1], B[2][1]);
+  CondNegSwap(c, V[1][1], V[2][1]);
+  CondNegSwap(c, B[1][2], B[2][2]);
+  CondNegSwap(c, V[1][2], V[2][2]);
+}
+// Implementation of Algorithm 4
+inline Givens QRGivensQuaternion(float a1, float a2) {
+  // a1 = pivot point on diagonal
+  // a2 = lower triangular entry we want to annihilate
+  float epsilon = (float)_SVD_EPSILON;
+  float rho = sqrt(fmaf(a1, a1, +a2 * a2));
+  Givens g{fabsf(a1) + fmaxf(rho, epsilon), rho > epsilon ? a2 : 0};
+  bool b = a1 < 0.f;
+  CondSwap(b, g.sh, g.ch);
+  float w = 1.f / sqrt(fmaf(g.ch, g.ch, g.sh * g.sh));
+  g.ch *= w;
+  g.sh *= w;
+  return g;
+}
+// Implements a QR decomposition of a Matrix, see Sec 4.2
+inline QR QRDecomposition(glm::mat3& B) {
+  glm::mat3 Q, R;
+  // first Givens rotation (ch,0,0,sh)
+  auto g1 = QRGivensQuaternion(B[0][0], B[0][1]);
+  auto a = fmaf(-2.f, g1.sh * g1.sh, 1.f);
+  auto b = 2.f * g1.ch * g1.sh;
+  // apply B = Q' * B
+  R[0][0] = fmaf(a, B[0][0], b * B[0][1]);
+  R[1][0] = fmaf(a, B[1][0], b * B[1][1]);
+  R[2][0] = fmaf(a, B[2][0], b * B[2][1]);
+  R[0][1] = fmaf(-b, B[0][0], a * B[0][1]);
+  R[1][1] = fmaf(-b, B[1][0], a * B[1][1]);
+  R[2][1] = fmaf(-b, B[2][0], a * B[2][1]);
+  R[0][2] = B[0][2];
+  R[1][2] = B[1][2];
+  R[2][2] = B[2][2];
+  // second Givens rotation (ch,0,-sh,0)
+  auto g2 = QRGivensQuaternion(R[0][0], R[0][2]);
+  a = fmaf(-2.f, g2.sh * g2.sh, 1.f);
+  b = 2.f * g2.ch * g2.sh;
+  // apply B = Q' * B;
+  B[0][0] = fmaf(a, R[0][0], b * R[0][2]);
+  B[1][0] = fmaf(a, R[1][0], b * R[1][2]);
+  B[2][0] = fmaf(a, R[2][0], b * R[2][2]);
+  B[0][1] = R[0][1];
+  B[1][1] = R[1][1];
+  B[2][1] = R[2][1];
+  B[0][2] = fmaf(-b, R[0][0], a * R[0][2]);
+  B[1][2] = fmaf(-b, R[1][0], a * R[1][2]);
+  B[2][2] = fmaf(-b, R[2][0], a * R[2][2]);
+  // third Givens rotation (ch,sh,0,0)
+  auto g3 = QRGivensQuaternion(B[1][1], B[1][2]);
+  a = fmaf(-2.f, g3.sh * g3.sh, 1.f);
+  b = 2.f * g3.ch * g3.sh;
+  // R is now set to desired value
+  R[0][0] = B[0][0];
+  R[1][0] = B[1][0];
+  R[2][0] = B[2][0];
+  R[0][1] = fmaf(a, B[0][1], b * B[0][2]);
+  R[1][1] = fmaf(a, B[1][1], b * B[1][2]);
+  R[2][1] = fmaf(a, B[2][1], b * B[2][2]);
+  R[0][2] = fmaf(-b, B[0][1], a * B[0][2]);
+  R[1][2] = fmaf(-b, B[1][1], a * B[1][2]);
+  R[2][2] = fmaf(-b, B[2][1], a * B[2][2]);
+  // construct the cumulative rotation Q=Q1 * Q2 * Q3
+  // the number of floating point operations for three quaternion
+  // multiplications is more or less comparable to the explicit form of the
+  // joined matrix. certainly more memory-efficient!
+  auto sh12 = 2.f * fmaf(g1.sh, g1.sh, -0.5f);
+  auto sh22 = 2.f * fmaf(g2.sh, g2.sh, -0.5f);
+  auto sh32 = 2.f * fmaf(g3.sh, g3.sh, -0.5f);
+  Q[0][0] = sh12 * sh22;
+  Q[1][0] = fmaf(4.f * g2.ch * g3.ch, sh12 * g2.sh * g3.sh,
+                 2.f * g1.ch * g1.sh * sh32);
+  Q[2][0] = fmaf(4.f * g1.ch * g3.ch, g1.sh * g3.sh,
+                 -2.f * g2.ch * sh12 * g2.sh * sh32);
+
+  Q[0][1] = -2.f * g1.ch * g1.sh * sh22;
+  Q[1][1] =
+      fmaf(-8.f * g1.ch * g2.ch * g3.ch, g1.sh * g2.sh * g3.sh, sh12 * sh32);
+  Q[2][1] = fmaf(
+      -2.f * g3.ch, g3.sh,
+      4.f * g1.sh * fmaf(g3.ch * g1.sh, g3.sh, g1.ch * g2.ch * g2.sh * sh32));
+
+  Q[0][2] = 2.f * g2.ch * g2.sh;
+  Q[1][2] = -2.f * g3.ch * sh22 * g3.sh;
+  Q[2][2] = sh22 * sh32;
+  return QR{Q, R};
+}
+}  // namespace
+
+namespace manifold {
+/** @addtogroup Connections
+ *  @{
+ */
+
+/**
+ * The three matrices of a Singular Value Decomposition.
+ */
+struct SVDSet {
+  glm::mat3 U, S, V;
+};
+
+/**
+ * Returns the Singular Value Decomposition of A: A = U * S * glm::transpose(V).
+ *
+ * @param A The matrix to decompose.
+ */
+inline SVDSet SVD(glm::mat3 A) {
+  glm::mat3 V = JacobiEigenAnalysis(glm::transpose(A) * A);
+  auto B = A * V;
+  SortSingularValues(B, V);
+  QR qr = QRDecomposition(B);
+  return SVDSet{qr.Q, qr.R, V};
+}
+
+/**
+ * Returns the largest singular value of A.
+ *
+ * @param A The matrix to measure.
+ */
+inline float SpectralNorm(glm::mat3 A) {
+  SVDSet usv = SVD(A);
+  return usv.S[0][0];
+}
+/** @} */
+}  // namespace manifold

src/utilities/include/utils.h

@@ -26,6 +26,18 @@
 #include "par.h"
 #include "vec.h"
 
+#if __has_include(<tracy/Tracy.hpp>)
+#include <tracy/Tracy.hpp>
+#else
+#define FrameMarkStart(x)
+#define FrameMarkEnd(x)
+// putting ZoneScoped in a function will instrument the function execution when
+// TRACY_ENABLE is set, which allows the profiler to record more accurate
+// timing.
+#define ZoneScoped
+#define ZoneScopedN(name)
+#endif
+
 namespace manifold {
 
 /** @defgroup Private

src/utilities/include/vec.h

@@ -24,6 +24,7 @@
 // #include "optional_assert.h"
 #include "par.h"
 #include "public.h"
+#include "vec_view.h"
 
 namespace manifold {
 
@@ -33,100 +34,6 @@ namespace manifold {
 template <typename T>
 class Vec;
 
-/**
- * View for Vec, can perform offset operation.
- * This will be invalidated when the original vector is dropped or changes
- * length.
- */
-template <typename T>
-class VecView {
- public:
-  using Iter = T *;
-  using IterC = const T *;
-
-  VecView(T *ptr_, int size_) : ptr_(ptr_), size_(size_) {}
-
-  VecView(const VecView &other) {
-    ptr_ = other.ptr_;
-    size_ = other.size_;
-  }
-
-  VecView &operator=(const VecView &other) {
-    ptr_ = other.ptr_;
-    size_ = other.size_;
-    return *this;
-  }
-
-  // allows conversion to a const VecView
-  operator VecView<const T>() const { return {ptr_, size_}; }
-
-  inline const T &operator[](int i) const {
-    if (i < 0 || i >= size_) throw std::out_of_range("Vec out of range");
-    return ptr_[i];
-  }
-
-  inline T &operator[](int i) {
-    if (i < 0 || i >= size_) throw std::out_of_range("Vec out of range");
-    return ptr_[i];
-  }
-
-  IterC cbegin() const { return ptr_; }
-  IterC cend() const { return ptr_ + size_; }
-
-  IterC begin() const { return cbegin(); }
-  IterC end() const { return cend(); }
-
-  Iter begin() { return ptr_; }
-  Iter end() { return ptr_ + size_; }
-
-  const T &front() const {
-    if (size_ == 0)
-      throw std::out_of_range("attempt to take the front of an empty vector");
-    return ptr_[0];
-  }
-
-  const T &back() const {
-    if (size_ == 0)
-      throw std::out_of_range("attempt to take the back of an empty vector");
-    return ptr_[size_ - 1];
-  }
-
-  T &front() {
-    if (size_ == 0)
-      throw std::out_of_range("attempt to take the front of an empty vector");
-    return ptr_[0];
-  }
-
-  T &back() {
-    if (size_ == 0)
-      throw std::out_of_range("attempt to take the back of an empty vector");
-    return ptr_[size_ - 1];
-  }
-
-  int size() const { return size_; }
-
-  bool empty() const { return size_ == 0; }
-
-#ifdef MANIFOLD_DEBUG
-  void Dump() {
-    std::cout << "Vec = " << std::endl;
-    for (int i = 0; i < size(); ++i) {
-      std::cout << i << ", " << ptr_[i] << ", " << std::endl;
-    }
-    std::cout << std::endl;
-  }
-#endif
-
- protected:
-  T *ptr_ = nullptr;
-  int size_ = 0;
-
-  VecView() = default;
-  friend class Vec<T>;
-  friend class Vec<typename std::remove_const<T>::type>;
-  friend class VecView<typename std::remove_const<T>::type>;
-};
-
 /*
  * Specialized vector implementation with multithreaded fill and uninitialized
  * memory optimizations.
@@ -157,7 +64,7 @@ class Vec : public VecView<T> {
     if (this->size_ != 0) {
       this->ptr_ = reinterpret_cast<T *>(malloc(this->size_ * sizeof(T)));
       if (this->ptr_ == nullptr) throw std::bad_alloc();
-      TracyAllocS(ptr_, size_ * sizeof(T), 3);
+      TracyAllocS(this->ptr_, this->size_ * sizeof(T), 3);
       uninitialized_copy(policy, vec.begin(), vec.end(), this->ptr_);
     }
   }
@@ -169,7 +76,7 @@ class Vec : public VecView<T> {
     if (this->size_ != 0) {
       this->ptr_ = reinterpret_cast<T *>(malloc(this->size_ * sizeof(T)));
       if (this->ptr_ == nullptr) throw std::bad_alloc();
-      TracyAllocS(ptr_, size_ * sizeof(T), 3);
+      TracyAllocS(this->ptr_, this->size_ * sizeof(T), 3);
       uninitialized_copy(policy, vec.begin(), vec.end(), this->ptr_);
     }
   }
@@ -188,7 +95,7 @@ class Vec : public VecView<T> {
 
   ~Vec() {
     if (this->ptr_ != nullptr) {
-      TracyFreeS(ptr_, 3);
+      TracyFreeS(this->ptr_, 3);
       free(this->ptr_);
     }
     this->ptr_ = nullptr;
@@ -208,7 +115,7 @@ class Vec : public VecView<T> {
     if (this->size_ != 0) {
       this->ptr_ = reinterpret_cast<T *>(malloc(this->size_ * sizeof(T)));
       if (this->ptr_ == nullptr) throw std::bad_alloc();
-      TracyAllocS(ptr_, size_ * sizeof(T), 3);
+      TracyAllocS(this->ptr_, this->size_ * sizeof(T), 3);
       uninitialized_copy(policy, other.begin(), other.end(), this->ptr_);
     }
     return *this;
@@ -217,7 +124,7 @@ class Vec : public VecView<T> {
   Vec<T> &operator=(Vec<T> &&other) {
     if (&other == this) return *this;
     if (this->ptr_ != nullptr) {
-      TracyFreeS(ptr_, 3);
+      TracyFreeS(this->ptr_, 3);
       free(this->ptr_);
     }
     this->size_ = other.size_;
@@ -257,7 +164,7 @@ class Vec : public VecView<T> {
         uninitialized_copy(autoPolicy(this->size_), this->ptr_,
                            this->ptr_ + this->size_, newBuffer);
       if (this->ptr_ != nullptr) {
-        TracyFreeS(ptr_, 3);
+        TracyFreeS(this->ptr_, 3);
         free(this->ptr_);
       }
       this->ptr_ = newBuffer;
@@ -281,12 +188,12 @@ class Vec : public VecView<T> {
     if (this->size_ > 0) {
       newBuffer = reinterpret_cast<T *>(malloc(this->size_ * sizeof(T)));
       if (newBuffer == nullptr) throw std::bad_alloc();
-      TracyAllocS(newBuffer, size_ * sizeof(T), 3);
+      TracyAllocS(newBuffer, this->size_ * sizeof(T), 3);
       uninitialized_copy(autoPolicy(this->size_), this->ptr_,
                          this->ptr_ + this->size_, newBuffer);
     }
     if (this->ptr_ != nullptr) {
-      TracyFreeS(ptr_, 3);
+      TracyFreeS(this->ptr_, 3);
       free(this->ptr_);
     }
     this->ptr_ = newBuffer;

src/utilities/include/vec_view.h

@@ -0,0 +1,112 @@
+// Copyright 2023 The Manifold Authors.
+//
+// Licensed under the Apache License, Version 2.0 (the "License");
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+//
+//      http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+
+#pragma once
+
+#include <stdexcept>
+
+namespace manifold {
+
+/**
+ * View for Vec, can perform offset operation.
+ * This will be invalidated when the original vector is dropped or changes
+ * length. Roughly equivalent to std::span<T> from c++20
+ */
+template <typename T>
+class VecView {
+ public:
+  using Iter = T *;
+  using IterC = const T *;
+
+  VecView(T *ptr_, int size_) : ptr_(ptr_), size_(size_) {}
+
+  VecView(const VecView &other) {
+    ptr_ = other.ptr_;
+    size_ = other.size_;
+  }
+
+  VecView &operator=(const VecView &other) {
+    ptr_ = other.ptr_;
+    size_ = other.size_;
+    return *this;
+  }
+
+  // allows conversion to a const VecView
+  operator VecView<const T>() const { return {ptr_, size_}; }
+
+  inline const T &operator[](int i) const {
+    if (i < 0 || i >= size_) throw std::out_of_range("Vec out of range");
+    return ptr_[i];
+  }
+
+  inline T &operator[](int i) {
+    if (i < 0 || i >= size_) throw std::out_of_range("Vec out of range");
+    return ptr_[i];
+  }
+
+  IterC cbegin() const { return ptr_; }
+  IterC cend() const { return ptr_ + size_; }
+
+  IterC begin() const { return cbegin(); }
+  IterC end() const { return cend(); }
+
+  Iter begin() { return ptr_; }
+  Iter end() { return ptr_ + size_; }
+
+  const T &front() const {
+    if (size_ == 0)
+      throw std::out_of_range("attempt to take the front of an empty vector");
+    return ptr_[0];
+  }
+
+  const T &back() const {
+    if (size_ == 0)
+      throw std::out_of_range("attempt to take the back of an empty vector");
+    return ptr_[size_ - 1];
+  }
+
+  T &front() {
+    if (size_ == 0)
+      throw std::out_of_range("attempt to take the front of an empty vector");
+    return ptr_[0];
+  }
+
+  T &back() {
+    if (size_ == 0)
+      throw std::out_of_range("attempt to take the back of an empty vector");
+    return ptr_[size_ - 1];
+  }
+
+  int size() const { return size_; }
+
+  bool empty() const { return size_ == 0; }
+
+#ifdef MANIFOLD_DEBUG
+  void Dump() {
+    std::cout << "Vec = " << std::endl;
+    for (int i = 0; i < size(); ++i) {
+      std::cout << i << ", " << ptr_[i] << ", " << std::endl;
+    }
+    std::cout << std::endl;
+  }
+#endif
+
+ protected:
+  T *ptr_ = nullptr;
+  int size_ = 0;
+
+  VecView() = default;
+};
+
+}  // namespace manifold

test/CMakeLists.txt

@@ -14,6 +14,7 @@
 
 project(manifold_test)
 
+set(CMAKE_CXX_STANDARD 17)
 enable_testing()
 
 set(SOURCE_FILES polygon_test.cpp cross_section_test.cpp manifold_test.cpp boolean_test.cpp sdf_test.cpp samples_test.cpp test_main.cpp)
@@ -44,6 +45,21 @@ endif()
 add_executable(${PROJECT_NAME} ${SOURCE_FILES})
 target_link_libraries(${PROJECT_NAME} polygon GTest::gtest_main manifold sdf samples cross_section)
 
+if(MANIFOLD_FUZZ)
+  FetchContent_Declare(fuzztest
+    GIT_REPOSITORY https://github.com/google/fuzztest.git
+    # note that if commit hash is used, it cannot be a shallow clone
+    GIT_TAG        2606e04a43e5a7730e437a849604a61f1cb0ff28
+    GIT_PROGRESS   TRUE
+  )
+  FetchContent_MakeAvailable(fuzztest)
+  fuzztest_setup_fuzzing_flags()
+  add_executable(polygon_fuzz polygon_fuzz.cpp)
+  target_link_libraries(polygon_fuzz PUBLIC polygon)
+  link_fuzztest(polygon_fuzz)
+  gtest_discover_tests(polygon_fuzz)
+endif()
+
 if(MANIFOLD_CBIND)
   target_link_libraries(${PROJECT_NAME} manifoldc)
 endif()

test/boolean_test.cpp

@@ -163,7 +163,7 @@ TEST(Boolean, Coplanar) {
   MeshGL cylinderGL = WithPositionColors(cylinder);
   cylinder = Manifold(cylinderGL);
 
-  Manifold cylinder2 = cylinder.Scale({0.8f, 0.8f, 1.0f}).Rotate(0, 0, 5);
+  Manifold cylinder2 = cylinder.Scale({0.8f, 0.8f, 1.0f}).Rotate(0, 0, 185);
   Manifold out = cylinder - cylinder2;
   ExpectMeshes(out, {{32, 64, 3, 48}});
   EXPECT_EQ(out.NumDegenerateTris(), 0);
@@ -187,7 +187,7 @@ TEST(Boolean, CoplanarProp) {
   MeshGL cylinderGL = WithIndexColors(cylinder.GetMesh());
   cylinder = Manifold(cylinderGL);
 
-  Manifold cylinder2 = cylinder.Scale({0.8f, 0.8f, 1.0f}).Rotate(0, 0, 5);
+  Manifold cylinder2 = cylinder.Scale({0.8f, 0.8f, 1.0f}).Rotate(0, 0, 185);
   Manifold out = cylinder - cylinder2;
   ExpectMeshes(out, {{52, 104, 3, 88}});
   EXPECT_EQ(out.NumDegenerateTris(), 0);
@@ -659,6 +659,23 @@ TEST(Boolean, TreeTransforms) {
   EXPECT_FLOAT_EQ((a + b).GetProperties().volume, 2);
 }
 
+TEST(Boolean, Spiral) {
+  ManifoldParams().deterministic = true;
+  const int d = 2;
+  std::function<Manifold(const int, const float, const float)> spiral =
+      [&](const int rec, const float r, const float add) {
+        const float rot = 360.0f / (glm::pi<float>() * r * 2) * d;
+        const float rNext = r + add / 360 * rot;
+        const Manifold cube =
+            Manifold::Cube(glm::vec3(1), true).Translate({0, r, 0});
+        if (rec > 0)
+          return spiral(rec - 1, rNext, add).Rotate(0, 0, rot) + cube;
+        return cube;
+      };
+  const Manifold result = spiral(120, 25, 2);
+  EXPECT_EQ(result.Genus(), -120);
+}
+
 TEST(Boolean, Sweep) {
   PolygonParams().processOverlaps = true;
 
@@ -902,3 +919,452 @@ TEST(Boolean, Sweep) {
 
   PolygonParams().processOverlaps = false;
 }
+
+TEST(Boolean, InterpolatedNormals) {
+  MeshGL a;
+  a.numProp = 8;
+  a.vertProperties = {
+      // 0
+      -409.0570983886719, -300, -198.83624267578125, 0, -1, 0,
+      590.9429321289062, 301.1637268066406,
+      // 1
+      -1000, -300, 500, 0, -1, 0, 0, 1000,
+      // 2
+      -1000, -300, -500, 0, -1, 0, 0, 0,
+      // 3
+      -1000, -300, -500, -1, 0, 0, 600, 0,
+      // 4
+      -1000, -300, 500, -1, 0, 0, 600, 1000,
+      // 5
+      -1000, 300, -500, -1, 0, 0, 0, 0,
+      // 6
+      7.179656982421875, -300, -330.03717041015625, 0, -1, 0,
+      1007.1796264648438, 169.9628448486328,
+      // 7
+      1000, 300, 500, 0, 0, 1, 2000, 600,
+      // 8
+      403.5837097167969, 300, 500, 0, 0, 1, 1403.583740234375, 600,
+      // 9
+      564.2904052734375, 21.64801025390625, 500, 0, 0, 1, 1564.29052734375,
+      321.64801025390625,
+      // 10
+      1000, -300, -500, 0, 0, -1, 2000, 600,
+      // 11
+      -1000, -300, -500, 0, 0, -1, 0, 600,
+      // 12
+      -1000, 300, -500, 0, 0, -1, 0, 0,
+      // 13
+      1000, 300, 500, 0, 1, 0, 0, 1000,
+      // 14
+      1000, 300, -500, 0, 1, 0, 0, 0,
+      // 15
+      724.5271606445312, 300, 398.83624267578125, 0, 1, 0, 275.47283935546875,
+      898.8362426757812,
+      // 16
+      -115.35255432128906, -300, 500, 0, -1, 0, 884.6475219726562,
+      1000.0001220703125,
+      // 17
+      -384.7195129394531, 166.55722045898438, 500, 0, 0, 1, 615.280517578125,
+      466.5572509765625,
+      // 18
+      -1000, -300, 500, 0, 0, 1, 0, 0,
+      // 19
+      -161.6136932373047, -219.87335205078125, 500, 0, 0, 1, 838.3862915039062,
+      80.12664794921875,
+      // 20
+      1000, -300, 500, 0, 0, 1, 2000, 0,
+      // 21
+      -115.35255432128906, -300, 500, 0, 0, 1, 884.6475219726562, 0,
+      // 22
+      1000, 300, 500, 1, 0, 0, 600, 1000,
+      // 23
+      1000, -300, 500, 1, 0, 0, 0, 1000,
+      // 24
+      1000, 300, -500, 1, 0, 0, 600, 0,
+      // 25
+      566.6257934570312, 300, 23.1280517578125, 0, 1, 0, 433.3742370605469,
+      523.1281127929688,
+      // 26
+      411.5867004394531, -66.51548767089844, -500, 0, 0, -1, 1411.586669921875,
+      366.5155029296875,
+      // 27
+      375.7498779296875, -4.444300651550293, -500, 0, 0, -1, 1375.7498779296875,
+      304.4443054199219,
+      // 28
+      346.7673034667969, 300, -500, 0, 1, 0, 653.2326049804688, 0,
+      // 29
+      -153.58984375, 300, -388.552490234375, 0, 1, 0, 1153.58984375,
+      111.447509765625,
+      // 30
+      199.9788818359375, 300, -500, 0, 1, 0, 800.0211791992188, 0,
+      // 31
+      -1000, 300, -500, 0, 1, 0, 2000, 0,
+      // 32
+      -153.58987426757812, 300, 44.22247314453125, 0, 1, 0, 1153.58984375,
+      544.2224731445312,
+      // 33
+      199.9788818359375, 300, -500, 0, 0, -1, 1199.9791259765625, 0,
+      // 34
+      521.6780395507812, -2.9542479515075684, -500, 0, 0, -1, 1521.677978515625,
+      302.9542541503906,
+      // 35
+      346.7673034667969, 300, -500, 0, 0, -1, 1346.767333984375, 0,
+      // 36
+      1000, 300, -500, 0, 0, -1, 2000, 0,
+      // 37
+      -1000, 300, 500, -1, 0, 0, 0, 1000,
+      // 38
+      -1000, 300, 500, 0, 0, 1, 0, 600,
+      // 39
+      -1000, 300, 500, 0, 1, 0, 2000, 1000,
+      // 40
+      -153.58985900878906, 300, 500, 0, 0, 1, 846.4102172851562, 600,
+      // 41
+      88.46627807617188, -253.06915283203125, 500, 0, 0, 1, 1088.4664306640625,
+      46.93084716796875,
+      // 42
+      -153.58985900878906, 300, 500, 0, 1, 0, 1153.58984375, 1000,
+      // 43
+      7.1796698570251465, -300, 500, 0, 0, 1, 1007.1797485351562, 0,
+      // 44
+      1000, -300, -500, 0, -1, 0, 2000, 0,
+      // 45
+      1000, -300, 500, 0, -1, 0, 2000, 1000,
+      // 46
+      7.1796698570251465, -300, 500, 0, -1, 0, 1007.1796264648438, 1000,
+      // 47
+      403.5837097167969, 300, 500, 0, 1, 0, 596.4163208007812, 1000,
+      // 48
+      1000, -300, -500, 1, 0, 0, 0, 0,
+      // 49
+      492.3005676269531, -19.915321350097656, -500, 0, 0, -1, 1492.300537109375,
+      319.91534423828125,
+      // 50
+      411.5867004394531, -66.51548767089844, -500, -0.5, 0.8660253882408142, 0,
+      880.5439453125, 0,
+      // 51
+      7.179656982421875, -300, -330.03717041015625, -0.5000000596046448,
+      0.866025447845459, 0, 383.6058654785156, 0,
+      // 52
+      492.3005676269531, -19.915321350097656, -500, -0.5, 0.8660253882408142, 0,
+      968.1235961914062, 31.876384735107422,
+      // 53
+      7.1796698570251465, -300, 500, -0.5000000596046448, 0.866025447845459, 0,
+      99.71644592285156, 779.979736328125,
+      // 54
+      88.46627807617188, -253.06915283203125, 500, -0.5, 0.8660253882408142, 0,
+      187.91758728027344, 812.0823974609375,
+      // 55
+      -153.58985900878906, 300, 500, 0.5, -0.866025447845459, 0,
+      749.2095947265625, 834.9661865234375,
+      // 56
+      -384.7195129394531, 166.55722045898438, 500, 0.5000000596046448,
+      -0.866025447845459, 0, 1000, 743.6859741210938,
+      // 57
+      -153.58987426757812, 300, 44.22247314453125, 0.5, -0.8660253882408142, 0,
+      593.3245239257812, 406.6754455566406,
+      // 58
+      564.2904052734375, 21.64801025390625, 500, -0.5, 0.866025447845459, 0,
+      704.217041015625, 1000.0000610351562,
+      // 59
+      -604.9979248046875, 39.37942886352539, -198.83624267578125, 0.5,
+      -0.8660253882408142, 0, 1000, 0,
+      // 60
+      199.9788818359375, 300, -500, 0.29619815945625305, 0.1710100919008255,
+      0.9396927356719971, 880.5438842773438, 176.7843475341797,
+      // 61
+      -153.58984375, 300, -388.552490234375, 0.29619815945625305,
+      0.1710100919008255, 0.9396927356719971, 554.6932373046875, 0,
+      // 62
+      375.7498779296875, -4.444300651550293, -500, 0.29619812965393066,
+      0.1710100919008255, 0.9396926760673523, 880.5438842773438,
+      528.3263549804688,
+      // 63
+      566.6257934570312, 300, 23.1280517578125, -0.8137977123260498,
+      -0.46984636783599854, 0.342020183801651, 239.89218139648438, 600.1796875,
+      // 64
+      346.7673034667969, 300, -500, -0.8137977719306946, -0.46984633803367615,
+      0.342020183801651, 349.8214111328125, 43.478458404541016,
+      // 65
+      521.6780395507812, -2.9542479515075684, -500, -0.8137977719306946,
+      -0.46984633803367615, 0.342020183801651, 0, 43.478458404541016,
+      // 66
+      804.9979248046875, 160.62057495117188, 398.83624267578125, -0.5,
+      0.8660253882408142, 0, 1000, 1000,
+      // 67
+      521.6780395507812, -2.9542479515075684, -500, -0.5, 0.8660253882408142, 0,
+      1000, 43.47837829589844,
+      // 68
+      -153.58984375, 300, -388.552490234375, 0.5, -0.8660253882408142, 0,
+      445.3067626953125, 0,
+      // 69
+      -604.9979248046875, 39.37942886352539, -198.83624267578125,
+      0.29619815945625305, 0.1710100919008255, 0.9396927356719971, 0, 0,
+      // 70
+      804.9979248046875, 160.62057495117188, 398.83624267578125,
+      -0.813797652721405, -0.46984630823135376, 0.3420201539993286, 0, 1000,
+      // 71
+      -161.6136932373047, -219.87335205078125, 500, 0.8137977123260498,
+      0.46984630823135376, -0.3420201539993286, 446.21160888671875,
+      743.68603515625,
+      // 72
+      -604.9979248046875, 39.37942886352539, -198.83624267578125,
+      0.813797652721405, 0.46984630823135376, -0.3420201539993286, 0, 0,
+      // 73
+      -384.7195129394531, 166.55722045898438, 500, 0.8137977123260498,
+      0.46984636783599854, -0.342020183801651, 0, 743.6859741210938,
+      // 74
+      -115.35255432128906, -300, 500, 0.813797652721405, 0.46984633803367615,
+      -0.3420201539993286, 538.73388671875, 743.68603515625,
+      // 75
+      -409.0570983886719, -300, -198.83624267578125, 0.813797652721405,
+      0.46984630823135376, -0.3420201539993286, 391.8816223144531, 0,
+      // 76
+      7.179656982421875, -300, -330.03717041015625, 0.29619815945625305,
+      0.1710100919008255, 0.9396927356719971, 383.6058654785156, 600,
+      // 77
+      564.2904052734375, 21.64801025390625, 500, -0.29619812965393066,
+      -0.1710100919008255, -0.9396926164627075, 704.2169189453125,
+      -0.000030517578125,
+      // 78
+      403.5837097167969, 300, 500, -0.29619812965393066, -0.1710100919008255,
+      -0.9396926164627075, 704.2169799804688, 321.4132385253906,
+      // 79
+      724.5271606445312, 300, 398.83624267578125, -0.29619815945625305,
+      -0.1710100919008255, -0.9396926760673523, 1000, 160.94149780273438,
+      // 80
+      804.9979248046875, 160.62057495117188, 398.83624267578125,
+      -0.29619815945625305, -0.1710100919008255, -0.9396927356719971, 1000, 0,
+      // 81
+      -409.0570983886719, -300, -198.83624267578125, 0.29619815945625305,
+      0.1710100919008255, 0.9396927356719971, 0, 391.88165283203125,
+      // 82
+      724.5271606445312, 300, 398.83624267578125, -0.813797652721405,
+      -0.46984630823135376, 0.342020183801651, 160.94149780273438,
+      1000.0000610351562,
+      // 83
+      411.5867004394531, -66.51548767089844, -500, 0.29619815945625305,
+      0.1710100769996643, 0.9396926164627075, 880.5440063476562,
+      600.0000610351562};
+  a.triVerts = {// 0
+                0, 1, 2,
+                // 1
+                3, 4, 5,
+                // 2
+                6, 0, 2,
+                // 3
+                7, 8, 9,
+                // 4
+                10, 11, 12,
+                // 5
+                13, 14, 15,
+                // 6
+                0, 16, 1,
+                // 7
+                17, 18, 19,
+                // 8
+                9, 20, 7,
+                // 9
+                18, 21, 19,
+                // 10
+                22, 23, 24,
+                // 11
+                14, 25, 15,
+                // 12
+                26, 12, 27,
+                // 13
+                14, 28, 25,
+                // 14
+                29, 30, 31,
+                // 15
+                29, 31, 32,
+                // 16
+                12, 33, 27,
+                // 17
+                34, 35, 36,
+                // 18
+                5, 4, 37,
+                // 19
+                17, 38, 18,
+                // 20
+                31, 39, 32,
+                // 21
+                40, 38, 17,
+                // 22
+                9, 41, 20,
+                // 23
+                39, 42, 32,
+                // 24
+                41, 43, 20,
+                // 25
+                6, 2, 44,
+                // 26
+                6, 45, 46,
+                // 27
+                26, 10, 12,
+                // 28
+                47, 13, 15,
+                // 29
+                48, 24, 23,
+                // 30
+                6, 44, 45,
+                // 31
+                26, 49, 10,
+                // 32
+                49, 34, 10,
+                // 33
+                34, 36, 10,
+                // 34
+                50, 51, 52,
+                // 35
+                51, 53, 54,
+                // 36
+                51, 54, 52,
+                // 37
+                55, 56, 57,
+                // 38
+                52, 54, 58,
+                // 39
+                59, 57, 56,
+                // 40
+                60, 61, 62,
+                // 41
+                63, 64, 65,
+                // 42
+                52, 66, 67,
+                // 43
+                59, 68, 57,
+                // 44
+                69, 62, 61,
+                // 45
+                65, 70, 63,
+                // 46
+                71, 72, 73,
+                // 47
+                52, 58, 66,
+                // 48
+                74, 72, 71,
+                // 49
+                74, 75, 72,
+                // 50
+                62, 69, 76,
+                // 51
+                77, 78, 79,
+                // 52
+                79, 80, 77,
+                // 53
+                69, 81, 76,
+                // 54
+                63, 70, 82,
+                // 55
+                76, 83, 62};
+  a.mergeFromVert = {3,  4,  11, 12, 13, 18, 21, 22, 23, 24, 31, 33, 35, 36,
+                     38, 39, 42, 44, 45, 46, 47, 48, 50, 51, 52, 53, 54, 55,
+                     56, 57, 58, 60, 61, 62, 63, 64, 65, 67, 68, 69, 70, 71,
+                     72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83};
+  a.mergeToVert = {2,  1,  2,  5,  7,  1,  16, 7,  20, 14, 5,  30, 28, 14,
+                   37, 37, 40, 10, 20, 43, 8,  10, 26, 6,  49, 43, 41, 40,
+                   17, 32, 9,  30, 29, 27, 25, 28, 34, 34, 29, 59, 66, 19,
+                   59, 17, 16, 0,  6,  9,  8,  15, 66, 0,  15, 26};
+
+  MeshGL b;
+  b.numProp = 8;
+  b.vertProperties = {// 0
+                      -1700, -600, -1000, -1, 0, 0, 1200, 0,
+                      // 1
+                      -1700, -600, 1000, -1, 0, 0, 1200, 2000,
+                      // 2
+                      -1700, 600, -1000, -1, 0, 0, 0, 0,
+                      // 3
+                      -1700, -600, -1000, 0, -1, 0, 0, 0,
+                      // 4
+                      300, -600, -1000, 0, -1, 0, 2000, 0,
+                      // 5
+                      -1700, -600, 1000, 0, -1, 0, 0, 2000,
+                      // 6
+                      -1700, -600, -1000, 0, 0, -1, 0, 1200,
+                      // 7
+                      -1700, 600, -1000, 0, 0, -1, 0, 0,
+                      // 8
+                      300, -600, -1000, 0, 0, -1, 2000, 1200,
+                      // 9
+                      -1700, -600, 1000, 0, 0, 1, 0, 0,
+                      // 10
+                      300, -600, 1000, 0, 0, 1, 2000, 0,
+                      // 11
+                      -1700, 600, 1000, 0, 0, 1, 0, 1200,
+                      // 12
+                      -1700, 600, 1000, -1, 0, 0, 0, 2000,
+                      // 13
+                      -1700, 600, -1000, 0, 1, 0, 2000, 0,
+                      // 14
+                      -1700, 600, 1000, 0, 1, 0, 2000, 2000,
+                      // 15
+                      300, 600, 1000, 0, 1, 0, 0, 2000,
+                      // 16
+                      300, -600, -1000, 1, 0, 0, 0, 0,
+                      // 17
+                      300, 600, -1000, 1, 0, 0, 1200, 0,
+                      // 18
+                      300, -600, 1000, 1, 0, 0, 0, 2000,
+                      // 19
+                      300, -600, 1000, 0, -1, 0, 2000, 2000,
+                      // 20
+                      300, 600, -1000, 0, 0, -1, 2000, 0,
+                      // 21
+                      300, 600, -1000, 0, 1, 0, 0, 0,
+                      // 22
+                      300, 600, 1000, 0, 0, 1, 2000, 1200,
+                      // 23
+                      300, 600, 1000, 1, 0, 0, 1200, 2000};
+  b.triVerts = {// 0
+                0, 1, 2,
+                // 1
+                3, 4, 5,
+                // 2
+                6, 7, 8,
+                // 3
+                9, 10, 11,
+                // 4
+                1, 12, 2,
+                // 5
+                13, 14, 15,
+                // 6
+                16, 17, 18,
+                // 7
+                4, 19, 5,
+                // 8
+                7, 20, 8,
+                // 9
+                21, 13, 15,
+                // 10
+                10, 22, 11,
+                // 11
+                17, 23, 18};
+  b.mergeFromVert = {3, 5, 6, 7, 8, 9, 12, 13, 14, 16, 18, 19, 20, 21, 22, 23};
+  b.mergeToVert = {0, 1, 0, 2, 4, 1, 11, 2, 11, 4, 10, 10, 17, 17, 15, 15};
+
+  a.runOriginalID = {Manifold::ReserveIDs(1)};
+  b.runOriginalID = {Manifold::ReserveIDs(1)};
+
+  Manifold aManifold(a);
+  Manifold bManifold(b);
+
+  auto aMinusB = aManifold - bManifold;
+
+  std::vector<MeshGL> meshList;
+  meshList.emplace_back(a);
+  meshList.emplace_back(b);
+
+  RelatedGL(aMinusB, meshList, false, false);
+}
+
+TEST(Boolean, CreatePropertiesSlow) {
+  Manifold a = Manifold::Sphere(10, 1024).SetProperties(
+      3, [](float* newprop, glm::vec3 pos, const float* old) {
+        for (int i = 0; i < 3; i++) newprop[i] = 0;
+      });
+  Manifold b = Manifold::Sphere(10, 1024).Translate({5, 0, 0});
+  Manifold result = a + b;
+  EXPECT_EQ(result.NumProp(), 3);
+}

test/manifold_test.cpp

@@ -232,6 +232,14 @@ TEST(Manifold, Revolve2) {
   EXPECT_NEAR(prop.surfaceArea, 96.0f * glm::pi<float>(), 1.0f);
 }
 
+TEST(Manifold, Revolve3) {
+  CrossSection circle = CrossSection::Circle(1, 32);
+  Manifold sphere = Manifold::Revolve(circle, 32);
+  auto prop = sphere.GetProperties();
+  EXPECT_NEAR(prop.volume, 4.0f / 3.0f * glm::pi<float>(), 0.1);
+  EXPECT_NEAR(prop.surfaceArea, 4 * glm::pi<float>(), 0.15);
+}
+
 TEST(Manifold, PartialRevolveOnYAxis) {
   Polygons polys = SquareHole(2.0f);
   Polygons offsetPolys = SquareHole(10.0f);
@@ -303,6 +311,23 @@ TEST(Manifold, Warp2) {
   EXPECT_NEAR(propBefore.volume, 321, 1);
 }
 
+TEST(Manifold, WarpBatch) {
+  Manifold shape1 =
+      Manifold::Cube({2, 3, 4}).Warp([](glm::vec3& v) { v.x += v.z * v.z; });
+  auto prop1 = shape1.GetProperties();
+
+  Manifold shape2 =
+      Manifold::Cube({2, 3, 4}).WarpBatch([](VecView<glm::vec3> vecs) {
+        for (glm::vec3& v : vecs) {
+          v.x += v.z * v.z;
+        }
+      });
+  auto prop2 = shape2.GetProperties();
+
+  EXPECT_EQ(prop1.volume, prop2.volume);
+  EXPECT_EQ(prop1.surfaceArea, prop2.surfaceArea);
+}
+
 TEST(Manifold, Smooth) {
   Manifold tet = Manifold::Tetrahedron();
   Manifold smooth = Manifold::Smooth(tet.GetMesh());
@@ -527,6 +552,14 @@ TEST(Manifold, Transform) {
   Identical(cube.GetMesh(), cube2.GetMesh());
 }
 
+TEST(Manifold, Slice) {
+  Manifold cube = Manifold::Cube();
+  CrossSection bottom = cube.Slice();
+  CrossSection top = cube.Slice(1);
+  EXPECT_EQ(bottom.Area(), 1);
+  EXPECT_EQ(top.Area(), 0);
+}
+
 TEST(Manifold, MeshRelation) {
   Mesh gyroidMesh = Gyroid();
   MeshGL gyroidMeshGL = WithIndexColors(gyroidMesh);

test/polygon_fuzz.cpp

@@ -0,0 +1,123 @@
+// Copyright 2023 The Manifold Authors.
+//
+// Licensed under the Apache License, Version 2.0 (the "License");
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+//
+//      http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+
+#include <atomic>
+#include <future>
+
+#include "fuzztest/fuzztest.h"
+#include "gtest/gtest.h"
+#include "polygon.h"
+
+using namespace fuzztest;
+
+void TriangulationNoCrash(
+    std::vector<std::vector<std::pair<float, float>>> input, float precision) {
+  if (precision < 0) precision = -1;
+  manifold::PolygonParams().intermediateChecks = true;
+  manifold::Polygons polys;
+  size_t size = 0;
+  for (const auto &simplePoly : input) {
+    polys.emplace_back();
+    for (const auto &p : simplePoly) {
+      polys.back().emplace_back(p.first, p.second);
+      size++;
+    }
+  }
+  std::atomic<pid_t> tid;
+  std::atomic<bool> faulted(true);
+  auto asyncFuture =
+      std::async(std::launch::async, [&polys, &faulted, &tid, precision]() {
+        tid.store(gettid());
+        try {
+          manifold::Triangulate(polys, precision);
+          faulted.store(false);
+        } catch (manifold::geometryErr e) {
+          // geometryErr is fine
+          faulted.store(false);
+        } catch (...) {
+          printf("got unexpected error\n");
+        }
+      });
+  if (asyncFuture.wait_for(std::chrono::milliseconds(
+          std::max<size_t>(size, 10000))) == std::future_status::timeout) {
+    printf("timeout after %ldms...\n", std::max<size_t>(size, 10000));
+    pthread_cancel(tid.load());
+  }
+
+  EXPECT_FALSE(faulted.load());
+}
+
+void TriangulationNoCrashRounded(
+    std::vector<std::vector<std::pair<float, float>>> input, float precision) {
+  TriangulationNoCrash(std::move(input), precision);
+}
+
+using Polygons = std::vector<std::vector<std::pair<float, float>>>;
+using TestCase = std::tuple<Polygons, float>;
+
+std::vector<TestCase> SeedProvider();
+
+auto PolygonDomain =
+    ContainerOf<Polygons>(ContainerOf<std::vector<std::pair<float, float>>>(
+                              PairOf(Finite<float>(), Finite<float>()))
+                              .WithMinSize(3)
+                              .WithMaxSize(10000))
+        .WithMinSize(1)
+        .WithMaxSize(10000);
+
+FUZZ_TEST(PolygonFuzz, TriangulationNoCrash)
+    .WithDomains(PolygonDomain, InRange<float>(-1.0, 0.1))
+    .WithSeeds(SeedProvider);
+
+FUZZ_TEST(PolygonFuzz, TriangulationNoCrashRounded)
+    .WithDomains(ReversibleMap(
+                     [](auto input) {
+                       for (auto &poly : input) {
+                         for (auto &pair : poly) {
+                           pair.first = std::round(pair.first);
+                           pair.second = std::round(pair.second);
+                         }
+                       }
+                       return input;
+                     },
+                     [](auto input) {
+                       return std::optional{std::tuple{input}};
+                     },
+                     PolygonDomain),
+                 InRange<float>(-1.0, 0.1))
+    .WithSeeds(SeedProvider);
+
+static std::vector<TestCase> TestCases;
+void TestPoly(Polygons polys, int _unused, float precision = -1.0f) {
+  TestCases.push_back({polys, precision});
+}
+
+// ----------------------------------------------------------------------------
+//                        UGLY HACK FOR SEED REUSE
+// ----------------------------------------------------------------------------
+
+#undef TEST
+#define TEST(_unused1, _unused2)
+
+std::vector<TestCase> SeedProvider() {
+  std::vector<TestCase> TestCases;
+  auto TestPoly = [&TestCases](Polygons polys, int _unused,
+                               float precision = -1.0f) {
+    TestCases.push_back({polys, precision});
+  };
+
+#include "polygon_corpus.cpp"
+
+  return TestCases;
+}

test/samples_test.cpp

@@ -182,9 +182,29 @@ TEST(Samples, Frame) {
 TEST(Samples, Bracelet) {
   Manifold bracelet = StretchyBracelet();
   CheckNormals(bracelet);
-  EXPECT_LE(bracelet.NumDegenerateTris(), 22);
+  EXPECT_EQ(bracelet.NumDegenerateTris(), 0);
   EXPECT_EQ(bracelet.Genus(), 1);
   CheckGL(bracelet);
+
+  CrossSection projection = bracelet.Project();
+  Rect rect = projection.Bounds();
+  Box box = bracelet.BoundingBox();
+  EXPECT_EQ(rect.min.x, box.min.x);
+  EXPECT_EQ(rect.min.y, box.min.y);
+  EXPECT_EQ(rect.max.x, box.max.x);
+  EXPECT_EQ(rect.max.y, box.max.y);
+  EXPECT_NEAR(projection.Area(), 649, 1);
+  EXPECT_EQ(projection.NumContour(), 2);
+  Manifold extrusion = Manifold::Extrude(projection, 1);
+  EXPECT_EQ(extrusion.NumDegenerateTris(), 0);
+  EXPECT_EQ(extrusion.Genus(), 1);
+
+  CrossSection slice = bracelet.Slice();
+  EXPECT_EQ(slice.NumContour(), 2);
+  EXPECT_NEAR(slice.Area(), 230.6, 0.1);
+  extrusion = Manifold::Extrude(slice, 1);
+  EXPECT_EQ(extrusion.Genus(), 1);
+
 #ifdef MANIFOLD_EXPORT
   if (options.exportModels) ExportMesh("bracelet.glb", bracelet.GetMesh(), {});
 #endif
@@ -202,9 +222,16 @@ TEST(Samples, GyroidModule) {
   const float precision = gyroid.Precision();
   EXPECT_NEAR(bounds.min.z, 0, precision);
   EXPECT_NEAR(bounds.max.z, size * glm::sqrt(2.0f), precision);
+
+  CrossSection slice = gyroid.Slice(5);
+  EXPECT_EQ(slice.NumContour(), 4);
+  EXPECT_NEAR(slice.Area(), 121.9, 0.1);
+  Manifold extrusion = Manifold::Extrude(slice, 1);
+  EXPECT_EQ(extrusion.Genus(), -3);
+
 #ifdef MANIFOLD_EXPORT
   if (options.exportModels)
-    ExportMesh("gyroidModule.gltf", gyroid.GetMesh(), {});
+    ExportMesh("gyroidModule.glb", gyroid.GetMesh(), {});
 #endif
 }
 
@@ -237,6 +264,18 @@ TEST(Samples, Sponge4) {
   EXPECT_EQ(cutSponge.first.Genus(), 13394);
   EXPECT_EQ(cutSponge.second.Genus(), 13394);
 
+  CrossSection projection = cutSponge.first.Project();
+  Rect rect = projection.Bounds();
+  Box box = cutSponge.first.BoundingBox();
+  EXPECT_EQ(rect.min.x, box.min.x);
+  EXPECT_EQ(rect.min.y, box.min.y);
+  EXPECT_EQ(rect.max.x, box.max.x);
+  EXPECT_EQ(rect.max.y, box.max.y);
+  EXPECT_NEAR(projection.Area(), 0.535, 0.001);
+  Manifold extrusion = Manifold::Extrude(projection, 1);
+  EXPECT_EQ(extrusion.NumDegenerateTris(), 0);
+  EXPECT_EQ(extrusion.Genus(), 502);
+
 #ifdef MANIFOLD_EXPORT
   if (options.exportModels) {
     ExportMesh("mengerHalf.glb", cutSponge.first.GetMesh(), {});
@@ -266,6 +305,28 @@ TEST(Samples, SelfIntersect) {
   res.GetMeshGL();  // test crash
   manifold::PolygonParams().processOverlaps = false;
 }
+
+TEST(Samples, GenericTwinBooleanTest7863) {
+  manifold::PolygonParams().processOverlaps = true;
+  std::string file = __FILE__;
+  std::string dir = file.substr(0, file.rfind('/'));
+  Manifold m1 = ImportMesh(dir + "/models/Generic_Twin_7863.1.t0_left.glb");
+  Manifold m2 = ImportMesh(dir + "/models/Generic_Twin_7863.1.t0_right.glb");
+  Manifold res = m1 + m2;  // Union
+  res.GetMeshGL();         // test crash
+  manifold::PolygonParams().processOverlaps = false;
+}
+
+TEST(Samples, Havocglass8Bool) {
+  manifold::PolygonParams().processOverlaps = true;
+  std::string file = __FILE__;
+  std::string dir = file.substr(0, file.rfind('/'));
+  Manifold m1 = ImportMesh(dir + "/models/Havocglass8_left.glb");
+  Manifold m2 = ImportMesh(dir + "/models/Havocglass8_right.glb");
+  Manifold res = m1 + m2;  // Union
+  res.GetMeshGL();         // test crash
+  manifold::PolygonParams().processOverlaps = false;
+}
 #endif
 
 TEST(Samples, CondensedMatter16) {

test/sdf_test.cpp

@@ -115,3 +115,21 @@ TEST(SDF, Resize) {
   EXPECT_EQ(layers.Status(), Manifold::Error::NoError);
   EXPECT_EQ(layers.Genus(), -8);
 }
+
+TEST(SDF, SineSurface) {
+  Mesh surface(LevelSet(
+      [](glm::vec3 p) {
+        float mid = glm::sin(p.x) + glm::sin(p.y);
+        return (p.z > mid - 0.5 && p.z < mid + 0.5) ? 1 : 0;
+      },
+      {glm::vec3(-4 * glm::pi<float>()), glm::vec3(4 * glm::pi<float>())}, 1));
+  Manifold smoothed = Manifold::Smooth(surface).Refine(2);
+
+  EXPECT_EQ(smoothed.Status(), Manifold::Error::NoError);
+  EXPECT_EQ(smoothed.Genus(), -2);
+
+#ifdef MANIFOLD_EXPORT
+  if (options.exportModels)
+    ExportMesh("sinesurface.glb", smoothed.GetMeshGL(), {});
+#endif
+}

test/test_main.cpp

@@ -373,7 +373,7 @@ void RelatedGL(const Manifold& out, const std::vector<MeshGL>& originals,
         glm::vec3 edges[3];
         for (int k : {0, 1, 2}) edges[k] = inTriPos[k] - outTriPos[j];
         const float volume = glm::dot(edges[0], glm::cross(edges[1], edges[2]));
-        ASSERT_LE(volume, area * 100 * out.Precision());
+        ASSERT_LE(volume, area * out.Precision());
 
         if (checkNormals) {
           glm::vec3 normal;
@@ -396,7 +396,7 @@ void RelatedGL(const Manifold& out, const std::vector<MeshGL>& originals,
             const float volumeP =
                 glm::dot(edgesP[0], glm::cross(edgesP[1], edgesP[2]));
 
-            ASSERT_LE(volumeP, area * 100 * out.Precision());
+            ASSERT_LE(volumeP, area * out.Precision());
           }
         }
       }