Add programme/solver/oracle + sizing experiments (negative result)

Adds the bottom-up ratio solver, programme parser, Perl-oracle bridge,
and two experiments. Headline finding: the "isolated size solver on a
frozen topology" hypothesis is NOT validated.

- resolve_ratios.py: re-solving candidate-002 from programme targets
  recovers areas accurately but scores below the original (introduces
  width/perpendicular/crinkliness failures the area objective ignores).
- refine_sweep.py: warm-start refine of all 34 evolved candidates
  regresses 34/34 (fails 124->297 perpendicular-tied; 124->626 area-only
  with free skew). Moving cuts to fix room area breaks the coupled
  adjacency/access/shape constraints those designs balanced.

Conclusion: sizing is not separable from the rest of Urb's fitness;
a geometry inner loop must optimise the full objective, not an area proxy.
Geometry port remains validated byte-identical to Urb.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>

experiments/refine_sweep.py

@@ -0,0 +1,79 @@
+"""Population sweep: warm-start refine every evolved candidate and tally results.
+
+For each real .dom in the example dir, score it, run the solver as a geometry
+optimiser (warm-start, no strip), and re-score. Reports how often bottom-up
+sizing improves vs regresses total fitness, plus aggregate fail-count change.
+
+This is a breadth check on the solver-as-optimiser role; raw fitness is still
+confounded by the 0.5^n failure cliff and any topological defects, so the
+fail-count and per-candidate detail matter as much as the win/loss tally.
+"""
+
+import shutil
+import sys
+from pathlib import Path
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1] / "src"))
+
+from homemaker import dom, oracle, programme, solver  # noqa: E402
+
+URB = Path("/home/bruno/src/urb")
+EX = URB / "examples/programme-house"
+
+
+def _is_candidate(p: Path) -> bool:
+    # real designs: 32-hex hashes or candidate-NNN; skip init and our scratch
+    name = p.stem
+    return name not in {"init", "original", "roundtrip", "solved", "refined"}
+
+
+def main() -> None:
+    scratch = Path(__file__).resolve().parents[1] / "scratch"
+    scratch.mkdir(exist_ok=True)
+    shutil.copy(EX / "patterns.config", scratch / "patterns.config")
+    targets = programme.load_programme(str(EX / "patterns.config"))
+
+    doms = sorted(p for p in EX.glob("*.dom") if _is_candidate(p))
+    win = loss = tie = 0
+    fails_before = fails_after = 0
+    rows = []
+    for src in doms:
+        try:
+            shutil.copy(src, scratch / "orig.dom")
+            s0 = oracle.score(scratch / "orig.dom", URB)
+            root = dom.load(str(src))
+            # gentlest refiner: nudge cut POSITIONS for programme-room area only,
+            # keep evolved cut angles and leave circulation/shape untouched.
+            solver.solve_ratios(
+                root, targets, strip=False, perpendicular=False,
+                weight_width=0.0, weight_proportion=0.0, min_width_generic=0.0,
+            )
+            dom.dump(root, str(scratch / "ref.dom"))
+            s1 = oracle.score(scratch / "ref.dom", URB)
+        except Exception as e:  # noqa: BLE001
+            rows.append(f"  {src.name:40s} ERROR {e}")
+            continue
+        fails_before += s0.n_fails
+        fails_after += s1.n_fails
+        if s1.fitness > s0.fitness * 1.001:
+            win += 1
+            mark = "＋"
+        elif s1.fitness < s0.fitness * 0.999:
+            loss += 1
+            mark = "－"
+        else:
+            tie += 1
+            mark = "="
+        rows.append(
+            f"  {mark} {src.name:40s} {s0.fitness:.4g} -> {s1.fitness:.4g}"
+            f"  fails {s0.n_fails}->{s1.n_fails}"
+        )
+
+    print("\n".join(rows))
+    n = win + loss + tie
+    print(f"\n{n} candidates: {win} improved, {loss} regressed, {tie} tied")
+    print(f"total fails: {fails_before} -> {fails_after}")
+
+
+if __name__ == "__main__":
+    main()

experiments/resolve_ratios.py

