fix Screw3D SDF accuracy: octree holes, taper slope, asymmetric profiles

[snowbldr]

Three bugs in Screw3D.Evaluate caused the octree marching cubes renderer to produce meshes with holes (non-manifold boundary edges). The uniform renderer was unaffected because it evaluates every cube unconditionally. One additional floating-point bug in Line2.IntersectLine was found along the way and fixed.

Bug 1: Helical stretch overestimates distance (all thread types)

Root cause. Screw3D maps 3D coordinates to a 2D thread profile via a helical unwrap: the angular position θ and axial position z are combined into a single linear coordinate along the pitch. This mapping is non-isometric — a unit step in 3D maps to a step of size √(1 + (lead/(2πr))²) in the unwrapped 2D space. The 2D thread SDF therefore returns a value larger than the true 3D distance.

The octree renderer's isEmpty check evaluates the SDF at cube centers and skips any cube where |d| ≥ half_diagonal. The overestimated distance causes it to skip cubes that actually contain the surface, producing holes. Steeper helices (high pitch-to-radius ratio, multi-start threads) have a larger stretch factor and are affected more.

Fix. Divide the final SDF value by the maximum singular value of the full Jacobian of the (r, θ, z) → (profileX, profileY) mapping:

J = | 0      k      1      |   where k = lead/(2πr)
    | 1      0      tan(t) |   and t = taper angle

JJᵀ = | 1+k²     tan(t)     |
      | tan(t)   1+tan²(t)  |

σ_max = √((a + c + √((a-c)² + 4·tan²(t))) / 2)
where a = 1+k², c = 1+tan²(t)

When taper=0 this simplifies to √(1+k²). Dividing by σ_max restores the Lipschitz-1 property: |SDF(p) - SDF(q)| ≤ |p - q|.

The correction is applied to max(d0, d1) (the intersection of thread profile and length constraint) rather than to d0 alone, so that end cap surfaces (where d=0) remain at d=0 after division.

The stretch factor varies with r (larger near the axis). The octree isEmpty check evaluates at cube centers, but the nearest surface may be at a smaller r where the stretch is larger. Using rEff = max(r - threadDepth, r/2) gives a conservative bound: the nearest thread surface can be at most threadDepth closer to the axis than the evaluation point.

Bug 2: Taper slope uses atan instead of tan

Root cause. The Evaluate method computed the taper slope as math.Atan(s.taper) instead of math.Tan(s.taper). For a taper angle in radians, tan(angle) gives the slope (rise/run), while atan(angle) gives the angle whose tangent is the input — a different operation entirely.

For small angles like NPT (~1.8°), atan(x) ≈ tan(x) ≈ x, so the error was invisible (0.06%). But at larger taper angles the error is significant:

taper=30°  tan=0.5774  atan=0.4823  (16.5% error)

This caused the thread crest to be at the wrong radius for tapered screws, and the incorrect slope fed into the Jacobian stretch correction, leaving residual holes in tapered screw meshes.

Fix. Replace math.Atan(s.taper) with math.Tan(s.taper). The value is precomputed in the constructor as tanTaper and tan2Taper fields.

Bug 3: SawTooth wrap discontinuity for asymmetric profiles (buttress)

Root cause. The helical coordinate is folded into one pitch period via SawTooth(z, pitch), mapping to [-pitch/2, +pitch/2]. For symmetric profiles (ISO, ACME), the SDF is continuous across the wrap boundary because the profile is mirror-symmetric: SDF(+pitch/2) = SDF(-pitch/2).

For asymmetric profiles (ANSI buttress, plastic buttress), the profile shape differs at the left vs right edge. The SDF has a discontinuity at the wrap boundary — for buttress threads, up to 0.76mm jump at the boundary vs 0.00mm for ISO. This discontinuity causes the octree to skip cubes straddling the boundary, producing holes.

Fix. Evaluate the 2D thread profile at the current period and both adjacent periods, taking the minimum (union semantics):

d0 = min(thread.Evaluate(p0),
         thread.Evaluate(p0 - pitch),
         thread.Evaluate(p0 + pitch))

Each thread tooth is geometrically identical, so the correct SDF for a periodic thread pattern is the union of all period copies. The minimum of signed distances gives union semantics: inside any copy → negative, outside all copies → positive with distance to the nearest surface.

