diff --git a/CONTEXT.md b/CONTEXT.md index 6d72724..6f73480 100644 --- a/CONTEXT.md +++ b/CONTEXT.md @@ -32,3 +32,30 @@ decimation/bottom-snap folds; the cross-module grid differences above are a suspected contributor. If unifying grids is ever attempted, it is a behaviour change — verify with the export→import round-trip, not the in-memory mesh. + +## Edge keys must be exact integers (`js/meshIndex.js`, `js/meshRepair.js`) + +Do **not** pack a vertex-id pair into one JS number as `a * 2**32 + b`. float64 +carries 53 bits of integer precision, so that form is exact only up to +`a = 2^21 = 2,097,152` — above it distinct edges collide onto one key, silently, +and only on meshes big enough that nobody verifies by hand. + +`meshRepair.js` used to do this in both `countEdgeDefects` and +`resolveTJunctions`, with different severities: + +* **countEdgeDefects** — colliding edges sum their incidence counts and trip the + `> 2` non-manifold test, so a *good* export is reported as broken. Measured on + a torus that is manifold by grid construction: 2.52 M vertices reported + 210,422 phantom non-manifold edges, 3.74 M reported 819,608. +* **resolveTJunctions** — worse, because it repairs rather than measures. A real + boundary edge (count 1) that collides reads as count 2 and its T-junction is + left unrepaired; and decoding the key back (`b = k % 4294967296`) returns + vertex ids that were never on that edge. + +Both now use `IntPairMap` (Int32 pair keys) over a dense edge table. Below the +2.1 M threshold the old keys were exact, so the change is a no-op there — which +is what the pipeline fingerprints confirm. + +`diag-edgekey-collision.mjs` reproduces the failure and is the regression test: +it builds meshes whose manifoldness is guaranteed by topology, not measured, so +any counter that disagrees is wrong by construction. diff --git a/diag-edgekey-collision.mjs b/diag-edgekey-collision.mjs new file mode 100644 index 0000000..61f8556 --- /dev/null +++ b/diag-edgekey-collision.mjs @@ -0,0 +1,91 @@ +/* + * Copyright (c) 2026 CNCKitchen (Stefan Hermann) and contributors + * SPDX-License-Identifier: AGPL-3.0-only + */ + +// countEdgeDefects / resolveTJunctions key an edge as `x * 4294967296 + y` +// (x * 2^32 + y) in a JS number. That is exact only while x * 2^32 + y stays +// within float64's 2^53 integer range, i.e. up to x = 2^21 = 2,097,152 +// vertices. Past that, DISTINCT edges collide onto one key, their incidence +// counts add up, and the sum trips the `> 2` test — so a perfectly manifold +// mesh is reported as having non-manifold edges. +// +// This builds a closed torus (manifold by construction: every edge has exactly +// two incident faces, guaranteed by the grid topology, not by measurement) at a +// range of sizes and compares the packed-float key against exact integer keys. +// +// node --max-old-space-size=8000 diag-edgekey-collision.mjs +import * as THREE from 'three'; +import { countEdgeDefects } from './js/meshRepair.js'; +import { IntPairMap } from './js/meshIndex.js'; + +// Exact reference counter: integer pair keys, no float packing, no Map cap. +function exactDefects(geometry, Q = 1e4) { + const p = geometry.attributes.position.array, n = p.length / 9; + const vmap = new IntPairMap(Math.ceil(n * 2)); + // Quantise to the same grid, then id via a 3-key weld built from two pairs. + const idOf = new Map(); + const id = new Int32Array(n * 3); + let next = 0; + for (let i = 0; i < n * 3; i++) { + const k = Math.round(p[i*3]*Q) + '/' + Math.round(p[i*3+1]*Q) + '/' + Math.round(p[i*3+2]*Q); + let v = idOf.get(k); + if (v === undefined) { v = next++; idOf.set(k, v); } + id[i] = v; + } + const edges = new IntPairMap(Math.ceil(n * 1.6)); + let nE = 0; + let counts = new Int32Array(Math.ceil(n * 1.8) + 16); + for (let t = 0; t < n; t++) { + const a = id[t*3], b = id[t*3+1], c = id[t*3+2]; + if (a === b || b === c || a === c) continue; + const tri = [a, b, c]; + for (let e = 0; e < 3; e++) { + const x = tri[e], y = tri[(e+1)%3]; + const lo = x < y ? x : y, hi = x < y ? y : x; + const s = edges.getOrSet(lo, hi, nE); + if (edges.inserted) { + if (nE >= counts.length) { const g = new Int32Array(counts.length*2); g.set(counts); counts = g; } + nE++; + } + counts[s]++; + } + } + let open = 0, nonManifold = 0; + for (let i = 0; i < nE; i++) { if (counts[i] === 1) open++; else if (counts[i] > 2) nonManifold++; } + return { open, nonManifold, verts: next }; +} + +// Closed torus: a full grid wrap in both directions, so EVERY edge borders +// exactly two faces. Manifold by topology, independent of any counter. +function torus(nu, nv, R = 100, r = 30) { + const pos = new Float32Array(nu * nv * 2 * 9); + let o = 0; + const P = (i, j) => { + const u = 2*Math.PI*(i % nu)/nu, v = 2*Math.PI*(j % nv)/nv; + return [(R + r*Math.cos(v))*Math.cos(u), (R + r*Math.cos(v))*Math.sin(u), r*Math.sin(v)]; + }; + for (let i = 0; i < nu; i++) for (let j = 0; j < nv; j++) { + const a = P(i,j), b = P(i+1,j), c = P(i+1,j+1), d = P(i,j+1); + for (const t of [[a,b,c],[a,c,d]]) for (const p of t) { pos[o++]=p[0]; pos[o++]=p[1]; pos[o++]=p[2]; } + } + const g = new THREE.BufferGeometry(); + g.setAttribute('position', new THREE.BufferAttribute(pos, 3)); + return g; +} + +console.log('Closed torus — every edge has exactly 2 incident faces by construction.\n'); +console.log(' vertices tris | upstream countEdgeDefects | exact integer keys'); +console.log(' ' + '-'.repeat(72)); +for (const [nu, nv] of [[400,300],[900,700],[1400,1100],[1800,1400],[2200,1700]]) { + const g = torus(nu, nv); + const tris = g.attributes.position.count / 3; + let up; + try { up = countEdgeDefects(g); } catch (e) { up = { err: e.message }; } + const ex = exactDefects(g); + const upStr = up.err ? `THREW: ${up.err}` : `open=${String(up.open).padStart(6)} nonManifold=${String(up.nonManifold).padStart(9)}`; + console.log(` ${String(ex.verts).padStart(8)} ${String(tris).padStart(8)} | ${upStr.padEnd(37)} | open=${ex.open} nonManifold=${ex.nonManifold}`); + g.dispose(); +} +console.log('\n float64 keeps x*2^32+y exact only to x = 2^21 = 2,097,152 vertices.'); +console.log(' Above that the upstream counter reports defects a manifold mesh does not have.'); diff --git a/js/meshIndex.js b/js/meshIndex.js index 9ec840d..7f87ce0 100644 --- a/js/meshIndex.js +++ b/js/meshIndex.js @@ -108,6 +108,98 @@ export class QuantizedPointMap { } } +/** + * IntPairMap — open-addressing map keyed by an ORDERED PAIR OF INTEGERS. + * + * Several places key on two vertex ids rather than a 3-D point — edge tables, + * above all. The tempting shortcut is to pack the pair into one JS number as + * `a * 2**32 + b` and use it in a Map. **That is only correct up to + * a = 2^21 = 2,097,152.** float64 carries 53 bits of integer precision, so + * beyond that the packed value rounds and distinct pairs collide — silently, + * and only on meshes large enough that nobody checks by hand. meshRepair.js + * did exactly this and reported hundreds of thousands of non-manifold edges on + * meshes that were provably manifold (see diag-edgekey-collision.mjs). + * + * This class keeps the two components as separate Int32s, so the key is exact + * for every id the pipeline can produce (< 2^31), and it has no equivalent of + * a JS Map's ~16.7 M entry cap. + * + * It is also 12 bytes per slot against QuantizedPointMap's 28 (3 × Float64 key + * + Int32 value), which matters wherever an integer-pair table is sized by + * triangle count. + * + * Values must be integers in [0, 2^31-1]; -1 is the "empty" sentinel and is + * what get() returns on a miss — same contract as QuantizedPointMap. + */ +export class IntPairMap { + /** + * @param {number} expected – expected number of unique pairs (sizing hint) + */ + constructor(expected = 256) { + /** true when the last getOrSet() inserted a new key */ + this.inserted = false; + this._size = 0; + let cap = 16; + const target = Math.max(16, Math.ceil(expected / 0.6)); + while (cap < target) cap *= 2; + this._alloc(cap); + } + + get size() { return this._size; } + + _alloc(cap) { + this._cap = cap; + this._mask = cap - 1; + this._k1 = new Int32Array(cap); + this._k2 = new Int32Array(cap); + this._val = new Int32Array(cap).fill(-1); + } + + _slot(a, b) { + let h = Math.imul(a, 0x9E3779B1) ^ Math.imul(b, 0x85EBCA77); + h ^= h >>> 15; + let i = h & this._mask; + const k1 = this._k1, k2 = this._k2, val = this._val, mask = this._mask; + while (val[i] !== -1) { + if (k1[i] === a && k2[i] === b) return i; + i = (i + 1) & mask; + } + return i; + } + + _grow() { + const ok1 = this._k1, ok2 = this._k2, oval = this._val, ocap = this._cap; + this._alloc(ocap * 2); + for (let i = 0; i < ocap; i++) { + if (oval[i] === -1) continue; + const s = this._slot(ok1[i], ok2[i]); + this._k1[s] = ok1[i]; this._k2[s] = ok2[i]; this._val[s] = oval[i]; + } + } + + /** Value stored for the pair (a,b), or -1 if absent. */ + get(a, b) { + return this._val[this._slot(a, b)]; + } + + /** + * Return the value already stored for (a,b); if absent, store `value` and + * return it. `this.inserted` tells which case occurred. + */ + getOrSet(a, b, value) { + const i = this._slot(a, b); + const existing = this._val[i]; + if (existing !== -1) { + this.inserted = false; + return existing; + } + this._k1[i] = a; this._k2[i] = b; this._val[i] = value; + this.inserted = true; + if (++this._size > this._cap * 0.7) this._grow(); + return value; + } +} + /** * Weld a non-indexed position buffer: assign each vertex the sequential id of * its quantised position (first occurrence wins). diff --git a/js/meshRepair.js b/js/meshRepair.js index 3dcf3d7..81fdb67 100644 --- a/js/meshRepair.js +++ b/js/meshRepair.js @@ -32,7 +32,7 @@ * @returns {THREE.BufferGeometry} repaired non-indexed geometry */ import { THREE } from './threeCompat.js'; -import { QuantizedPointMap } from './meshIndex.js'; +import { QuantizedPointMap, IntPairMap } from './meshIndex.js'; /** * Count open (1-face) and non-manifold (3+-face) edges of a non-indexed geometry, @@ -46,19 +46,40 @@ export function countEdgeDefects(geometry, Q = 1e4) { for (let i = 0; i < n * 3; i++) { id[i] = vmap.getOrSet(p[i*3], p[i*3+1], p[i*3+2], vmap.size); } - const ec = new Map(); + // Exact integer pair keys over typed arrays. The previous implementation + // packed each edge as `x * 4294967296 + y` into a JS number and counted in a + // Map, which failed twice on large exports: + // + // * float64 holds x * 2^32 + y exactly only up to x = 2^21 = 2,097,152 + // vertices. Beyond that distinct edges share a key and their incidence + // counts add, so the `> 2` test fires on a perfectly manifold mesh — + // measured on a torus that is manifold by construction, 2.52 M vertices + // reported 210,422 phantom non-manifold edges and 3.74 M reported + // 819,608. Exports were being declared broken when they were fine. + // * a JS Map throws past V8's ~16.7 M entry cap, i.e. above ~11 M + // triangles, which the pipeline can now reach. + const edges = new IntPairMap(Math.max(16, n * 2)); + let nE = 0; + let counts = new Int32Array(n * 3 + 16); for (let t = 0; t < n; t++) { const a = id[t*3], b = id[t*3+1], c = id[t*3+2]; if (a === b || b === c || a === c) continue; - const tri = [a, b, c]; for (let e = 0; e < 3; e++) { - const x = tri[e], y = tri[(e+1)%3]; - const key = x < y ? x * 4294967296 + y : y * 4294967296 + x; - ec.set(key, (ec.get(key) || 0) + 1); + const x = e === 0 ? a : e === 1 ? b : c; + const y = e === 0 ? b : e === 1 ? c : a; + const lo = x < y ? x : y, hi = x < y ? y : x; + const s = edges.getOrSet(lo, hi, nE); + if (edges.inserted) { + if (nE >= counts.length) { + const g = new Int32Array(counts.length * 2); g.set(counts); counts = g; + } + nE++; + } + counts[s]++; } } let open = 0, nonManifold = 0; - for (const c of ec.values()) { if (c === 1) open++; else if (c > 2) nonManifold++; } + for (let i = 0; i < nE; i++) { if (counts[i] === 1) open++; else if (counts[i] > 2) nonManifold++; } return { open, nonManifold, tris: n }; } @@ -133,21 +154,50 @@ export function resolveTJunctions(geometry, opts = {}) { faces.push([a, b, c]); } - const ekey = (a, b) => (a < b ? a * 4294967296 + b : b * 4294967296 + a); - for (let iter = 0; iter < maxIters; iter++) { - // Edge → adjacent-face count + sample. - const eCount = new Map(); + // Edge → adjacent-face count, in a dense table keyed on the EXACT integer + // pair. The previous key packed the pair into one JS number as + // `a * 4294967296 + b`, which float64 represents exactly only while the + // product stays inside 2^53 — i.e. up to a = 2^21 = 2,097,152 vertices. + // Above that distinct edges land on the same key and their counts merge, + // which here is a CORRECTNESS bug and not merely a reporting one: a real + // boundary edge (count 1) colliding with another reads as count 2, so its + // T-junction is never repaired, and decoding the key back into (a, b) — + // `b = k % 4294967296` — yields vertex ids that were never on that edge. + // Measured on a torus that is manifold by construction: 2.52 M vertices + // reported 210,422 phantom non-manifold edges, 3.74 M reported 819,608. + const eMap = new IntPairMap(Math.max(16, faces.length * 2)); + let nE = 0; + let eLo = new Int32Array(faces.length * 3 + 16); + let eHi = new Int32Array(eLo.length); + let eCnt = new Int32Array(eLo.length); + const edgeSlot = (a, b) => { + const lo = a < b ? a : b, hi = a < b ? b : a; + const s = eMap.getOrSet(lo, hi, nE); + if (eMap.inserted) { + if (nE >= eLo.length) { + const grow = (o) => { const g = new Int32Array(o.length * 2); g.set(o); return g; }; + eLo = grow(eLo); eHi = grow(eHi); eCnt = grow(eCnt); + } + eLo[nE] = lo; eHi[nE] = hi; nE++; + } + return s; + }; + const edgeCount = (a, b) => { + const lo = a < b ? a : b, hi = a < b ? b : a; + const s = eMap.get(lo, hi); + return s === -1 ? 0 : eCnt[s]; + }; + for (let fi = 0; fi < faces.length; fi++) { const f = faces[fi]; - for (let e = 0; e < 3; e++) eCount.set(ekey(f[e], f[(e+1)%3]), (eCount.get(ekey(f[e], f[(e+1)%3])) || 0) + 1); + for (let e = 0; e < 3; e++) eCnt[edgeSlot(f[e], f[(e+1)%3])]++; } // Boundary edges (exactly one face) and the set of boundary vertices. const bverts = new Set(); - for (const [k, c] of eCount) { - if (c !== 1) continue; - const b = k % 4294967296, a = (k - b) / 4294967296; - bverts.add(a); bverts.add(b); + for (let s = 0; s < nE; s++) { + if (eCnt[s] !== 1) continue; + bverts.add(eLo[s]); bverts.add(eHi[s]); } if (bverts.size === 0) break; const bvArr = [...bverts]; @@ -159,7 +209,7 @@ export function resolveTJunctions(geometry, opts = {}) { const f = faces[fi]; for (let e = 0; e < 3; e++) { const a = f[e], b = f[(e+1)%3]; - if ((eCount.get(ekey(a, b)) || 0) !== 1) continue; // only boundary edges + if (edgeCount(a, b) !== 1) continue; // only boundary edges const ax = vx[a], ay = vy[a], az = vz[a]; const ex = vx[b]-ax, ey = vy[b]-ay, ez = vz[b]-az; const elen2 = ex*ex + ey*ey + ez*ez;