@@ -0,0 +1,90 @@
+"""Go/no-go experiment for bottom-up sizing.
+
+Take a real evolved .dom, throw away its division ratios, and re-solve them from
+the programme's target sizes alone. Score three versions through the Perl oracle:
+
+  original  -- the evolved .dom as-is (baseline)
+  roundtrip -- loaded and re-emitted unmodified (checks dump fidelity)
+  solved    -- ratios stripped to 0.5 then solved from programme targets
+
+If `solved` scores >= `original`, sizing can be recovered from the programme
+without the evolved geometry, and the EA only needs to search topology.
+
+Usage:
+  python experiments/resolve_ratios.py [source.dom] [--urb /path/to/urb]
+"""
+
+import argparse
+import shutil
+import sys
+from pathlib import Path
+
+sys.path.insert(0, str(Path(__file__).resolve().parents[1] / "src"))
+
+from homemaker import dom, oracle, programme, solver  # noqa: E402
+
+DEFAULT_URB = Path("/home/bruno/src/urb")
+DEFAULT_SRC = DEFAULT_URB / "examples/programme-house/candidate-002.dom"
+
+
+def _score_in(scratch: Path, name: str, root: dom.Node, urb: Path) -> oracle.Score:
+    path = scratch / name
+    dom.dump(root, str(path))
+    return oracle.score(path, urb)
+
+
+def main() -> None:
+    ap = argparse.ArgumentParser()
+    ap.add_argument("source", nargs="?", default=str(DEFAULT_SRC))
+    ap.add_argument("--urb", default=str(DEFAULT_URB))
+    args = ap.parse_args()
+
+    src = Path(args.source).resolve()
+    urb = Path(args.urb).resolve()
+    config = src.parent / "patterns.config"
+
+    scratch = Path(__file__).resolve().parents[1] / "scratch"
+    scratch.mkdir(exist_ok=True)
+    shutil.copy(config, scratch / "patterns.config")
+
+    targets = programme.load_programme(str(config))
+
+    # baseline: score the original file as-is
+    orig_copy = scratch / "original.dom"
+    shutil.copy(src, orig_copy)
+    s_orig = oracle.score(orig_copy, urb)
+
+    # roundtrip: load + re-emit unmodified
+    s_round = _score_in(scratch, "roundtrip.dom", dom.load(str(src)), urb)
+
+    # solved: strip ratios and re-solve from programme targets
+    root = dom.load(str(src))
+    print("--- programme leaves BEFORE solve (ratios intact) ---")
+    print(solver.area_report(root, targets))
+    res = solver.solve_ratios(root, targets, strip=True)
+    print("\n--- programme leaves AFTER solve (from targets, ratios stripped) ---")
+    print(solver.area_report(root, targets))
+    print(f"\nsolver: cost={res.cost:.4f} nfev={res.nfev} success={res.success}")
+    s_solved = _score_in(scratch, "solved.dom", root, urb)
+
+    # refined: warm-start from the evolved ratios (solver as geometry optimiser)
+    root2 = dom.load(str(src))
+    solver.solve_ratios(root2, targets, strip=False)
+    s_refined = _score_in(scratch, "refined.dom", root2, urb)
+
+    print("\n=== FITNESS (via urb-fitness.pl oracle) ===")
+    for label, s in (
+        ("original", s_orig),
+        ("roundtrip", s_round),
+        ("solved", s_solved),
+        ("refined", s_refined),
+    ):
+        print(f"  {label:9s} fitness={s.fitness:.10g}  fails={s.n_fails}")
+
+    print("\nVERDICT:")
+    print(f"  solved  (strip, frozen topology): {'>=' if s_solved.fitness >= s_orig.fitness else '<'} original")
+    print(f"  refined (warm-start optimiser):   {'>=' if s_refined.fitness >= s_orig.fitness else '<'} original")
+
+
+if __name__ == "__main__":
+    main()

src/homemaker/dom.py

@@ -156,6 +156,7 @@ def load(path: str) -> Node:
     if root.node is not None:
         root.node_file = [list(p) for p in root.node]
         root.node = geometry.offset_quad(root.node, -root.wall_outer)
+    geometry.clear_cache()  # fresh tree: drop any stale coordinates
     return root
 
 

src/homemaker/geometry.py

@@ -19,6 +19,17 @@ from .dom import Node
 
 Point = list[float]
 
+# Memoisation of derived coordinates. The pull-based recursion mirrors Urb but,
+# uncached, re-derives ancestor/below corners exponentially with depth. Urb
+# itself caches in the node and clears via Clean_Cache(); we do the same with a
+# module cache keyed by node identity. Callers that mutate divisions (the
+# solver) must call clear_cache(); dom.load() clears it for a fresh tree.
+_cache: dict = {}
+
+
+def clear_cache() -> None:
+    _cache.clear()
+
 
 def _interp(a: Point, b: Point, t: float) -> Point:
     return [a[0] * (1 - t) + b[0] * t, a[1] * (1 - t) + b[1] * t]