For symmetric profiles this is a no-op: the adjacent-period values are always >= the current-period value at the boundary. For asymmetric profiles it eliminates the discontinuity by ensuring the SDF always reflects the nearest thread tooth regardless of which period SawTooth maps to.

Bug 4: Line2.IntersectLine floating point comparison

Root cause. The collinear intersection path in IntersectLine compared interval endpoints with exact equality (x[0] == x[1]) to decide whether the overlap is a single point or a range. Floating point arithmetic can produce x[0]=0, x[1]=1e-16 for what should be a single-point intersection, causing IntersectLine to return 2 nearly-identical points instead of 1.

Fix. Replace exact equality with tolerance comparison: math.Abs(x[0]-x[1]) <= tolerance (1e-9).

New files

render/mesh.go — Mesh analysis utilities:

• CollectTriangles: renders an SDF3 and returns all triangles
• CountBoundaryEdges: counts edges shared by only one triangle
• IsWatertight: convenience check for zero boundary edges
• MaxZ: finds maximum Z vertex coordinate

sdf/screw_test.go — Unit test for taper slope (tan vs atan validation)

render/screw_test.go — 82 subtests across 3 test functions:

• Test_Screw_Watertight_Straight (33 configs)
• Test_Screw_Watertight_Tapered (16 configs)
• Test_Screw_EndCap_Position (33 configs)

examples/screw_assortment/ — Example rendering 44 screw configs at 300 cells (the dominant cell count across sdfx examples) covering all thread profiles, start counts, handedness, taper angles, and extreme dimensions. Outputs screws.stl with SHA1SUM verification. Layout spaces screws by max-neighbor radius so oversized configs don't crowd small ones. Low-resolution stress configs (25/50 cells) live in render/screw_test.go rather than here.

tools/stldiff/ — Cross-branch STL regression checker. main.go parses binary STLs and reports IDENTICAL (canonical hash match), MINOR (tiny float drift), or MATERIAL (real geometry change). run.sh creates detached worktrees at BASE and HEAD refs, renders every example on both, compares every produced STL, and prints a summary. Invocable via make stldiff from the repo root; override refs with make stldiff BASE=x HEAD=y.

Validation for this PR (./tools/stldiff/run.sh master HEAD, 197 STLs across 76 examples): 183 IDENTICAL, 6 MINOR (same triangle counts, sub-micron bbox drift), 8 MATERIAL — all 8 MATERIAL entries are in screw-using examples (3dp_nutbolt, bolt_container, gas_cap, nutsandbolts, tapers, test), as expected from the screw SDF fix. No unexpected geometry changes in any non-screw example.

Regenerated all example SHA1SUMs. Most were stale pre-existing drift (byte-level noise from sha1tool.py file ordering and renderer nondeterminism); only the 8 MATERIAL entries above reflect genuine screw SDF output changes.

[deadsy]

Looks like some good work. The screw distance issue has been bothering me for some time, I just didn't know how to fix it. :-) Instead I just took the hit on rendering performance by using the uniform renderer, but it's always felt a bit dirty.

[snowbldr]

Yeah honestly I had a lot of help from claude opus, it's actually the one that figured out the jacobian math bit. I had discovered a rudimentary fix by wrapping the SDF3, overwriting the evaluate function and multiplying by 0.7 or so, and that tricked the octree renderer into not skipping the cubes that contained what was becoming holes. Then I was reading through issues and PRs in the repo and saw you had mentioned that there was a bug in the screw logic, which is what lead me down this path to try and fix it.

I really wanted to make sure it was as correct as possible too hence all the testing and the stldiff tool to help verify correctness. The make target for that should help ensure that any future changes are stable as well.

I've been working on a modeling project with cadquery for the last year or so and finally just gave up after I was trying to get two helixes to intersect, that apparently makes open cascade very angry. This SDF based system you've got here is much nicer.

Thanks for taking a look at this, and creating this tool, I really appreciate it.

[deadsy]

Some regression issues: Probably with the buttress thread code.

examples/gas_cap

The internal threads have gone away (maybe they are buried in the wall...) after the PR has been applied.

[deadsy]

test issues:

jasonh@crun:~/personal/sdfx/render$ go test 
--- FAIL: Test_Screw_Watertight_Tapered (4.36s)
    --- FAIL: Test_Screw_Watertight_Tapered/taper_30 (0.28s)
        screw_test.go:141: octree mesh has 23 boundary edges (want 0 for watertight)
        screw_test.go:143: 84787 tris, 23 boundary edges
    --- FAIL: Test_Screw_Watertight_Tapered/taper_30_coarse (0.23s)
        screw_test.go:141: octree mesh has 54 boundary edges (want 0 for watertight)
        screw_test.go:143: 71136 tris, 54 boundary edges
    --- FAIL: Test_Screw_Watertight_Tapered/taper_30_left (0.46s)
        screw_test.go:141: octree mesh has 23 boundary edges (want 0 for watertight)
        screw_test.go:143: 84787 tris, 23 boundary edges
    --- FAIL: Test_Screw_Watertight_Tapered/taper_30_buttress (0.24s)
        screw_test.go:141: octree mesh has 4 boundary edges (want 0 for watertight)
        screw_test.go:143: 75594 tris, 4 boundary edges
    --- FAIL: Test_Screw_Watertight_Tapered/taper_30_left_buttress (0.26s)
        screw_test.go:141: octree mesh has 4 boundary edges (want 0 for watertight)
        screw_test.go:143: 75594 tris, 4 boundary edges
FAIL
exit status 1
FAIL    github.com/deadsy/sdfx/render   28.511s

[snowbldr]

Ah, thanks for double checking! I'll get it fixed an push an update.

[snowbldr]

gas cap is fixed

Tracked down the regression. The buttress threads disappear because the original "Bug 3" fix had an unintended side-effect on asymmetric profiles.

Quick background. Thread profiles describe one repeating unit of the thread in 2D. Screw3D sweeps that shape into a helix, wrapping the position around every pitch units. For the wrap to be smooth, what's just past +pitch/2 must match what's just past -pitch/2. Symmetric profiles (ISO, ACME) get this for
free — they mirror around the center. The original buttress profile is asymmetric and cheated: it padded both sides of the central tooth with flat crest, even though the real periodic shape has a slanted flank descending into the next valley.

  // before — flat-crest filler on both sides of the central tooth
  tp.Add(pitch, 0)
  tp.Add(pitch, radius)               // ← flat at full radius
  tp.Add(hp-((h0-h1)*t1), radius)
  tp.Add(t0*h0-hp, radius-h1).Smooth(0.0714*pitch, 5)
  tp.Add((h0-h1)*t0-hp, radius)
  tp.Add(-pitch, radius)              // ← flat at full radius
  tp.Add(-pitch, 0)

That's harmless for the uniform renderer (which only cares where the surface is). It produces holes on the octree renderer, which skips cubes that look "deep inside" the material — and the flat-crest mismatch makes points near the wrap look deeper inside than they actually are.

Bug 3 in the original commit tried to fix the holes by overlaying three shifted copies of the profile and keeping the closer one:

  // the workaround that broke buttress
  d0 := s.thread.Evaluate(p0)
  dL := s.thread.Evaluate(v2.Vec{X: p0.X - s.pitch, Y: p0.Y})
  dR := s.thread.Evaluate(v2.Vec{X: p0.X + s.pitch, Y: p0.Y})
  d0 = math.Min(d0, math.Min(dL, dR))

For symmetric profiles this changed nothing. For buttress, the flat-crest filler of one copy landed directly on the central valley of the next copy and filled it in — that's the gas cap's missing threads.

The actual fix. Replace the flat-crest filler with the actual flanks the periodic shape has there. The wrap matches automatically and the three-copy workaround in Evaluate isn't needed (and is removed):

  // after — both sides describe the real next/previous tooth
  xV  := t0*h0 - hp           // valley root
  x7  := hp - (h0-h1)*t1      // 7° flank top
  x45 := (h0-h1)*t0 - hp      // 45° flank top
  yEdge := radius + x45       // y on 45° flank where it crosses ±pitch

  tp.Add(pitch, 0)
  tp.Add(pitch, yEdge)                                  // flank crosses the wrap
  tp.Add(x45+pitch, radius)                             // next period's crest
  tp.Add(x7, radius)
  tp.Add(xV, radius-h1).Smooth(0.0714*pitch, 5)         // central valley
  tp.Add(x45, radius)
  tp.Add(x7-pitch, radius)
  tp.Add(xV-pitch, radius-h1).Smooth(0.0714*pitch, 5)   // previous period's valley
  tp.Add(-pitch, yEdge)                                 // mirror flank crosses the wrap
  tp.Add(-pitch, 0)

