From f4badebd1a9a3b462e9800031f112b610450f7da Mon Sep 17 00:00:00 2001 From: Revar Desmera Date: Thu, 27 Aug 2026 10:58:29 -0700 Subject: [PATCH] levelset: cut cleanly where a 3D surface meets the box Reported from a real model: a gyroid built with levelset() came out visibly staircased along every edge where the sheet meets the bounding box, next to a smooth one from BOSL2's isosurface() of the same field. Manifold closes the mesh where the surface reaches the edge of the box it is given, and that closure follows the sample lattice. Anything touching the boundary gets a voxel staircase; a gyroid touches it everywhere. Fixed the way the 2D path already was: sample a box padded by 2*edge, let the ragged closure happen out there, then intersect with a cube at the bounds actually asked for, so the cut is a plane. The grid sampler now reads "outside" beyond its data rather than clamping, which would smear the boundary values outward and extend the surface instead of closing it. The tests did not catch this because every shape in them sat inside its box. The numeric form of the bug is sharp once looked for: a half-space should be exactly half the box at any resolution, and now is -- 32000 to the digit at edge 2, 1 and 0.5. Worth recording: the 2D path got this exact fix first, and 3D was not checked at the time. A boundary bug found in one dimension is worth looking for in the other. Three regression tests. 1078 passing under both engines. --- pyproject.toml | 2 +- src/builtins/topology.cpp | 31 ++++++++++++++++++++++++-- tests/test_booleans.cpp | 46 +++++++++++++++++++++++++++++++++++++++ 3 files changed, 76 insertions(+), 3 deletions(-) diff --git a/pyproject.toml b/pyproject.toml index d6692b7..59c9231 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -4,7 +4,7 @@ build-backend = "scikit_build_core.build" [project] name = "openscad_cpp_evaluator" -version = "0.51.0" +version = "0.51.1" description = "C++ OpenSCAD evaluator with Python bindings" readme = "README.md" requires-python = ">=3.12" diff --git a/src/builtins/topology.cpp b/src/builtins/topology.cpp index d620da7..f78f5b6 100644 --- a/src/builtins/topology.cpp +++ b/src/builtins/topology.cpp @@ -970,6 +970,19 @@ std::vector generateLevelSet(Evaluator& ev, const CSGParams& params if (i0[a] > n[a] - 2) i0[a] = n[a] - 2; t[a] = g[a] - static_cast(i0[a]); } + // Outside the sampled block the answer is "outside", falling away + // with distance. Clamping instead (which is what the interpolation + // above does on its own) would smear the boundary values outward and + // extend the surface rather than close it. + double outside = 0.0; + for (int a = 0; a < 3; ++a) { + const double maxg = static_cast(n[a] - 1); + const double raw = (pos[a] - origin[a]) / spacing[a]; + if (raw < 0.0) outside = std::max(outside, -raw); + else if (raw > maxg) outside = std::max(outside, raw - maxg); + } + if (outside > 0.0) return invert ? (1.0 + outside) : -(1.0 + outside); + double acc = 0.0; for (int c = 0; c < 8; ++c) { const size_t di = static_cast(c & 1), dj = static_cast((c >> 1) & 1), @@ -980,8 +993,18 @@ std::vector generateLevelSet(Evaluator& ev, const CSGParams& params return bandDistance(acc, isoLo, isoHi, invert); }; - manifold::Box bounds(manifold::vec3(origin[0], origin[1], origin[2]), - manifold::vec3((*hi)[0], (*hi)[1], (*hi)[2])); + // Manifold closes the mesh where the surface meets the edge of the box it + // is given, and that closure follows the sample lattice -- a visible + // staircase on anything that reaches the boundary, which a gyroid does + // everywhere. So sample a PADDED box, let the ragged closure happen out + // there, and cut back to the box actually asked for. The cut is then a + // clean plane. + // + // Exactly the fix already applied to the 2D path (clipToBounds); it was + // not applied here at the time because the 2D error was the one measured. + const double pad = 2.0 * edge; + manifold::Box bounds(manifold::vec3(origin[0] - pad, origin[1] - pad, origin[2] - pad), + manifold::vec3((*hi)[0] + pad, (*hi)[1] + pad, (*hi)[2] + pad)); // tolerance -1: a positive value makes Manifold do EXTRA evaluations per // output vertex to snap nearer the true surface. Against a fixed grid // those only re-interpolate data already used -- cost, no information. @@ -996,6 +1019,10 @@ std::vector generateLevelSet(Evaluator& ev, const CSGParams& params : manifold::Manifold::LevelSet(sampleGrid, bounds, edge, 0.0, -1.0, /*canParallel=*/true); if (solid.IsEmpty()) return {}; + solid = solid ^ manifold::Manifold::Cube( + manifold::vec3((*hi)[0] - origin[0], (*hi)[1] - origin[1], (*hi)[2] - origin[2])) + .Translate(manifold::vec3(origin[0], origin[1], origin[2])); + if (solid.IsEmpty()) return {}; ColoredBody b; b.body = std::move(solid); std::vector out; diff --git a/tests/test_booleans.cpp b/tests/test_booleans.cpp index c4bc9e6..86cd861 100644 --- a/tests/test_booleans.cpp +++ b/tests/test_booleans.cpp @@ -1143,3 +1143,49 @@ TEST(LevelSet2d, A3dArrayWithA2dBoundsWarns) { EXPECT_FALSE(levelsetWarnings("levelset([[[1,2],[3,4]],[[5,6],[7,8]]], bounds=[[0,0],[1,1]]);") .empty()); } + +TEST(LevelSet, ASurfaceMeetingTheBoxIsCutCleanly) { + // Manifold closes the mesh where the surface reaches the edge of the box + // it is given, and that closure follows the SAMPLE LATTICE -- a visible + // staircase on anything that touches the boundary. A gyroid touches it + // everywhere, which is how this was found: the shape looked ragged next + // to BOSL2's isosurface() of the same field. + // + // Sampling a padded box and cutting back makes the cut a plane. The + // numeric form of that: a half-space is EXACTLY half the box, at every + // resolution. Before the fix it was resolution-dependent. + // + // The 2D path had this fixed already (clipToBounds); 3D did not, because + // 2D was the case that got measured at the time. + for (const char* edge : {"2", "1", "0.5"}) { + Evaluated e = evalSrc(std::string("levelset(function(x,y,z) x, " + "bounds=[[-20,-20,-20],[20,20,20]], " + "isovalue=[-1e18,0], edge=") + + edge + ");"); + ASSERT_EQ(e.bodies.size(), 1u) << "edge=" << edge; + EXPECT_NEAR(soleBody(e).Volume(), 32000.0, 1e-6) << "edge=" << edge; + } +} + +TEST(LevelSet, PaddingDoesNotDisturbAShapeInsideTheBox) { + // The padding is sampled outside the requested bounds, so a shape that + // never reaches the boundary must come out exactly as before. + Evaluated e = evalSrc(sphereFieldSrc(50, 30) + + "levelset(f, bounds=[[-30,-30,-30],[30,30,30]], isovalue=20);"); + ASSERT_EQ(e.bodies.size(), 1u); + const double analytic = 4.0 / 3.0 * 3.14159265358979 * 8000; + EXPECT_NEAR(soleBody(e).Volume(), analytic, 0.01 * analytic); + EXPECT_EQ(soleBody(e).Genus(), 0); +} + +TEST(LevelSet, AGridFieldIsAlsoCutCleanlyAtTheBox) { + // The grid intake cannot sample beyond its data, so the sampler reads + // "outside" past the block rather than clamping and smearing the edge + // values outward. Same clean cut as the function form. + Evaluated e = evalSrc( + "N = 41;\nfunction co(t) = -20 + 40*t/(N-1);\n" + "f = [for (i=[0:N-1]) [for (j=[0:N-1]) [for (k=[0:N-1]) co(i) ]]];\n" + "levelset(f, bounds=[[-20,-20,-20],[20,20,20]], isovalue=[-1e18,0]);"); + ASSERT_EQ(e.bodies.size(), 1u); + EXPECT_NEAR(soleBody(e).Volume(), 32000.0, 40.0); // one sample layer of slack +}