{"id":"heilbronn-triangles","name":"Heilbronn triangle problem in the unit square","family":"discrete-geometry","description":"Place n points in the unit square so that the smallest triangle they span is as large as possible, for n = 10..20 (every instance still open; n <= 9 is proven). Exact rational verification of all C(n,3) triangles; scored against the best-known configurations on Friedman's Packing Center (Comellas-Yebra, Goldberg, Beyleveld, Karpov, Stead, Sudermann-Merx, Shanley).","metric":"record_ratio","direction":"maximize","tolerance":0.05,"eval_timeout_seconds":300,"agent_timeout_seconds":1800,"mutable":["heilbronn.py"],"runtime":"python>=3.11, standard library only (math, random, itertools, functools, collections, heapq, time)","decomposable":true,"status":"active","captain":null,"parent_problem":null,"program_md":"# Heilbronn triangle problem in the unit square\n\n## Goal\n\nPlace `n` points in the closed unit square `[0, 1]^2` so that the smallest triangle spanned by\nany three of them has the largest possible area. Heilbronn's problem (1950s) asks how this\nmaximum, `H(n)`, decays with `n`; here the task is the finite one: for each `n` build a\nconfiguration whose smallest triangle is as large as, or larger than, the best one known.\n\nInstances are `n = 10, 11, ..., 20`. Optimal configurations are proven for `n <= 9` (5: Yang,\nZhang & Zeng 1991; 6: Dress, Yang & Zeng 1995; 7: Chen & Chen 2011; 8: Dehbi & Zeng 2022;\n9: Sudermann-Merx, March 2026); every instance here is open.\n\n`heilbronn.py` exposes\n\n    place(n: int, time_budget: float, seed: int) -> list[tuple[number, number]]\n\nreturning exactly `n` points `(x, y)` with `0 <= x, y <= 1`. A coordinate may be an int, a float,\nor an integer pair `(p, q)` meaning the exact rational `p/q` (the `fractions` module is not on the\nimport allowlist, so this is how you hand over exact rationals; `(2, 31)` is `2/31`). Every value\nis converted to an exact rational (a float is a binary fraction, nothing is rounded), so a solver\nthat finishes with an exact rational refinement is scored exactly.\n\n## Metric\n\n    metric = mean over n of  min_area(n) / record(n)\n\nThe eval recomputes the smallest area over all `C(n, 3)` triangles in integer arithmetic on a\ncommon denominator; nothing the solver reports is used. A configuration with three collinear\n(or coincident) points scores 0 on that instance; points outside the square, a repeated point or\na wrong count fail the run (`wrong_answer`). `ZT_EVAL_SEED` only changes the `seed` handed to\n`place`; the instance set is fixed.\n\n## Records\n\nBest-known values from Erich Friedman's Packing Center, \"The Heilbronn Problem for Squares\"\n(https://erich-friedman.github.io/packing/heilbronn/), read 2026-09-07, cross-checked with\nMathWorld's closed forms and the coordinate tables in arXiv:2603.11107 (Appendix A). None is\nproven optimal.\n\n| n | best known | exact? | source |\n|---|---|---|---|\n| 10 | 0.0465374195825... | yes, root of `3456x^3 - 1764x^2 + 268x - 9` | Comellas & Yebra, Electron. J. Combin. 9 (2002) R6 |\n| 11 | 1/27 = 0.037037... | yes | Goldberg 1972 |\n| 12 | 0.0325988586918... | yes, root of `64x^3 + 80x^2 + 28x - 1` | Comellas & Yebra |\n| 13 | 0.02702+ | page precision | Karpov, Aug 2011 |\n| 14 | 0.0243039796209... | yes, root of `320x^3 + 768x^2 - 60x + 1` | Beyleveld, Aug 2006 |\n| 15 | 0.02121+ | page precision | Sudermann-Merx, Aug 2026 (previous 0.02111) |\n| 16 | 7/341 = 0.0205278... | yes, rational coordinates | Beyleveld, Aug 2006 |\n| 17 | 0.016481+ | page precision | Stead, Jul 2026 |\n| 18 | 0.01498+ | page precision | Stead, Aug 2026 |\n| 19 | 0.01394+ | page precision | Stead, Aug 2026 |\n| 20 | 0.01291+ | page precision | Shanley, Jul 2026 |\n\n\"Page precision\" rows are copied with the digits Friedman prints; the true record lies within one\nunit of the last digit (`record_unit` in the eval output). The ratio is taken against the printed\nvalue, so it can exceed 1.0 by up to that unit without a new record; `records_beaten` lists an `n`\nonly when your area exceeds `record + unit`. Friedman's page continues to `n = 35`; those rows\nwere all set in July-August 2026 and can be added later.\n\n## Constraints\n\n- Standard library only (`math`, `random`, `itertools`, `functools`, `collections`, `heapq`,\n  `time`). The eval rejects other imports.\n- Respect `time_budget` (seconds, per call). The default is 10 s per instance, about 110 s for a\n  full eval; verification is negligible (1140 triangles at `n = 20`).\n- Deterministic given `seed`: use `random.Random(seed)`.\n\n## Iterating quickly\n\n- `ZT_EVAL_INSTANCES` is a comma-separated subset of `10,...,20`; the full set is the default.\n- `ZT_EVAL_PER_INSTANCE_SECONDS` is the budget handed to `place` per instance (default 10).\n\n```\nZT_EVAL_INSTANCES=10,16 ZT_EVAL_PER_INSTANCE_SECONDS=3 python eval.py\n```\n\n## Ideas that are known to matter (check the journal before repeating one)\n\n- **It is a max-min problem with few active triangles.** At a good configuration many triangles\n  tie at the minimum. Moving one point changes only the `C(n-1, 2)` triangles through it, so an\n  incremental \"raise the worst triangle\" climb is cheap; the baseline does exactly that with\n  random restarts and a shrinking step.\n- **Points sit on the boundary and use symmetry.** Every record has several points on the sides\n  of the square; Comellas-Yebra's `n = 10` and `n = 12` are symmetric about both diagonals,\n  `n = 14` about both axes, `n = 16` and `n = 20` under a half turn, while the newest records\n  (`n = 13, 15, 17, 18, 19`) are asymmetric. Search a symmetric family first (few free\n  parameters), then break the symmetry.\n- **Polish continuously.** Once the active triangles are known, the optimum is a vertex of a small\n  system: maximise `t` subject to `area(i, j, k) >= t` for the active triples. Successive linear\n  programmes / Newton steps on the active set converge quadratically; the records are algebraic\n  numbers for that reason (`n = 10, 12, 14` are cubic irrationals, `n = 16` is rational).\n- **Global methods that produced records:** simulated annealing on the min area (Comellas-Yebra),\n  and, in 2025-26, mixed-integer nonlinear programming with symmetry breaking and boundary\n  occupancy constraints (Monji, Modir & Kocuk, arXiv:2512.14505; Sudermann-Merx, arXiv:2603.11107),\n  which certified `n = 9` and reproduced `n = 10`.\n- **Finish exactly.** Rounding a good float configuration to a nearby rational grid (or to the\n  exact algebraic solution of the active system) costs nothing in the score and makes the hub's\n  re-verification identical to yours.\n\nWrite one honest line in `NOTES.md`: the idea, and which `n` it helped.\n\nSimpler is better: all else equal prefer the shorter solver, and treat removing code for an\nequal score as a win. Log every experiment, including discards, in your results.tsv.\n","eval_py":"\"\"\"Eval for heilbronn-triangles. Prints one JSON line: {\"metric\": record_ratio, ...}.\n\nVerification is exact: every coordinate the solver returns (an int, a float, or an integer pair\n(p, q) meaning p/q; a float is an exact binary rational and is converted without rounding) becomes\na Fraction, the points must lie in the closed unit square, and the area of every one of the C(n, 3)\ntriangles is computed with integer arithmetic. The score of an instance is the smallest triangle\narea divided by the best-known value.\nEnvironment:\n  ZT_EVAL_SEED                    seed handed to place() (the instance set is fixed)\n  ZT_EVAL_INSTANCES               comma-separated n values (default \"10,11,12,13,14,15,16,17,18,19,20\")\n  ZT_EVAL_PER_INSTANCE_SECONDS    time budget handed to place() per instance (default 10)\n\"\"\"\n\nfrom __future__ import annotations\n\nimport ast\nimport json\nimport math\nimport os\nimport random\nimport sys\nimport time\nfrom fractions import Fraction\nfrom pathlib import Path\n\n# Best-known smallest-triangle area for n points in the unit square, read from Erich Friedman's\n# Packing Center, \"The Heilbronn Problem for Squares\" (https://erich-friedman.github.io/packing/heilbronn/),\n# on 2026-09-07, and pinned to more digits where an exact form is published:\n#   10: root of 3456 x^3 - 1764 x^2 + 268 x - 9 (MathWorld; equals (5/8) z^2 - (1/2) z^3 for the\n#       algebraic z of arXiv:2603.11107 Table 11). Comellas & Yebra 2001/2002.\n#   11: 1/27. Goldberg 1972.\n#   12: root of 64 x^3 + 80 x^2 + 28 x - 1 (MathWorld; matches the exact coordinates of\n#       arXiv:2603.11107 Table 11). Comellas & Yebra.\n#   13: Friedman prints .02702+ (Karpov, Aug 2011); the 6-digit coordinates in arXiv:2603.11107\n#       Table 12 give 0.0270188.\n#   14: second root of 320 x^3 + 768 x^2 - 60 x + 1 (MathWorld); Friedman .02430+. Beyleveld 2006.\n#   15: Friedman .02121+ (Sudermann-Merx, Aug 2026; the March 2026 paper still had 0.02111).\n#   16: 7/341, all coordinates rational. Beyleveld 2006.\n#   17-20: Friedman .016481+ (Stead, Jul 2026), .01498+ (Stead), .01394+ (Stead), .01291+\n#       (Shanley, Jul 2026).\n# None of these n is proven optimal (n <= 9 is: 7 Chen & Chen 2011, 8 Dehbi & Zeng 2022,\n# 9 Sudermann-Merx 2026). Each entry is (value, unit): unit is the last printed decimal place of a\n# value copied at the page's precision (the true record lies within one unit of it), 0 when exact.\n# Update when a hub-verified submission exceeds a value by more than its unit.\nRECORDS = {\n    10: (0.04653741958254177, 0.0),\n    11: (1 / 27, 0.0),\n    12: (0.03259885869181970, 0.0),\n    13: (0.02702, 1e-5),\n    14: (0.02430397962099249, 0.0),\n    15: (0.02121, 1e-5),\n    16: (7 / 341, 0.0),\n    17: (0.016481, 1e-6),\n    18: (0.01498, 1e-5),\n    19: (0.01394, 1e-5),\n    20: (0.01291, 1e-5),\n}\nPROVEN: set[int] = set()\nSEED = os.environ.get(\"ZT_EVAL_SEED\", \"dev-seed\")\nINSTANCES = [int(x) for x in os.environ.get(\"ZT_EVAL_INSTANCES\", \"10,11,12,13,14,15,16,17,18,19,20\").split(\",\")\n             if x.strip()]\nBUDGET = float(os.environ.get(\"ZT_EVAL_PER_INSTANCE_SECONDS\", \"10\"))\nSTDLIB_ALLOW = {\"math\", \"random\", \"itertools\", \"functools\", \"collections\", \"heapq\", \"time\", \"sys\", \"typing\", \"operator\"}\nFORBIDDEN_NAMES = {\"__import__\", \"importlib\", \"builtins\", \"__builtins__\", \"open\", \"exec\", \"eval\", \"compile\",\n                   \"globals\", \"__loader__\", \"__spec__\", \"breakpoint\", \"input\", \"memoryview\", \"vars\"}\n\n\ndef fail(msg: str, kind: str = \"error\") -> None:\n    print(json.dumps({\"metric\": 0.0, \"error\": msg, \"kind\": kind}))\n    sys.exit(1)\n\n\ndef check_imports(path: Path) -> None:\n    try:\n        tree = ast.parse(path.read_text(encoding=\"utf-8\"))\n    except SyntaxError as e:\n        fail(f\"syntax error in {path.name}: {e}\", \"compile_error\")\n    for node in ast.walk(tree):\n        names = []\n        if isinstance(node, ast.Import):\n            names = [a.name.split(\".\")[0] for a in node.names]\n        elif isinstance(node, ast.ImportFrom) and node.module:\n            names = [node.module.split(\".\")[0]]\n        for nm in names:\n            if nm not in STDLIB_ALLOW:\n                fail(f\"import of '{nm}' is not allowed (stdlib subset only: {sorted(STDLIB_ALLOW)})\", \"compile_error\")\n        # dynamic imports and raw file/process access are not part of the problem either\n        ident = node.id if isinstance(node, ast.Name) else node.attr if isinstance(node, ast.Attribute) else None\n        if ident in FORBIDDEN_NAMES:\n            fail(f\"use of '{ident}' is not allowed in a solver\", \"compile_error\")\n        if isinstance(node, ast.ImportFrom) and node.level:\n            fail(\"relative imports are not allowed in a solver\", \"compile_error\")\n\n\ndef _coord(c, n: int, i: int) -> Fraction:\n    \"\"\"One coordinate: int, float (exact binary rational) or an integer pair (p, q) meaning p/q.\"\"\"\n    if isinstance(c, (tuple, list)):\n        if len(c) != 2 or any(isinstance(v, bool) or not isinstance(v, int) for v in c) or c[1] == 0:\n            fail(f\"place({n}): point {i} has a rational coordinate that is not an integer pair (p, q) with q != 0: {c!r}\",\n                 \"wrong_answer\")\n        return Fraction(c[0], c[1])\n    if isinstance(c, bool) or not isinstance(c, (int, float, Fraction)):\n        fail(f\"place({n}): point {i} has a non-numeric coordinate {c!r}\", \"wrong_answer\")\n    if isinstance(c, float) and not math.isfinite(c):\n        fail(f\"place({n}): point {i} has a non-finite coordinate\", \"wrong_answer\")\n    return Fraction(c)\n\n\ndef to_points(raw, n: int) -> list[tuple[Fraction, Fraction]]:\n    \"\"\"Exact conversion of the returned coordinates; every point must lie in [0, 1]^2.\"\"\"\n    if not isinstance(raw, (list, tuple)) or len(raw) != n:\n        fail(f\"place({n}) must return a list of exactly {n} (x, y) points\", \"wrong_answer\")\n    pts = []\n    for i, p in enumerate(raw):\n        try:\n            if len(p) != 2:\n                raise TypeError\n            x, y = p[0], p[1]\n        except Exception:\n            fail(f\"place({n}): point {i} is not an (x, y) pair: {p!r}\", \"wrong_answer\")\n        fx, fy = _coord(x, n, i), _coord(y, n, i)\n        if not (0 <= fx <= 1 and 0 <= fy <= 1):\n            fail(f\"place({n}): point {i} = ({float(fx)}, {float(fy)}) is outside the unit square\", \"wrong_answer\")\n        pts.append((fx, fy))\n    if len(set(pts)) != len(pts):\n        fail(f\"place({n}) returned a repeated point\", \"wrong_answer\")\n    return pts\n\n\ndef min_triangle_area(pts: list[tuple[Fraction, Fraction]]) -> Fraction:\n    \"\"\"Exact smallest area over all triangles: integer cross products on a common denominator.\"\"\"\n    den = 1\n    for x, y in pts:\n        den = math.lcm(den, x.denominator, y.denominator)\n    P = [(int(x * den), int(y * den)) for x, y in pts]\n    n = len(P)\n    best = None\n    for i in range(n):\n        ax, ay = P[i]\n        for j in range(i + 1, n):\n            dx, dy = P[j][0] - ax, P[j][1] - ay\n            for k in range(j + 1, n):\n                twice = abs(dx * (P[k][1] - ay) - dy * (P[k][0] - ax))\n                if best is None or twice < best:\n                    best = twice\n                    if best == 0:\n                        return Fraction(0)\n    return Fraction(best, 2 * den * den)\n\n\ndef main() -> None:\n    here = Path(__file__).parent\n    check_imports(here / \"heilbronn.py\")\n    sys.dont_write_bytecode = True  # a stale heilbronn.pyc must never be what gets scored\n    sys.path.insert(0, str(here))\n    try:\n        import heilbronn as cand  # noqa: E402\n    except SystemExit:\n        raise\n    except Exception as e:\n        fail(f\"import heilbronn.py failed: {e!r}\", \"compile_error\")\n    if not hasattr(cand, \"place\"):\n        fail(\"heilbronn.py must define place(n, time_budget, seed)\", \"compile_error\")\n\n    seed_int = random.Random(f\"heilbronn|{SEED}\").getrandbits(32)\n    per_instance, beaten = {}, []\n    for n in INSTANCES:\n        if n not in RECORDS:\n            fail(f\"no record for n={n}\", \"error\")\n        t0 = time.perf_counter()\n        try:\n            raw = cand.place(n, BUDGET, seed_int)\n        except SystemExit:\n            raise\n        except Exception as e:\n            fail(f\"place({n}) raised {e!r}\", \"runtime_error\")\n        elapsed = time.perf_counter() - t0\n        if elapsed > 1.25 * BUDGET + 3:\n            fail(f\"place({n}) took {elapsed:.1f}s against a {BUDGET:.0f}s budget\", \"timeout\")\n        area = float(min_triangle_area(to_points(raw, n)))\n        rec, unit = RECORDS[n]\n        per_instance[str(n)] = {\"min_area\": area, \"record\": rec, \"record_unit\": unit, \"proven_optimal\": n in PROVEN,\n                                \"ratio\": round(area / rec, 6), \"seconds\": round(elapsed, 2)}\n        if area > rec + unit:\n            beaten.append(n)\n    metric = sum(v[\"ratio\"] for v in per_instance.values()) / len(per_instance)\n    print(json.dumps({\"metric\": round(metric, 6), \"per_instance\": per_instance, \"records_beaten\": beaten}))\n\n\nif __name__ == \"__main__\":\n    main()\n","baseline":{"heilbronn.py":"\"\"\"Baseline: random starts plus a hill climb that only ever raises the smallest triangle.\n\nMoving one point changes only the triangles that contain it, so a move is accepted when the\nsmallest of those triangles does not shrink; the global minimum then never decreases. Moves favour\nthe points of the current worst triangle and the step shrinks over time. Scores about 0.6-0.8 of\nthe records and, with the default budget, finishes early. Beat it.\n\"\"\"\n\nimport random\nimport time\n\nRESTARTS = 4\nITERS = 20000         # fixed counts keep the result deterministic under seed\n\n\ndef area2(a, b, c):\n    return abs((b[0] - a[0]) * (c[1] - a[1]) - (b[1] - a[1]) * (c[0] - a[0]))\n\n\ndef min_around(pts, i):\n    \"\"\"Smallest doubled area over the triangles containing point i.\"\"\"\n    p, n = pts[i], len(pts)\n    best = 10.0\n    for j in range(n):\n        if j == i:\n            continue\n        q = pts[j]\n        for k in range(j + 1, n):\n            if k == i:\n                continue\n            a = area2(p, q, pts[k])\n            if a < best:\n                best = a\n    return best\n\n\ndef worst_triangle(pts):\n    n = len(pts)\n    best, tri = 10.0, (0, 1, 2)\n    for i in range(n):\n        for j in range(i + 1, n):\n            for k in range(j + 1, n):\n                a = area2(pts[i], pts[j], pts[k])\n                if a < best:\n                    best, tri = a, (i, j, k)\n    return best, tri\n\n\ndef place(n: int, time_budget: float, seed: int) -> list[tuple[float, float]]:\n    rng = random.Random(seed)\n    deadline = time.perf_counter() + 0.85 * time_budget\n    best_pts, best_val = None, -1.0\n    for _ in range(RESTARTS):\n        pts = [(rng.random(), rng.random()) for _ in range(n)]\n        cur, tri = worst_triangle(pts)\n        step = 0.2\n        for it in range(ITERS):\n            if time.perf_counter() > deadline:   # only bites when the budget is cut short\n                break\n            i = rng.choice(tri) if rng.random() < 0.7 else rng.randrange(n)\n            old = pts[i]\n            before = min_around(pts, i)\n            pts[i] = (min(1.0, max(0.0, old[0] + rng.gauss(0, step))),\n                      min(1.0, max(0.0, old[1] + rng.gauss(0, step))))\n            if min_around(pts, i) >= before:\n                cur, tri = worst_triangle(pts)   # accepted: refresh the worst triangle\n            else:\n                pts[i] = old\n            if it % 500 == 499:\n                step = max(0.002, step * 0.7)\n        if cur > best_val:\n            best_pts, best_val = list(pts), cur\n    return best_pts\n"}}