The original buttress polygon was correct for the renderer that existed when it was written (uniform marching cubes — only the surface position matters). The octree renderer added later introduced a new constraint — what's just past +pitch/2 must match what's just past -pitch/2 — that the old polygon
happens not to satisfy. This update brings buttress in line with that constraint.

New runtime check. To prevent anyone (us or a future contributor writing a custom asymmetric profile) from re-introducing this bug, Screw3D now validates the contract at construction time and fails fast with a useful error instead of producing silent octree holes:

if err := checkThreadWrapContinuous(thread, pitch); err != nil {
return nil, err
}

A bad profile gets rejected with a message pointing at the canonical example to copy from:

thread profile is discontinuous at the SawTooth wrap boundary x=±pitch/2
(pitch=2): SDF(+pitch/2,4.5)=-0.301 vs SDF(-pitch/2,4.5)=-0.5, |Δ|=0.199.
Asymmetric profiles must extend past ±pitch/2 with the actual periodic
continuation; see ANSIButtressThread for an example.

Tests:

• Test_Buttress_WrapContinuity — both buttress profiles satisfy the contract.
• Test_Screw3D_RejectsDiscontinuousProfile — a deliberately bad asymmetric polygon (the original buttress shape) gets refused at construction.
• Test_Screw_Watertight_* (existing) — every buttress config renders watertight on the octree at 25/100/300 cells.

Untouched: ISO and ACME profiles. Their symmetry already satisfies the wrap contract.

Visual: gas_cap (uniform, 200 cells) shows proper threads again; screw_assortment buttress configs render watertight on the octree renderer.

[deadsy]

Re: The SHA checks...

I've noticed that different CPUs create different results. Specifically amd64 vs aarch64. I would have thought that IEEE754 would have led to consistency, but maybe some different choices are made wrt rounding mode.

Bottom line: If you run this stuff on a Mac with an ARM you get different results than with an x86 PC.

[deadsy]

A few fails...