@@ -26,30 +37,52 @@ def _interp(a: Point, b: Point, t: float) -> Point:
 
 def coordinate(n: Node, idx: int) -> Point:
     """Corner ``idx`` (0..3) of ``n``; mirrors ``Urb::Quad::Coordinate``."""
+    key = (id(n), idx)
+    hit = _cache.get(key)
+    if hit is not None:
+        return hit
     if n.below is not None:  # upper storey inherits geometry from below
-        return coordinate(n.below, idx)
-    rid = (idx + n.rotation) % 4
-    if n.parent is None:  # level root: stored, rotation-adjusted corner
-        return list(n.node[rid])
-    p = n.parent
-    if n.position == "l":
-        return {0: coordinate(p, 0), 1: coord_a(p), 2: coord_b(p), 3: coordinate(p, 3)}[rid]
-    # position == 'r'
-    return {0: coord_a(p), 1: coordinate(p, 1), 2: coordinate(p, 2), 3: coord_b(p)}[rid]
+        result = coordinate(n.below, idx)
+    else:
+        rid = (idx + n.rotation) % 4
+        if n.parent is None:  # level root: stored, rotation-adjusted corner
+            result = list(n.node[rid])
+        else:
+            p = n.parent
+            if n.position == "l":
+                result = {0: coordinate(p, 0), 1: coord_a(p), 2: coord_b(p), 3: coordinate(p, 3)}[rid]
+            else:  # 'r'
+                result = {0: coord_a(p), 1: coordinate(p, 1), 2: coordinate(p, 2), 3: coord_b(p)}[rid]
+    _cache[key] = result
+    return result
 
 
 def coord_a(n: Node) -> Point:
     """End 'a' of the division line; mirrors ``Urb::Quad::Coordinate_a``."""
+    key = (id(n), "a")
+    hit = _cache.get(key)
+    if hit is not None:
+        return hit
     if n.below is not None and n.below.divided:
-        return coord_a(n.below)
-    return _interp(coordinate(n, 0), coordinate(n, 1), n.division[0])
+        result = coord_a(n.below)
+    else:
+        result = _interp(coordinate(n, 0), coordinate(n, 1), n.division[0])
+    _cache[key] = result
+    return result
 
 
 def coord_b(n: Node) -> Point:
     """End 'b' of the division line; mirrors ``Urb::Quad::Coordinate_b``."""
+    key = (id(n), "b")
+    hit = _cache.get(key)
+    if hit is not None:
+        return hit
     if n.below is not None and n.below.divided:
-        return coord_b(n.below)
-    return _interp(coordinate(n, 3), coordinate(n, 2), n.division[1])
+        result = coord_b(n.below)
+    else:
+        result = _interp(coordinate(n, 3), coordinate(n, 2), n.division[1])
+    _cache[key] = result
+    return result
 
 
 def _dist(a: Point, b: Point) -> float:

src/homemaker/oracle.py

@@ -0,0 +1,52 @@
+"""Phase-1 fitness oracle: score a ``.dom`` via Urb's ``urb-fitness.pl``.
+
+This is the only throwaway component. It shells out to the Perl evaluator so we
+can validate the Python search core against the trusted fitness before porting
+fitness to Python (Phase 2). ``urb-fitness.pl`` reads ``patterns.config`` from
+its working directory, so the ``.dom`` must live beside the programme config.
+"""
+
+from __future__ import annotations
+
+import os
+import subprocess
+from dataclasses import dataclass
+from pathlib import Path
+
+DEFAULT_URB_ROOT = Path("/home/bruno/src/urb")
+
+
+@dataclass
+class Score:
+    fitness: float
+    fails: str  # raw .fails content (YAML and/or plain lines)
+
+    @property
+    def n_fails(self) -> int:
+        return sum(1 for line in self.fails.splitlines() if line.strip() and line.strip() != "---")
+
+
+def score(dom_path: str | Path, urb_root: str | Path = DEFAULT_URB_ROOT) -> Score:
+    dom_path = Path(dom_path).resolve()
+    urb_root = Path(urb_root).resolve()
+    score_file = Path(f"{dom_path}.score")
+    fails_file = Path(f"{dom_path}.fails")
+    for f in (score_file, fails_file):
+        f.unlink(missing_ok=True)
+
+    env = {**os.environ, "DEBUG": "1"}
+    proc = subprocess.run(
+        ["perl", f"-I{urb_root}/lib", str(urb_root / "bin" / "urb-fitness.pl"), dom_path.name],
+        cwd=dom_path.parent,
+        env=env,
+        capture_output=True,
+        text=True,
+    )
+    if not score_file.exists():
+        raise RuntimeError(
+            f"urb-fitness.pl produced no score for {dom_path}\n"
+            f"stdout:\n{proc.stdout}\nstderr:\n{proc.stderr}"
+        )
+    fitness = float(score_file.read_text().strip())
+    fails = fails_file.read_text() if fails_file.exists() else ""
+    return Score(fitness=fitness, fails=fails)

