TheRON/mechanical-compiler
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

stl: export a member as a sealed mesh

Browse Source 
prism_mesh takes a region and a length, not a Result. A member, its support, or
the support alone are then three calls with different regions rather than three
special cases inside one function -- support is a designed part of an artifact
where no printable orientation exists, not something a slicer adds, and
flattening everything into one mesh now would be expensive to undo later.

Named and documented as a MEMBER exporter. ROADMAP section 2 has three artifact
classes: members are prismatic and exist, nodes are non-prismatic and do not.
A prism sweep will never produce a node and no amount of extending this should
be attempted. A module called stl.py would invite exactly that.

No CAD kernel. Caps are a constrained Delaunay triangulation of the section
with its holes, walls are a quad strip per ring. That is why STL costs no new
dependency while STEP would.

Winding is forced, never assumed. region_parts returns each part outer-clockwise
with holes counter-clockwise because region_area depends on that to make holes
subtract; a prism swept along +Z needs the opposite, since traversing the outer
boundary counter-clockwise puts material on the left and the wall normal on the
right, which is outward. Both are reversed here.

The identifiers in the 80-byte header are PROVENANCE, NOT A CHECKSUM. Principle
4: mesh bytes are not reproducible across toolchain versions, so two exports
carrying the same build_id may differ byte for byte and both be correct -- a
GEOS release can move a triangulation without moving the geometry, which is the
F-034 mechanism applied to tessellation. The docstring says so where someone
diffing two exports will find it. Facet normals are written as zeros: every
consumer recomputes from vertex order, and a stored normal would be a second
source of truth for which way a face points.

A section failing is_region_simple is refused rather than exported. Its own
docstring names the case: an outline that measures perfectly, touches itself at
a point, and cannot be sealed. Emitting it would produce a file that opens in a
slicer and prints wrong.

member_stl defaults to the report's LENGTH_MM, the value the record's
VOLUME_MM3 and MASS_G were computed against. Any other sweep makes the record
describe a different object than the file beside it.

61 tests across all eleven profiles in both families, four properties each,
every negative assertion with a positive control. Edge pairing catches a missing
wall and inconsistent winding; Euler catches a filled or invented hole; area
overlap catches a bore paved over; volume against region_area times length
catches inverted normals and a wrong length together.

Two test defects found by mutation and fixed before landing, both worth
recording because both passed while being wrong.

The hole check originally tested triangle centroids. Ignoring holes at
triangulation time slipped through it -- two triangles can pave a bore with
both centroids outside it. Replaced with intersection area, which catches it.

The sweep-length test built at the family defaults, where preview_length_mm is
100, so a mutation hardcoding 100.0 passed the entire suite: the assertion was
comparing a value against the constant that had replaced it. Now parametrised
over 37.5 mm and 3048 mm, and proven against hardcoded mutations at both
values, each caught by the other case.

Mutation results against the real sections: holes not reversed fails 23,
exterior-only triangulation fails 45, outer ring not forced CCW fails 25,
simplicity precondition removed fails exactly 1 -- the refusal test written for
it. Restoration verified by checksum rather than by a green run, after a
same-length mutation plus a cached .pyc produced a false clean reading earlier
in the session. PYTHONDONTWRITEBYTECODE=1 throughout, per HANDOFF section 8.

