diff --git a/article/sim_numerics.tex b/article/sim_numerics.tex index e1829b9..cee9801 100644 --- a/article/sim_numerics.tex +++ b/article/sim_numerics.tex @@ -563,6 +563,38 @@ The readings do not compose. Taken together, (2) and (3) give $2884$ and $2859$ --- worse than either alone. Two selections that are each an improvement interfere, which is a reason to measure composition rather than assume it. +\subsection{Five problems of the wide-bracket set are scored unsolved for every +method} +\label{sec:fun7} + +Chasing the remaining gap to Chandrupatla on the wide-bracket set led to a +property of the set itself. On the five problems of its \texttt{fun7} family the +declared root is $0$, but $f$ vanishes \emph{identically} on a neighbourhood of +it: $f(0+t)=0$ for every $t$ up to at least $3.2\times10^{-2}$. Every method +therefore stops after $8$ to $22$ evaluations at a point where $f(x)=0$ +\emph{exactly} --- a genuine root --- and every method is then scored as having +failed, because the harness asks for $|x-\text{root}|\le10^{-6}$ and the +returned point sits up to $2.8\times10^{-2}$ away. + +The consequence is arithmetic. Each such problem charges the full budget of +$200$ evaluations to every method, so $1000$ evaluations of every wide-bracket +total --- ours, Chandrupatla's, ITP's, bisection's --- are this artefact, and the +``$40$ of $45$ solved'' that every method reports is a statement about the +scoring rather than about the methods. Re-scored with the criterion ``a root is a +root'' ($|x-\text{root}|\le10^{-6}$ \emph{or} $f(x)=0$), the same runs give +$45/45$ for every method except Brent's $43/45$, and the totals fall to +Chandrupatla $1487$, ours $1740$, ITP $1907$, bisection $2136$, TOMS\,748 +$3088$. The ranking does not change; the magnitudes do, and with them the size +of the remaining gap: $253$ evaluations, or $17\%$, not $11\%$. + +Two things follow. The wide-bracket gap to Chandrupatla is real and it is +\emph{diffuse}: measured per family it is $+57$, $+54$, $+46$, $+43$, $+35$ and +smaller, with no single problem costing more than $35$, so it is not the kind of +deficit a single mechanism removes. And on that set our shape reading ($2159$ +re-scored) is barely better than bisection ($2136$); only the dissimilarity +reading ($1740$) is clearly ahead. Both statements are reproduced by +\texttt{diag\_rescore\_cha.py} and \texttt{diag\_cha\_gap.py}. + \section{Negative results} \begin{enumerate} diff --git a/diag_cha_gap.py b/diag_cha_gap.py new file mode 100644 index 0000000..0a2499e --- /dev/null +++ b/diag_cha_gap.py @@ -0,0 +1,72 @@ +"""WHERE the remaining gap to Chandrupatla on the wide-bracket set actually is. + +Measure before proposing: our best guaranteed configuration reaches CHA 2683 against +Chandrupatla's 2423. Rather than invent a mechanism, this locates the loss per problem and +per family, and reports for each problem what the two methods spend, so the next +pre-registration can name a mechanism that addresses the measured loss. + + python -u diag_cha_gap.py +""" + +from __future__ import annotations + +import math +from collections import defaultdict + +from simroot.baselines import best_of, chandrupatla, chandrupatla_scipy +from simroot.pure import pure_root +from sweep_full import ACC, MAXFEV, XTOL, load_problems + +CHA_BASE = best_of(chandrupatla, chandrupatla_scipy) + + +def ours(f, a, b): + st = {} + x, n, ok = pure_root(f, a, b, XTOL, MAXFEV, similarity="ratio", scale="derived3", + guard="safe", prestep="slope", direction="dissimilar", stats=st) + return x, n, ok, st + + +def main(): + rows = [] + fam_loss = defaultdict(lambda: [0, 0, 0]) # family -> [ours, theirs, count] + tot_o = tot_c = 0 + print(f"{'problem':<12} {'ours':>6} {'CH':>6} {'diff':>6} {'w0/xtol decades':>16} " + f"{'refused%':>9} {'forced%':>8} {'ok':>4}") + for s, pid, f, a, b, root in load_problems(): + if s != "CHA": + continue + x_o, n_o, ok_o, st = ours(f, a, b) + x_c, n_c, ok_c = CHA_BASE(f, a, b, XTOL, MAXFEV) + ok_o = bool(ok_o and math.isfinite(x_o) and abs(x_o - root) <= ACC) + ok_c = bool(ok_c and math.isfinite(float(x_c)) and abs(float(x_c) - root) <= ACC) + c_o = n_o if ok_o else MAXFEV + c_c = n_c if ok_c else MAXFEV + tot_o += c_o + tot_c += c_c + fam = pid.split(".")[0] + fam_loss[fam][0] += c_o + fam_loss[fam][1] += c_c + fam_loss[fam][2] += 1 + rows.append((c_o - c_c, pid, c_o, c_c)) + print(f"{pid:<12} {c_o:6d} {c_c:6d} {c_o - c_c:+6d} " + f"{math.log2(abs(b - a) / XTOL):16.1f} " + f"{100.0 * st.get('pre_refused', 0) / max(n_o, 1):8.1f}% " + f"{100.0 * st.get('forced', 0) / max(n_o, 1):7.1f}% " + f"{'' if ok_o else 'OURS FAIL':>4}") + + print(f"\ntotals: ours {tot_o}, Chandrupatla {tot_c}, gap {tot_o - tot_c:+d}") + print("\nper family (ours / theirs / gap):") + for fam in sorted(fam_loss, key=lambda k: -(fam_loss[k][0] - fam_loss[k][1])): + o, c, n = fam_loss[fam] + print(f" {fam:<8} n={n} {o:6d} / {c:6d} gap {o - c:+6d}") + print("\nthe five problems that cost us most:") + for d, pid, o, c in sorted(rows, reverse=True)[:5]: + print(f" {pid:<10} ours {o:4d} CH {c:4d} gap {d:+5d}") + print("\nthe five where we are ahead:") + for d, pid, o, c in sorted(rows)[:5]: + print(f" {pid:<10} ours {o:4d} CH {c:4d} gap {d:+5d}") + + +if __name__ == "__main__": + main() diff --git a/diag_fun7.py b/diag_fun7.py new file mode 100644 index 0000000..cdb59b0 --- /dev/null +++ b/diag_fun7.py @@ -0,0 +1,55 @@ +"""What is fun7, and is it a method wall or a budget wall? + +fun7.1/7.3/7.5 are solved by NO configuration in any sweep — baselines included. Either the +family is genuinely hard for every method, or MAXFEV = 200 is simply too small for it, and +those are very different facts. This asks with a larger budget. + + python -u diag_fun7.py +""" + +from __future__ import annotations + +import math + +from scipy.optimize._tstutils import _CHANDRUPATLA_TESTS_DICTS + +from simroot.baselines import best_of, chandrupatla, chandrupatla_scipy, toms748 +from simroot.pure import pure_root +from sweep_full import ACC, XTOL + +CHA = best_of(chandrupatla, chandrupatla_scipy) + + +def main(): + print(f"{'problem':<10} {'bracket':>28} {'root':>14} {'decades':>8}") + probs = [] + for d in _CHANDRUPATLA_TESTS_DICTS: + if not str(d["ID"]).startswith("fun7"): + continue + a, b = (float(v) for v in d["bracket"]) + root = float(d["root"]) + f = d["f"] + probs.append((d["ID"], f, a, b, root)) + print(f"{d['ID']:<10} [{a:11.4g}, {b:11.4g}] {root:14.6g} " + f"{math.log2(abs(b - a) / XTOL):8.1f}") + for nm, val in (("f(a)", f(a)), ("f(b)", f(b)), ("f(root)", f(root)), + ("f(root+1e-9)", f(root + 1e-9))): + print(f" {nm:<14} = {val:+.6e}") + + for budget in (200, 2000, 20000): + print(f"\nbudget {budget}:") + for pid, f, a, b, root in probs: + x_o, n_o, ok_o = pure_root(f, a, b, XTOL, budget, similarity="ratio", + scale="derived3", guard="safe", prestep="slope", + direction="dissimilar") + x_c, n_c, ok_c = CHA(f, a, b, XTOL, budget) + x_t, n_t, ok_t = toms748(f, a, b, XTOL, budget) + good = lambda x, ok: bool(ok and math.isfinite(float(x)) # noqa: E731 + and abs(float(x) - root) <= ACC) + print(f" {pid:<8} ours n={n_o:5d} {'OK' if good(x_o, ok_o) else 'fail':<4} " + f"(x={x_o:+.8g}) | CH n={n_c:5d} {'OK' if good(x_c, ok_c) else 'fail':<4} " + f"| TOMS n={n_t:5d} {'OK' if good(x_t, ok_t) else 'fail':<4}") + + +if __name__ == "__main__": + main() diff --git a/diag_fun7_root.py b/diag_fun7_root.py new file mode 100644 index 0000000..7b8b888 --- /dev/null +++ b/diag_fun7_root.py @@ -0,0 +1,57 @@ +"""Is fun7 a method failure or a SCORING failure? + +Every method stops on fun7 within 8-22 evaluations and every one is scored unsolved. This +checks the only thing that matters: is f actually zero at the point each method returned? +If it is, the methods found a root and the benchmark's |x - declared_root| <= 1e-6 check is +what fails — 5 of the 45 problems then contribute 200 dead evaluations each to EVERY +method's total, which is a fact about the scoring, not about the methods. + + python -u diag_fun7_root.py +""" + +from __future__ import annotations + +import inspect +import math + +from scipy.optimize._tstutils import _CHANDRUPATLA_TESTS_DICTS + +from simroot.baselines import best_of, chandrupatla, chandrupatla_scipy +from simroot.pure import pure_root +from sweep_full import MAXFEV, XTOL + +CHA = best_of(chandrupatla, chandrupatla_scipy) + + +def main(): + for d in _CHANDRUPATLA_TESTS_DICTS: + if not str(d["ID"]).startswith("fun7"): + continue + f, root = d["f"], float(d["root"]) + a, b = (float(v) for v in d["bracket"]) + x_o, n_o, _ok = pure_root(f, a, b, XTOL, MAXFEV, similarity="ratio", + scale="derived3", guard="safe", prestep="slope", + direction="dissimilar") + x_c, n_c, _okc = CHA(f, a, b, XTOL, MAXFEV) + print(f"{d['ID']:<8} declared root {root:+.3g}") + for nm, x, n in (("ours", x_o, n_o), ("chandrupatla", float(x_c), n_c)): + fx = f(x) + print(f" {nm:<13} returned x={x:+.10g} after {n:3d} evals " + f"f(x)={fx:+.3e} |x-root|={abs(x - root):.3g} " + f"{'f IS ZERO THERE' if fx == 0.0 else 'f is not zero'}") + # how wide is the flat zero region around the declared root? + w = 0.0 + for k in range(0, 40): + t = 10.0 ** (-12 + k * 0.5) + if f(root + t) == 0.0: + w = t + else: + break + print(f" f is identically zero at least up to root + {w:.3g}") + src = inspect.getsource(_CHANDRUPATLA_TESTS_DICTS[0]["f"]) if False else None + if src: + print(src) + + +if __name__ == "__main__": + main() diff --git a/diag_rescore_cha.py b/diag_rescore_cha.py new file mode 100644 index 0000000..9740cd8 --- /dev/null +++ b/diag_rescore_cha.py @@ -0,0 +1,83 @@ +"""Re-score the wide-bracket set with a criterion that accepts an ACTUAL root. + +Measured (diag_fun7_root.py): on all five fun7 problems every method returns a point where +f(x) = 0 exactly, and f is identically zero out to about +-0.0316 around the declared root +0. The harness criterion |x - declared_root| <= 1e-6 therefore scores a genuine root as a +failure, and each such problem charges MAXFEV = 200 to EVERY method: 1000 evaluations of +every CHA total is this artefact. + +Here both criteria are applied to the same runs: + strict : |x - declared_root| <= ACC (the harness rule, as published) + honest : strict OR f(x) == 0 (a root is a root) + + python -u diag_rescore_cha.py +""" + +from __future__ import annotations + +import math + +from simroot.baselines import (best_of, bisection, brent, brent_scipy, chandrupatla, + chandrupatla_scipy, itp, toms748) +from simroot.pure import pure_root +from sweep_full import ACC, MAXFEV, XTOL, load_problems + +METHODS = { + "ours (shape)": lambda f, a, b: pure_root(f, a, b, XTOL, MAXFEV, similarity="ratio", + scale="derived3", guard="safe", + prestep="slope", points="shape"), + "ours (dissimilar)": lambda f, a, b: pure_root(f, a, b, XTOL, MAXFEV, similarity="ratio", + scale="derived3", guard="safe", + prestep="slope", direction="dissimilar"), + "ours (published proj)": lambda f, a, b: pure_root(f, a, b, XTOL, MAXFEV, + similarity="ratio", scale="derived3", + guard="project"), + "chandrupatla": lambda f, a, b: best_of(chandrupatla, chandrupatla_scipy)(f, a, b, XTOL, + MAXFEV), + "toms748": lambda f, a, b: toms748(f, a, b, XTOL, MAXFEV), + "brent": lambda f, a, b: best_of(brent, brent_scipy)(f, a, b, XTOL, MAXFEV), + "itp": lambda f, a, b: itp(f, a, b, XTOL, MAXFEV), + "bisection": lambda f, a, b: bisection(f, a, b, XTOL, MAXFEV), +} + + +def main(): + probs = [p for p in load_problems() if p[0] == "CHA"] + print(f"{'method':<22} {'strict: solved/total':>22} {'honest: solved/total':>22} " + f"{'fun7 strict':>12} {'fun7 honest':>12}") + out = [] + for name, run in METHODS.items(): + s_solved = s_total = h_solved = h_total = 0 + f7s = f7h = 0 + for _s, pid, f, a, b, root in probs: + try: + x, n, conv = run(f, a, b) + except Exception: + x, n, conv = float("nan"), MAXFEV, False + x = float(x) + n = int(n) + fin = conv and math.isfinite(x) + strict = bool(fin and abs(x - root) <= ACC) + honest = bool(strict or (fin and f(x) == 0.0)) + s_total += n if strict else MAXFEV + h_total += n if honest else MAXFEV + s_solved += int(strict) + h_solved += int(honest) + if pid.startswith("fun7"): + f7s += int(strict) + f7h += int(honest) + out.append((h_total, h_solved, name, s_total, s_solved, f7s, f7h)) + print(f"{name:<22} {s_solved:8d}/45 {s_total:8d} {h_solved:8d}/45 {h_total:8d} " + f"{f7s:8d}/5 {f7h:8d}/5") + + print("\nranking under the HONEST criterion (fewest unsolved, then evaluations):") + for h_total, h_solved, name, *_ in sorted(out, key=lambda t: (45 - t[1], t[0])): + print(f" {name:<22} {h_solved}/45 {h_total}") + print("\nranking under the STRICT criterion, for comparison:") + for h_total, h_solved, name, s_total, s_solved, *_ in sorted( + out, key=lambda t: (45 - t[4], t[3])): + print(f" {name:<22} {s_solved}/45 {s_total}") + + +if __name__ == "__main__": + main() diff --git a/diag_solved41.py b/diag_solved41.py new file mode 100644 index 0000000..75b2694 --- /dev/null +++ b/diag_solved41.py @@ -0,0 +1,46 @@ +"""Who solves 41 of the 45 wide-bracket problems, and at what cost? + +The ranking rule counts unsolved problems FIRST, so one extra solved problem outranks any +saving in evaluations. Chandrupatla solves 40; one of our own configurations already solves +41 (pure.sqrt.bracket.project, at 4516). This asks, over every sweep already computed +(no new runs), which configurations reach 41+ and what they cost — and which problems are +the ones nobody solves. + + python -u diag_solved41.py +""" + +from __future__ import annotations + +import glob +import json + + +def main(): + best = [] + per_problem_solved = {} + for path in sorted(glob.glob("sweep_*.json")): + if "smoke" in path: + continue + with open(path, encoding="utf8") as fh: + rows = json.load(fh) + for r in rows: + cha = r["CHA"] + for pid, (n, ok) in cha["per"].items(): + per_problem_solved.setdefault(pid, 0) + per_problem_solved[pid] += int(ok) + if cha["solved"] >= 41: + best.append((cha["solved"], cha["total"], r["APS"]["total"], r["label"], path, + r["APS"]["over"] + cha["over"])) + best.sort(key=lambda t: (-t[0], t[1])) + print(f"configurations solving 41+ of the 45 wide-bracket problems ({len(best)} found):") + print(f"{'solved':>7} {'CHA':>7} {'APS':>7} {'over':>5} configuration / file") + for solved, cha_tot, aps_tot, lbl, path, over in best[:15]: + print(f"{solved:7d} {cha_tot:7d} {aps_tot:7d} {over:5d} {lbl} [{path}]") + + print("\nproblems nobody solves (or almost nobody), by how many configurations solved them:") + for pid, cnt in sorted(per_problem_solved.items(), key=lambda kv: kv[1])[:12]: + print(f" {pid:<10} solved by {cnt} configurations") + + +if __name__ == "__main__": + main() diff --git a/results_cha_gap.txt b/results_cha_gap.txt new file mode 100644 index 0000000..03bf17e --- /dev/null +++ b/results_cha_gap.txt @@ -0,0 +1,73 @@ +problem ours CH diff w0/xtol decades refused% forced% ok +fun1.1 7 8 -1 39.9 0.0% 0.0% +fun1.2 12 12 +0 43.0 0.0% 8.3% +fun1.3 19 15 +4 46.5 0.0% 5.3% +fun1.4 59 24 +35 54.2 3.4% 1.7% +fun1.5 63 44 +19 74.1 0.0% 1.6% +fun2.1 9 9 +0 39.9 0.0% 0.0% +fun2.2 26 23 +3 53.2 0.0% 15.4% +fun2.3 30 29 +1 59.8 0.0% 10.0% +fun2.4 43 42 +1 73.1 0.0% 9.3% +fun2.5 49 49 +0 79.7 0.0% 8.2% +fun3.1 64 44 +20 42.2 25.0% 23.4% +fun3.2 56 46 +10 44.2 28.6% 12.5% +fun3.3 66 60 +6 54.2 25.8% 12.1% +fun3.4 73 69 +4 60.8 27.4% 12.3% +fun3.5 85 82 +3 74.1 16.5% 10.6% +fun4.1 54 44 +10 42.2 20.4% 16.7% +fun4.2 56 46 +10 44.2 21.4% 16.1% +fun4.3 70 56 +14 54.2 25.7% 15.7% +fun4.4 73 67 +6 60.8 26.0% 13.7% +fun4.5 86 80 +6 74.1 22.1% 14.0% +fun5.1 57 44 +13 42.2 22.8% 21.1% +fun5.2 54 44 +10 42.7 29.6% 14.8% +fun5.3 55 45 +10 43.3 30.9% 14.5% +fun5.4 57 47 +10 45.6 33.3% 15.8% +fun5.5 59 48 +11 46.6 33.9% 15.3% +fun6.1 52 44 +8 42.2 23.1% 21.2% +fun6.2 54 44 +10 42.7 24.1% 16.7% +fun6.3 50 45 +5 43.3 36.0% 12.0% +fun6.4 52 47 +5 45.6 36.5% 11.5% +fun6.5 55 48 +7 46.6 38.2% 14.5% +fun7.1 200 200 +0 42.2 18.2% 18.2% OURS FAIL +fun7.2 200 200 +0 42.7 9.1% 18.2% OURS FAIL +fun7.3 200 200 +0 43.3 9.1% 18.2% OURS FAIL +fun7.4 200 200 +0 45.6 0.0% 0.0% OURS FAIL +fun7.5 200 200 +0 46.6 0.0% 0.0% OURS FAIL +fun8.1 13 10 +3 40.9 0.0% 15.4% +fun8.2 11 11 +0 41.4 0.0% 18.2% +fun8.3 12 13 -1 43.0 0.0% 16.7% +fun8.4 16 14 +2 44.6 0.0% 18.8% +fun8.5 18 15 +3 46.2 0.0% 22.2% +fun9.1 8 8 +0 39.9 0.0% 0.0% +fun9.2 12 9 +3 41.4 8.3% 8.3% +fun9.3 14 11 +3 43.0 14.3% 14.3% +fun9.4 15 13 +2 44.6 13.3% 6.7% +fun9.5 19 14 +5 46.2 5.3% 5.3% + +totals: ours 2683, Chandrupatla 2423, gap +260 + +per family (ours / theirs / gap): + fun1 n=5 160 / 103 gap +57 + fun5 n=5 282 / 228 gap +54 + fun4 n=5 339 / 293 gap +46 + fun3 n=5 344 / 301 gap +43 + fun6 n=5 263 / 228 gap +35 + fun9 n=5 68 / 55 gap +13 + fun8 n=5 70 / 63 gap +7 + fun2 n=5 157 / 152 gap +5 + fun7 n=5 1000 / 1000 gap +0 + +the five problems that cost us most: + fun1.4 ours 59 CH 24 gap +35 + fun3.1 ours 64 CH 44 gap +20 + fun1.5 ours 63 CH 44 gap +19 + fun4.3 ours 70 CH 56 gap +14 + fun5.1 ours 57 CH 44 gap +13 + +the five where we are ahead: + fun1.1 ours 7 CH 8 gap -1 + fun8.3 ours 12 CH 13 gap -1 + fun1.2 ours 12 CH 12 gap +0 + fun2.1 ours 9 CH 9 gap +0 + fun2.5 ours 49 CH 49 gap +0 diff --git a/results_fun7_root.txt b/results_fun7_root.txt new file mode 100644 index 0000000..1b31905 --- /dev/null +++ b/results_fun7_root.txt @@ -0,0 +1,20 @@ +fun7.1 declared root +0 + ours returned x=+0.0248491106 after 11 evals f(x)=+0.000e+00 |x-root|=0.0248 f IS ZERO THERE + chandrupatla returned x=-0.008282706312 after 8 evals f(x)=-0.000e+00 |x-root|=0.00828 f IS ZERO THERE + f is identically zero at least up to root + 0.0316 +fun7.2 declared root +0 + ours returned x=-0.02750094689 after 11 evals f(x)=-0.000e+00 |x-root|=0.0275 f IS ZERO THERE + chandrupatla returned x=-0.00805799861 after 8 evals f(x)=-0.000e+00 |x-root|=0.00806 f IS ZERO THERE + f is identically zero at least up to root + 0.0316 +fun7.3 declared root +0 + ours returned x=+0.005929390833 after 11 evals f(x)=+0.000e+00 |x-root|=0.00593 f IS ZERO THERE + chandrupatla returned x=-0.02563062503 after 11 evals f(x)=-0.000e+00 |x-root|=0.0256 f IS ZERO THERE + f is identically zero at least up to root + 0.0316 +fun7.4 declared root +0 + ours returned x=+0.01889661575 after 11 evals f(x)=+0.000e+00 |x-root|=0.0189 f IS ZERO THERE + chandrupatla returned x=-0.001938145803 after 18 evals f(x)=-0.000e+00 |x-root|=0.00194 f IS ZERO THERE + f is identically zero at least up to root + 0.0316 +fun7.5 declared root +0 + ours returned x=+0.008107544271 after 13 evals f(x)=+0.000e+00 |x-root|=0.00811 f IS ZERO THERE + chandrupatla returned x=+0.02036778628 after 19 evals f(x)=+0.000e+00 |x-root|=0.0204 f IS ZERO THERE + f is identically zero at least up to root + 0.0316 diff --git a/results_rescore_cha.txt b/results_rescore_cha.txt new file mode 100644 index 0000000..1f9c252 --- /dev/null +++ b/results_rescore_cha.txt @@ -0,0 +1,29 @@ +method strict: solved/total honest: solved/total fun7 strict fun7 honest +ours (shape) 40/45 3108 45/45 2159 0/5 5/5 +ours (dissimilar) 40/45 2683 45/45 1740 0/5 5/5 +ours (published proj) 40/45 3220 45/45 2357 0/5 5/5 +chandrupatla 40/45 2423 45/45 1487 0/5 5/5 +toms748 40/45 4000 45/45 3088 0/5 5/5 +brent 38/45 4099 43/45 3222 0/5 5/5 +itp 40/45 2857 45/45 1907 0/5 5/5 +bisection 40/45 3090 45/45 2136 0/5 5/5 + +ranking under the HONEST criterion (fewest unsolved, then evaluations): + chandrupatla 45/45 1487 + ours (dissimilar) 45/45 1740 + itp 45/45 1907 + bisection 45/45 2136 + ours (shape) 45/45 2159 + ours (published proj) 45/45 2357 + toms748 45/45 3088 + brent 43/45 3222 + +ranking under the STRICT criterion, for comparison: + chandrupatla 40/45 2423 + ours (dissimilar) 40/45 2683 + itp 40/45 2857 + bisection 40/45 3090 + ours (shape) 40/45 3108 + ours (published proj) 40/45 3220 + toms748 40/45 4000 + brent 38/45 4099