src/homemaker/programme.py

@@ -0,0 +1,64 @@
+"""Parse a ``patterns.config`` programme into per-code space requirements.
+
+Only the ``spaces:`` section is read here. Generic codes (c/o/s) carry no
+explicit targets and are left unconstrained by the solver (they absorb slack).
+"""
+
+from __future__ import annotations
+
+from dataclasses import dataclass, field
+
+import yaml
+
+# Urb::Dom::Fitness defaults for optional params (ProgrammeDriven.default_params).
+_DEFAULT_WIDTH = (4.0, 1.0)
+_DEFAULT_PROPORTION = (1.5, 0.5)
+
+
+@dataclass
+class SpaceReq:
+    code: str
+    name: str = ""
+    size: float = 0.0  # target floor area, m^2
+    size_sigma: float = 1.0
+    width: float = _DEFAULT_WIDTH[0]
+    width_sigma: float = _DEFAULT_WIDTH[1]
+    proportion: float = _DEFAULT_PROPORTION[0]  # max length/width ratio
+    proportion_sigma: float = _DEFAULT_PROPORTION[1]
+    adjacency: list[str] = field(default_factory=list)
+    level: int | None = None
+    requires_below: str | None = None
+    count: int = 1
+
+
+def _pair(d: dict, key: str, default: tuple[float, float]) -> tuple[float, float]:
+    v = d.get(key)
+    if v is None:
+        return default
+    return float(v[0]), float(v[1])
+
+
+def load_programme(path: str) -> dict[str, SpaceReq]:
+    with open(path) as fh:
+        conf = yaml.safe_load(fh)
+    spaces = conf.get("spaces") or {}
+    out: dict[str, SpaceReq] = {}
+    for code, c in spaces.items():
+        size = _pair(c, "size", (0.0, 1.0))
+        width = _pair(c, "width", _DEFAULT_WIDTH)
+        prop = _pair(c, "proportion", _DEFAULT_PROPORTION)
+        out[code] = SpaceReq(
+            code=code,
+            name=c.get("name", ""),
+            size=size[0],
+            size_sigma=size[1],
+            width=width[0],
+            width_sigma=width[1],
+            proportion=prop[0],
+            proportion_sigma=prop[1],
+            adjacency=list(c.get("adjacency") or []),
+            level=c.get("level"),
+            requires_below=c.get("requires_below"),
+            count=int(c.get("count") or 1),
+        )
+    return out

src/homemaker/solver.py