jasonh@eccles:~/work/go/src/github.com/deadsy/sdfx/examples/screw_assortment$ ./screw_assortment 
M10x2                        r=  5.0 p=2.0 taper=  0.0° starts=  1 → 671832 tris,    0 boundary edges [PASS]
M5x3                         r=  2.5 p=3.0 taper=  0.0° starts=  1 → 494108 tris,    0 boundary edges [PASS]
coarse_M10x5                 r=  5.0 p=5.0 taper=  0.0° starts=  1 → 533924 tris,    0 boundary edges [PASS]
steep_M3x3                   r=  1.5 p=3.0 taper=  0.0° starts=  1 → 216632 tris,    0 boundary edges [PASS]
extreme_M2x3                 r=  1.0 p=3.0 taper=  0.0° starts=  1 → 192664 tris,    0 boundary edges [PASS]
short_1pitch                 r=  5.0 p=2.0 taper=  0.0° starts=  1 → 423712 tris,    0 boundary edges [PASS]
dual_start                   r=  5.0 p=2.0 taper=  0.0° starts=  2 → 676328 tris,    0 boundary edges [PASS]
triple_start                 r=  5.0 p=2.0 taper=  0.0° starts=  3 → 684572 tris,    0 boundary edges [PASS]
multi8                       r=  5.0 p=2.0 taper=  0.0° starts=  8 → 757872 tris,    0 boundary edges [PASS]
multi16                      r=  5.0 p=2.0 taper=  0.0° starts= 16 → 932728 tris,    0 boundary edges [PASS]
multi8_coarse                r=  5.0 p=5.0 taper=  0.0° starts=  8 → 918836 tris,    0 boundary edges [PASS]
left_M10x2                   r=  5.0 p=2.0 taper=  0.0° starts= -1 → 671832 tris,    0 boundary edges [PASS]
left_8start                  r=  5.0 p=2.0 taper=  0.0° starts= -8 → 757872 tris,    0 boundary edges [PASS]
left_steep_M3x3              r=  1.5 p=3.0 taper=  0.0° starts= -1 → 216632 tris,    0 boundary edges [PASS]
internal_M10x2               r=  5.0 p=2.0 taper=  0.0° starts=  1 → 664756 tris,    0 boundary edges [PASS]
internal_M5x3                r=  2.5 p=3.0 taper=  0.0° starts=  1 → 507728 tris,    0 boundary edges [PASS]
acme_M10x2                   r=  5.0 p=2.0 taper=  0.0° starts=  1 → 638476 tris,    0 boundary edges [PASS]
acme_steep_M5x3              r=  2.5 p=3.0 taper=  0.0° starts=  1 → 519000 tris,    0 boundary edges [PASS]
buttress_M10x2               r=  5.0 p=2.0 taper=  0.0° starts=  1 → 698468 tris,    0 boundary edges [PASS]
buttress_steep               r=  2.5 p=3.0 taper=  0.0° starts=  1 → 523804 tris,    0 boundary edges [PASS]
buttress_left                r=  5.0 p=2.0 taper=  0.0° starts= -1 → 698320 tris,    0 boundary edges [PASS]
buttress_multi4              r=  5.0 p=2.0 taper=  0.0° starts=  4 → 740128 tris,    0 boundary edges [PASS]
plastic_butt_M10             r=  5.0 p=2.0 taper=  0.0° starts=  1 → 661164 tris,    0 boundary edges [PASS]
plastic_butt_left            r=  5.0 p=2.0 taper=  0.0° starts= -1 → 661184 tris,    0 boundary edges [PASS]
plastic_butt_multi4          r=  5.0 p=2.0 taper=  0.0° starts=  4 → 692616 tris,    0 boundary edges [PASS]
fine_M20x0.5                 r= 10.0 p=0.5 taper=  0.0° starts=  1 → 1629472 tris,    0 boundary edges [PASS]
sub_pitch_len                r=  5.0 p=3.0 taper=  0.0° starts=  1 → 368632 tris,    0 boundary edges [PASS]
taper_1.8_NPT                r=  5.0 p=2.0 taper=  1.8° starts=  1 → 679344 tris,    0 boundary edges [PASS]
taper_5                      r=  5.0 p=2.0 taper=  5.0° starts=  1 → 693400 tris,    0 boundary edges [PASS]
taper_15                     r=  5.0 p=2.0 taper= 15.0° starts=  1 → 751496 tris,    0 boundary edges [PASS]
taper_30                     r=  5.0 p=2.0 taper= 30.0° starts=  1 → 749916 tris,    0 boundary edges [PASS]
taper_45                     r=  5.0 p=2.0 taper= 45.0° starts=  1 → 549978 tris,   32 boundary edges [FAIL]
taper_30_4start              r=  5.0 p=2.0 taper= 30.0° starts=  4 → 770060 tris,    0 boundary edges [PASS]
taper_30_coarse              r=  5.0 p=5.0 taper= 30.0° starts=  1 → 612952 tris,  142 boundary edges [FAIL]
taper_15_steep               r=  2.5 p=3.0 taper= 15.0° starts=  1 → 578316 tris,    0 boundary edges [PASS]
taper_30_left                r=  5.0 p=2.0 taper= 30.0° starts= -1 → 749916 tris,    0 boundary edges [PASS]
taper_15_internal            r=  5.0 p=2.0 taper= 15.0° starts=  1 → 747948 tris,    0 boundary edges [PASS]
taper_15_acme                r=  5.0 p=2.0 taper= 15.0° starts=  1 → 753232 tris,    0 boundary edges [PASS]
taper_30_buttress            r=  5.0 p=2.0 taper= 30.0° starts=  1 → 871262 tris,   68 boundary edges [FAIL]
taper_15_plastic_butt        r=  5.0 p=2.0 taper= 15.0° starts=  1 → 767624 tris,    0 boundary edges [PASS]
taper_30_left_buttress       r=  5.0 p=2.0 taper= 30.0° starts= -1 → 871262 tris,   68 boundary edges [FAIL]
taper_15_multi4_buttress     r=  5.0 p=2.0 taper= 15.0° starts=  4 → 857288 tris,    0 boundary edges [PASS]
large_M64x6                  r= 32.0 p=6.0 taper=  0.0° starts=  1 → 1436496 tris,    0 boundary edges [PASS]

[snowbldr]

Should be fixed in #87 .

Interestingly these all did pass on my mac, so hopefully that fixes it.