Suite 627 passed. None of the 566 moved.
This commit is contained in:
TheRON
2026-09-12 10:40:35 -05:00
parent ecac644204
commit 0545b79674
2 changed files with 625 additions and 0 deletions
Download Patch File Download Diff File Expand all files Collapse all files
src/mechcomp/stl.py
+291
View File
@@ -0,0 +1,291 @@
"""
Export a member as a sealed triangle mesh.
WHAT THIS IS, AND WHAT IT IS NOT
**A member exporter.** ``ROADMAP.md`` §2 names three artifact classes:
members are prismatic and exist, nodes are non-prismatic and do not, panels
are sheet and do not. What is here sweeps a cross-section along +Z and
caps both ends. It will never produce a node, and no amount of extending it
should be attempted -- a node is a different representation, not a harder
prism.
The module is deliberately named for the artifact rather than the file
format for that reason. ``stl.py`` would invite someone to add nodes to it.
IDENTIFIERS IN THE OUTPUT ARE PROVENANCE, NOT A CHECKSUM
``ROADMAP.md`` §7 principle 4: mesh bytes are not reproducible across
toolchain versions. Two exports carrying the same ``build_id`` may differ
byte for byte and both be correct -- a Shapely or GEOS release can move a
triangulation without moving the geometry, which is the F-034 mechanism
applied to tessellation instead of to booleans.
So the identifiers written into the file say *what produced this*. They do
not say *these bytes*. Anyone diffing two exports of the same design and
concluding the compiler is broken has misread them, and the header says so
in the file itself.
WHY NO CAD KERNEL
A member is a prism. Caps are a constrained Delaunay triangulation of the
section, walls are a quad strip per ring, and both are decided by the
section's own rings. Nothing here needs a solid modeller, which is why STL
export costs no new dependency while STEP would.
THE PRECONDITION IS THE SECTION'S, NOT OURS
``region.is_region_simple`` exists because a section that touches itself at
a point is a valid 2D outline that cannot be tessellated: it measures
perfectly and then fails to extrude into a sealed solid. That is checked
here and refused, rather than emitted as a mesh that opens in a slicer and
prints wrong.
"""
from __future__ import annotations
import struct
from typing import Dict, List, NamedTuple, Sequence, Tuple
Point = Tuple[float, float]
Path = Sequence[Point]
Region = Sequence[Path]
Vertex = Tuple[float, float, float]
Triangle = Tuple[int, int, int]
class ExportRefused(Exception):
"""
The section cannot become a sealed solid.
A refusal rather than an error, in the same sense as ``ProfileRejected``:
the compiler declining to produce something it cannot vouch for is the
compiler working.
"""
class Mesh(NamedTuple):
vertices: List[Vertex]
triangles: List[Triangle]
@property
def volume(self) -> float:
"""
Signed volume by the divergence theorem.
Positive for outward-facing normals. Comparing this against
``region_area(section) * length`` checks the winding, the cap
orientation and the sweep length in a single number -- three defects
that are individually hard to see and jointly impossible to miss.
"""
total = 0.0
v = self.vertices
for i, j, k in self.triangles:
ax, ay, az = v[i]
bx, by, bz = v[j]
cx, cy, cz = v[k]
total += (ax * (by * cz - bz * cy)
- ay * (bx * cz - bz * cx)
+ az * (bx * cy - by * cx))
return total / 6.0
# ---------------------------------------------------------------------------
# Building the mesh
# ---------------------------------------------------------------------------
def _signed_area(ring: Path) -> float:
total = 0.0
n = len(ring)
for i in range(n):
x0, y0 = ring[i]
x1, y1 = ring[(i + 1) % n]
total += x0 * y1 - x1 * y0
return total / 2.0
def _ccw(ring: Path) -> List[Point]:
return list(ring) if _signed_area(ring) >= 0 else list(reversed(ring))
def _cw(ring: Path) -> List[Point]:
return list(ring) if _signed_area(ring) < 0 else list(reversed(ring))
class _Vertices:
"""Deduplicating vertex table, keyed on exact coordinates."""
def __init__(self) -> None:
self.items: List[Vertex] = []
self._index: Dict[Vertex, int] = {}
def add(self, x: float, y: float, z: float) -> int:
key = (x, y, z)
found = self._index.get(key)
if found is None:
found = len(self.items)
self._index[key] = found
self.items.append(key)
return found
def prism_mesh(section: Region, length: float) -> Mesh:
"""
Sweep ``section`` from z=0 to z=``length`` and seal both ends.
Takes a region and a length rather than a ``Result``, which is the whole
reason the signature looks like this. A member, its support, or the support
alone are then three calls with different regions instead of three special
cases inside one function -- and support is a designed part of an artifact
when no printable orientation exists, not something a slicer adds.
WINDING
``region_parts`` returns each part as an outer ring clockwise with its
holes counter-clockwise, because ``region_area`` depends on that
convention to make holes subtract. A prism swept along +Z needs the
opposite: traversing the outer boundary counter-clockwise puts the
material on the left and the wall normal on the right, which is
outward. Both are therefore reversed here. The convention is not
assumed from the input -- every ring is forced.
"""
from .geom.region import is_region_simple, region_parts, to_shapely
if length <= 0:
raise ExportRefused("length must be greater than zero, got %r" % (length,))
if not section:
raise ExportRefused("the section is empty; there is nothing to sweep")
if not is_region_simple(section):
raise ExportRefused(
"the section touches or crosses itself, so it cannot be sealed "
"into a solid. It is a valid outline and measures correctly; it "
"is not extrudable. See geom.region.is_region_simple.")
try:
from shapely import constrained_delaunay_triangles
except ImportError as exc: # pragma: no cover
raise ExportRefused(
"constrained_delaunay_triangles requires Shapely 2.1 or newer"
) from exc
verts = _Vertices()
tris: List[Triangle] = []
# -- caps ------------------------------------------------------------
#
# The triangulator is given the polygon with its holes, so no triangle
# lands inside a bore. That is asserted in the tests rather than trusted:
# a triangulator that filled the bores would produce a mesh that looks
# right in a slicer and prints solid through the conduit.
tessellation = constrained_delaunay_triangles(to_shapely(section))
faces = getattr(tessellation, "geoms", [tessellation])
cap_count = 0
for face in faces:
if face.is_empty:
continue
ring = [(x, y) for x, y in list(face.exterior.coords)[:-1]]
if len(ring) != 3:
continue
a, b, c = _ccw(ring)
cap_count += 1
# Top at z=length: counter-clockwise seen from +Z faces +Z, outward.
tris.append((verts.add(a[0], a[1], length),
verts.add(b[0], b[1], length),
verts.add(c[0], c[1], length)))
# Bottom at z=0: the same triangle reversed faces -Z, also outward.
tris.append((verts.add(c[0], c[1], 0.0),
verts.add(b[0], b[1], 0.0),
verts.add(a[0], a[1], 0.0)))
if cap_count == 0:
raise ExportRefused(
"the section produced no triangles; it encloses no area")
# -- walls -----------------------------------------------------------
for part in region_parts(section):
rings = [_ccw(part[0])] + [_cw(hole) for hole in part[1:]]
for ring in rings:
n = len(ring)
for i in range(n):
x0, y0 = ring[i]
x1, y1 = ring[(i + 1) % n]
b0 = verts.add(x0, y0, 0.0)
b1 = verts.add(x1, y1, 0.0)
t0 = verts.add(x0, y0, length)
t1 = verts.add(x1, y1, length)
tris.append((b0, b1, t1))
tris.append((b0, t1, t0))
return Mesh(verts.items, tris)
# ---------------------------------------------------------------------------
# Writing it out
# ---------------------------------------------------------------------------
HEADER_BYTES = 80
def stl_header(input_id: str = "", build_id: str = "") -> bytes:
"""
The binary format's 80-byte header, carrying provenance.
A file separated from its design record should still say what it is. Eighty
bytes is not much, so it holds the two identifiers and nothing else --
about fifty-four characters of the eighty.
These are NOT a checksum of the bytes that follow. See the module
docstring: mesh bytes are not reproducible across toolchain versions, and
two files with the same ``build_id`` may legitimately differ.
"""
text = "mechcomp"
if input_id:
text += " input=%s" % input_id
if build_id:
text += " build=%s" % build_id
raw = text.encode("ascii", "replace")[:HEADER_BYTES]
return raw + b"\0" * (HEADER_BYTES - len(raw))
def binary_stl(mesh: Mesh, input_id: str = "", build_id: str = "") -> bytes:
"""
The mesh as a binary STL.
Facet normals are written as zeros. The format has a field for them and
every consumer recomputes from the vertex order anyway; writing a stored
normal that could disagree with the winding would create a second source
of truth for which way a face points.
"""
out = [stl_header(input_id, build_id),
struct.pack("<I", len(mesh.triangles))]
v = mesh.vertices
for i, j, k in mesh.triangles:
out.append(struct.pack("<3f", 0.0, 0.0, 0.0))
for idx in (i, j, k):
out.append(struct.pack("<3f", *v[idx]))
out.append(struct.pack("<H", 0))
return b"".join(out)
def member_stl(result, length: float = None) -> bytes:
"""
Export whatever ``build()`` returned.
``length`` defaults to the report's ``LENGTH_MM``, which is the value the
record's ``VOLUME_MM3`` and ``MASS_G`` were computed against. Sweeping any
other length would make the record describe a different object than the
file shipped beside it, so the default is the only sensible one and an
override exists for callers that know what they are doing.
"""
if length is None:
length = result.report["LENGTH_MM"]
mesh = prism_mesh(result.section, length)
record = getattr(result, "record", None)
return binary_stl(mesh,
record.input_id if record else "",
record.build_id if record else "")
def stl_filename(result) -> str:
"""A name that identifies the design without needing the record beside it."""
record = getattr(result, "record", None)
ident = record.input_id if record else "unidentified"
family = result.report.get("FAMILY", "member")
profile = str(result.report.get("PROFILE", "")).replace(" ", "-").lower()
return "mechcomp-%s-%s-%s.stl" % (family, profile, ident)
tests/test_stl.py
+334
View File
@@ -0,0 +1,334 @@
"""
The member mesh: sealed, correctly oriented, and the right size.
WHY THESE FOUR CHECKS AND NOT A VISUAL ONE
A bad mesh opens in a slicer. That is the whole problem: inverted normals,
a hole filled in, a wall missing or the wrong sweep length all produce a
file that looks plausible and prints wrong, and the failure is discovered
in plastic hours later.
So each property is asserted arithmetically, and each is chosen to fail
loudly for a different defect:
edge pairing a missing or duplicated wall, and inconsistent winding
Euler a hole that was filled, or one invented
centroid-in-hole a bore tessellated over -- prints solid through the conduit
volume inverted normals and a wrong length, in one number
Every negative assertion has a positive control, per ROADMAP principle 11.
"""
from __future__ import annotations
import struct
import pytest
stl = pytest.importorskip("mechcomp.stl")
profiles = pytest.importorskip("mechcomp.profiles")
region = pytest.importorskip("mechcomp.geom.region")
LENGTH = 40.0
def every_profile():
"""Each catalogue entry of each family, at its family defaults."""
for name, family in profiles.FAMILIES.items():
for profile in sorted(family.catalogue):
yield name, profile
def built(family="3x", profile="Y", **params):
return profiles.build(family, profile, params)
ALL = list(every_profile())
def ids(pair):
return "%s/%s" % pair
# ---------------------------------------------------------------------------
# Sealed and consistently wound
# ---------------------------------------------------------------------------
@pytest.mark.parametrize("pair", ALL, ids=ids)
def test_every_directed_edge_appears_exactly_once(pair):
"""
Manifoldness and consistent orientation in one assertion.
In a closed, consistently oriented surface each directed edge (a, b) occurs
once and its reverse (b, a) occurs once. A duplicated directed edge means
two faces wind the same way across a shared edge; a missing reverse means a
hole in the surface.
"""
mesh = stl.prism_mesh(built(*pair).section, LENGTH)
seen = set()
for i, j, k in mesh.triangles:
for edge in ((i, j), (j, k), (k, i)):
assert edge not in seen, "directed edge %r twice" % (edge,)
seen.add(edge)
for a, b in seen:
assert (b, a) in seen, "edge %r has no opposite" % ((a, b),)
def test_the_edge_check_can_fail():
"""
The positive control. Without it, a mesh with no triangles at all would
satisfy the test above vacuously.
"""
mesh = stl.prism_mesh(built().section, LENGTH)
assert len(mesh.triangles) > 0
broken = stl.Mesh(mesh.vertices, mesh.triangles[:-1])
seen = set()
for i, j, k in broken.triangles:
for edge in ((i, j), (j, k), (k, i)):
seen.add(edge)
assert any((b, a) not in seen for a, b in seen)
# ---------------------------------------------------------------------------
# The right number of holes, still holes
# ---------------------------------------------------------------------------
@pytest.mark.parametrize("pair", ALL, ids=ids)
def test_euler_characteristic_matches_the_section(pair):
"""
V - E + F = 2P - 2H.
A prism over a part with H holes is a genus-H handlebody, so its boundary
surface has characteristic 2 - 2H, and P disjoint parts sum. Filling a bore
would raise the left side; inventing a tunnel would lower it.
"""
result = built(*pair)
section = result.section
parts = region.region_parts(section)
p = len(parts)
h = sum(len(part) - 1 for part in parts)
mesh = stl.prism_mesh(section, LENGTH)
edges = {frozenset((a, b))
for i, j, k in mesh.triangles
for a, b in ((i, j), (j, k), (k, i))}
chi = len(mesh.vertices) - len(edges) + len(mesh.triangles)
assert chi == 2 * p - 2 * h, (
"chi=%d for %d part(s) and %d hole(s)" % (chi, p, h))
def test_at_least_one_profile_actually_has_a_hole():
"""
The positive control for the Euler test and the one below it: if no
section in the catalogue had holes, both would be asserting the easy case
and a triangulator that ignored holes would pass everything.
"""
with_holes = 0
for pair in ALL:
parts = region.region_parts(built(*pair).section)
with_holes += sum(len(part) - 1 for part in parts)
assert with_holes > 0
# ---------------------------------------------------------------------------
# Nothing tessellated across a bore
# ---------------------------------------------------------------------------
@pytest.mark.parametrize("pair", ALL, ids=ids)
def test_no_cap_triangle_sits_inside_a_hole(pair):
"""
The check that matters most in practice.
A triangulator that filled the bores produces a mesh that is sealed,
correctly wound, passes Euler if the holes vanish consistently -- and
prints solid where the strap has to pass. Only the volume test and this one
would notice.
"""
result = built(*pair)
section = result.section
mesh = stl.prism_mesh(section, LENGTH)
holes = [hole for part in region.region_parts(section) for hole in part[1:]]
if not holes:
pytest.skip("this section has no holes")
caps = [t for t in mesh.triangles
if len({mesh.vertices[i][2] for i in t}) == 1]
assert caps, "no cap triangles found"
from shapely.geometry import Polygon
# Overlapping AREA, not the centroid. A centroid test passes by luck on a
# coarse triangulation: two triangles covering a square with a central bore
# can both have their centroids outside it while jointly paving it over.
# Found by mutation -- ignoring the holes at triangulation time slipped
# through a centroid check and is caught by this one.
hole_polys = [Polygon(h) for h in holes]
for i, j, k in caps:
tri = Polygon([(mesh.vertices[v][0], mesh.vertices[v][1])
for v in (i, j, k)])
if tri.area <= 0:
continue
for hp in hole_polys:
assert tri.intersection(hp).area <= 1e-9 * tri.area, (
"a cap triangle paves over a hole -- the part would print "
"solid where the stock has to pass")
# ---------------------------------------------------------------------------
# The right size, the right way out
# ---------------------------------------------------------------------------
@pytest.mark.parametrize("pair", ALL, ids=ids)
def test_volume_equals_area_times_length(pair):
"""
Catches inverted normals and a wrong sweep length together.
Compared against ``region_area`` at full precision rather than against the
report's ``SECTION_AREA_MM2``, which is rounded to six significant figures
and would cap this test's resolution at the oracle's rather than the
geometry's.
"""
result = built(*pair)
expected = region.area(result.section) * LENGTH
actual = stl.prism_mesh(result.section, LENGTH).volume
assert actual > 0, "negative volume means the normals face inward"
assert abs(actual - expected) <= 1e-9 * abs(expected)
def test_volume_agrees_with_the_report_too():
"""
The same quantity through the published numbers, at the precision they are
published to. Six significant figures, so the bound is loose on purpose.
"""
result = built()
length = result.report["LENGTH_MM"]
expected = result.report["SECTION_AREA_MM2"] * length
actual = stl.prism_mesh(result.section, length).volume
assert abs(actual - expected) <= 1e-5 * abs(expected)
def test_volume_scales_with_length():
"""Positive control: a volume that ignored its length would be constant."""
section = built().section
assert stl.prism_mesh(section, 80.0).volume == pytest.approx(
2.0 * stl.prism_mesh(section, 40.0).volume)
def test_reversing_a_wall_would_be_caught():
"""
Positive control for the sign: a mesh with every triangle flipped has
negative volume, so ``actual > 0`` above is a real assertion.
"""
mesh = stl.prism_mesh(built().section, LENGTH)
flipped = stl.Mesh(mesh.vertices, [(k, j, i) for i, j, k in mesh.triangles])
assert flipped.volume < 0
# ---------------------------------------------------------------------------
# The length that ships
# ---------------------------------------------------------------------------
@pytest.mark.parametrize("params,expected", [
({"preview_length_mm": 37.5}, 37.5),
({"length_view": "Full Length", "member_length_ft": 10}, 3048.0),
])
def test_member_stl_sweeps_the_reported_length(params, expected):
"""
The record's VOLUME_MM3 and MASS_G are computed against LENGTH_MM. Sweeping
anything else would make the record describe a different object than the
file beside it.
NEITHER CASE IS 100 mm, AND THAT IS THE POINT.
An earlier version of this test built at the family defaults, where
``preview_length_mm`` is 100. A mutation replacing the lookup with a
hardcoded 100.0 passed the whole suite -- the assertion compared a
value against the constant that had replaced it. Both cases here differ
from 100, and the second differs by two orders of magnitude, so no
plausible fixed length satisfies them.
"""
result = built(**params)
assert result.report["LENGTH_MM"] == expected
data = stl.member_stl(result)
count = struct.unpack("<I", data[80:84])[0]
zs = set()
for n in range(count):
base = 84 + n * 50 + 12
for v in range(3):
zs.add(round(struct.unpack("<3f", data[base + v * 12:base + v * 12 + 12])[2], 4))
assert max(zs) == pytest.approx(expected, rel=1e-6)
assert min(zs) == 0.0
def test_full_length_export_is_not_a_preview():
"""The defect that motivated the length_view control, asserted end to end."""
preview = built(length_view="Preview")
full = built(length_view="Full Length", member_length_ft=10)
assert full.report["LENGTH_MM"] == 3048
assert (stl.prism_mesh(full.section, full.report["LENGTH_MM"]).volume
> 30 * stl.prism_mesh(preview.section, preview.report["LENGTH_MM"]).volume)
# ---------------------------------------------------------------------------
# The file
# ---------------------------------------------------------------------------
def test_binary_stl_is_exactly_the_right_size():
mesh = stl.prism_mesh(built().section, LENGTH)
data = stl.binary_stl(mesh)
assert len(data) == 84 + 50 * len(mesh.triangles)
assert struct.unpack("<I", data[80:84])[0] == len(mesh.triangles)
def test_the_header_is_eighty_bytes_and_carries_both_ids():
result = built()
head = stl.stl_header(result.record.input_id, result.record.build_id)
assert len(head) == stl.HEADER_BYTES
text = head.rstrip(b"\0").decode("ascii")
assert result.record.input_id in text
assert result.record.build_id in text
def test_the_header_survives_having_no_record():
assert len(stl.stl_header()) == stl.HEADER_BYTES
def test_the_filename_identifies_the_design():
result = built()
name = stl.stl_filename(result)
assert result.record.input_id in name
assert name.endswith(".stl")
# ---------------------------------------------------------------------------
# Refusals
# ---------------------------------------------------------------------------
def test_an_empty_section_is_refused():
with pytest.raises(stl.ExportRefused):
stl.prism_mesh([], LENGTH)
@pytest.mark.parametrize("bad", [0.0, -1.0])
def test_a_non_positive_length_is_refused(bad):
with pytest.raises(stl.ExportRefused):
stl.prism_mesh(built().section, bad)
def test_a_section_that_touches_itself_is_refused():
"""
Two squares meeting at a single corner: a valid outline, correct area, and
not extrudable. ``is_region_simple`` exists for exactly this, and the
exporter must refuse rather than emit a mesh that opens in a slicer.
"""
figure_eight = [
[(0.0, 0.0), (10.0, 0.0), (10.0, 10.0), (0.0, 10.0)],
[(10.0, 10.0), (20.0, 10.0), (20.0, 20.0), (10.0, 20.0)],
]
assert not region.is_region_simple(figure_eight)
with pytest.raises(stl.ExportRefused):
stl.prism_mesh(figure_eight, LENGTH)
def test_a_well_formed_section_is_not_refused():
"""Positive control: the refusal is not simply rejecting everything."""
assert stl.prism_mesh(built().section, LENGTH).triangles