@@ -0,0 +1,150 @@
+"""Bottom-up division-ratio solver.
+
+Given a *fixed* slicing topology (types, rotations, tree shape), solve every
+free division ratio so that each programme leaf best meets its target area —
+with soft, one-sided penalties for being too narrow or too elongated. This is
+the inversion of Urb's top-down sizing: rooms declare targets, geometry follows.
+
+Only generic leaves (circulation/outside/storage) and unconstrained types are
+left to absorb the residual area, exactly as a real plan lets corridors flex.
+
+A division is *free* only at the lowest storey where its tree path is divided;
+higher storeys inherit that cut via Below-inheritance (see geometry.coordinate),
+so their stored ratios are dead variables and must not be optimised.
+"""
+
+from __future__ import annotations
+
+import numpy as np
+from scipy.optimize import least_squares
+
+from . import geometry
+from .dom import Node, levels
+from .programme import SpaceReq
+
+_EPS = 0.02  # keep cuts off the edges
+
+
+def _branches(n: Node) -> list[Node]:
+    if not n.divided:
+        return []
+    return [n] + _branches(n.left) + _branches(n.right)
+
+
+def free_branches(root: Node) -> list[Node]:
+    """Branches whose division actually drives geometry (not inherited)."""
+    out: list[Node] = []
+    for lvl in levels(root):
+        for b in _branches(lvl):
+            if b.below is None or not b.below.divided:
+                out.append(b)
+    return out
+
+
+def _width(leaf: Node) -> float:
+    l0 = geometry.edge_length(leaf, 0)
+    l1 = geometry.edge_length(leaf, 1)
+    l2 = geometry.edge_length(leaf, 2)
+    l3 = geometry.edge_length(leaf, 3)
+    return min((l0 + l2) / 2, (l1 + l3) / 2)
+
+
+def _aspect(leaf: Node) -> float:
+    l0 = geometry.edge_length(leaf, 0)
+    l1 = geometry.edge_length(leaf, 1)
+    l2 = geometry.edge_length(leaf, 2)
+    l3 = geometry.edge_length(leaf, 3)
+    a = (l0 + l2) / 2
+    b = (l1 + l3) / 2
+    if a <= 0 or b <= 0:
+        return 1.0
+    return max(a, b) / min(a, b)
+
+
+def solve_ratios(
+    root: Node,
+    targets: dict[str, SpaceReq],
+    *,
+    strip: bool = True,
+    perpendicular: bool = True,
+    weight_width: float = 1.0,
+    weight_proportion: float = 0.3,
+    min_width_generic: float = 1.2,
+    max_nfev: int = 4000,
+):
+    """Solve free division ratios in place. Returns the scipy result object.
+
+    ``strip=True`` discards the existing ratios first (start from 0.5) — the
+    honest test that sizes are recoverable from the programme alone.
+
+    ``perpendicular=True`` ties the two ends of each cut (``a == b``), one DOF
+    per branch, so cuts stay perpendicular to their walls (matches Urb's
+    ``perpendicular`` quality and the slicing-tree model). ``min_width_generic``
+    keeps unconstrained circulation/outside leaves from collapsing into slivers.
+    """
+    free = free_branches(root)
+    if not free:
+        return None
+
+    if strip:
+        for b in free:
+            b.division = [0.5, 0.5]
+
+    per = 1 if perpendicular else 2
+    x0 = np.array(
+        [b.division[0] for b in free] if perpendicular
+        else [v for b in free for v in b.division],
+        dtype=float,
+    )
+
+    all_leaves = [leaf for lvl in levels(root) for leaf in lvl.leaves()]
+
+    def apply(x: np.ndarray) -> None:
+        for j, b in enumerate(free):
+            if perpendicular:
+                b.division = [float(x[j]), float(x[j])]
+            else:
+                b.division = [float(x[2 * j]), float(x[2 * j + 1])]
+        geometry.clear_cache()  # divisions changed; invalidate derived coords
+
+    def residuals(x: np.ndarray) -> list[float]:
+        apply(x)
+        r: list[float] = []
+        for leaf in all_leaves:
+            req = targets.get(leaf.type)
+            if req is not None:
+                area = geometry.area(leaf)
+                r.append((area - req.size) / req.size)
+                if weight_width:
+                    w = _width(leaf)
+                    r.append(weight_width * min(0.0, (w - req.width) / req.width))
+                if weight_proportion:
+                    asp = _aspect(leaf)
+                    r.append(weight_proportion * max(0.0, (asp - req.proportion) / req.proportion))
+            elif min_width_generic:
+                # keep circulation/outside from collapsing to slivers
+                w = _width(leaf)
+                r.append(min(0.0, (w - min_width_generic) / min_width_generic))
+        return r
+
+    res = least_squares(
+        residuals, x0, bounds=(_EPS, 1 - _EPS), max_nfev=max_nfev, xtol=1e-10, ftol=1e-10
+    )
+    apply(res.x)
+    return res
+
+
+def area_report(root: Node, targets: dict[str, SpaceReq]) -> str:
+    """Human-readable per-programme-leaf area vs target (for experiment output)."""
+    rows = []
+    for lvl_idx, lvl in enumerate(levels(root)):
+        for leaf in lvl.leaves():
+            if leaf.type in targets:
+                req = targets[leaf.type]
+                a = geometry.area(leaf)
+                rows.append(
+                    f"  {lvl_idx}/{leaf.id:6s} {leaf.type:4s} area={a:6.2f} "
+                    f"target={req.size:6.2f} err={(a - req.size):+6.2f} "
+                    f"w={_width(leaf):.2f}/{req.width:.2f} asp={_aspect(leaf):.2f}/{req.proportion:.2f}"
+                )
+    return "\n".join(rows)