order_b.py
Cited on: SIZE.
To run it: save order_b.py, with erring_cells.py, erring_label.py, erring_lower.py, erring_order.py, erring_widths.py, erring_zones.py beside it; needs numpy.
"""
order_b.py -- does putting the widths out of order buy the gap sum B a
constant at large m, at a fixed multiset of widths.
QUESTION. R(P) = sum_a T_aa + B/2 <= E(P) at every P, B = sum_j (S0_j +
o_j) (erring_order.py's DERIVATION (3)). The route to lim kappa_m log m
= 1/2 needs B >= (1 - o(1))/log m at every P once max w and the prefix
are spent. B reads ratios of widths only and never S_m, so at a fixed
multiset S_m/8 is fixed and order moves R alone. At the multiset of a
sorted near-optimum, does some order lower R by a ratio that stays away
from 1 as m grows?
DESIGN. The exponential profile w(t) ~ e^(Kt) cut into m intervals,
c_j = -log(1 - j (1 - e^-K)/m)/K, K the one minimizing max(R(sort P),
S_m/8) over 40 values in (0, 3 log m] (B falls and S_m rises as K
grows). Orders: (O1) sorted; (O2) every move of the interval at
position j to position i, j and i on a 16 x 16 grid of positions (i < j
plants a spike earlier, i > j a dip later); (O3) a greedy descent from
sorted, each round taking the best move of the O2 grid, 10 rounds or
until none lowers R; (O4) sawtooths: the ranks = r mod N in sorted
order, r = 0..N-1, concatenated, N teeth each spanning the whole ratio,
N = 2, 4, 8; (O5) front-loading: remove the n intervals of ranks
m - t - n .. m - t - 1, insert them in increasing width at positions
round(i P m / n), i = 0..n-1, of the rest (where two positions round
alike the later insert lands first; no reported optimum collides), n in {1, 2, 4, 8, 16, 32, 64}
(n <= m/4), t/m in {0, 0.02, 0.05, 0.1}, P in {0.05, 0.1, 0.2, 0.3,
0.5}, rho_F the least R over them / R(sort). At m = 250, 500, 1000,
2000. Ingredient print at the best O5 order, never a kill clause: the
width moved, and over the intervals with c_(j-1) >= delta = 1/log m the
count of H (B_j < theta w_j, theta = 2/log m) and Lambda_H / log m,
Lambda_H = sum over H of w_j^2/B_j. (EDGE, added after the first
run) the least R over O5 past its grid: P in {0, 0.01, 0.02, 0.05}, n
in {16, 32, 64, 128, 256} (n <= m/4), t/m in {0, 0.01, 0.02}, at m =
500 and 2000. (E, added then) the multiset and the O5 grid at m = 120,
250, 500: E1, E at R's best O5 order over E(sort); E2, the least E over
the O5 grid over E(sort). The frozen design named a 40 x 40
grid and 30 rounds; cut to 16 x 16 and 10 before the first run, for
wall-clock (each R at m = 2000 costs about 0.05 s). C3 reads the move
back directly, not through the grid. The frozen family was geometric in
the index, w_j ~ q^j; a rehearsal at m = 150 showed it is a linear
profile in position (w ~ t), where B stays a constant and the balance
ran to the grid's edge (S_m/8 never bound), so the exponential profile,
the sorted proof's own corner, replaced it before the first full run. A
move is taken when it lowers R by more than 1e-12 R, not 1e-15: the
rehearsal's descent moved microscopic intervals at rounding level. The
rehearsal (m = 150, 400) found O3 never stalls in 10 rounds, each round
moving an interval of rank about 0.93m into the first 0.07m to 0.35m of
positions (R / R(sort) 0.958 and 0.974 after 10), so a capped descent
cannot read the limit; O5 and its kill were frozen then, before the
full run, and the frozen kill on O3 and on R_delta (R over the
intervals with c_(j-1) >= 1/log m) was retired with it (a moved
interval leaves R_delta's set with its share).
PREDICTIONS (written before the run). (Q1) No single move (O2) lowers R
by more than 1% at any m read. (Q2) O3's descent stalls within 2% of
sorted. (Q3) O4's sawtooth raises R by a factor growing with N. (Q4)
The moves that do lower R put a narrow interval first or move one
within the first 1/log m of [0, 1). (Q5, at the re-freeze) rho_F rises
toward 1 with m, a prefix effect. (Q6, before the E arm) E1 and E2 at
least 0.995 at every m read: the gain is R's slack. (Q7, before EDGE)
(1 - rho) log m moves by at most 0.03 from O5's at m = 500 and 2000.
KILL (of "order buys B no constant at large m"; re-frozen after the
rehearsal). rho_F <= 0.95 at each of m = 500, 1000, 2000, rising across
them by no more than 0.01. (E) of "front-loading buys E nothing": E2
at most 0.98 at m = 250 and 500. (EDGE) of "the gain is about c/log
m": (1 - rho) log m at m = 2000 above that at m = 500 by more than
0.02.
CONTROLS (read before any verdict). C1 R <= LB2' <= E at 200 random
partitions, m = 2..12, and the fast diagonal against tents()' trace.
C2 B by hand at the multiset {0.4, 0.6}: 0.18 sorted, 0.08 reversed.
C3 the move machinery: at m = 200, sorted with the interval of rank 180
moved to position 1, the move back gives R(sort) to 1e-15. C4 (with
the E arm) E against E(H_n) = (n^2 + 3n - 1)/(8 n^2), n = 50 and 200.
FINDINGS (copied from the printed run; m = 250 to 2000). Controls 3/3:
C1 least margin 1.37e-07, diagonal against tents 5.55e-17; C2 B 0.18
and 0.08, R 0.150 sorted and 0.140 reversed; C3 back to R(sort)
exactly. KILL missed: rho_F 0.9404, 0.9469, 0.9527 at m = 500, 1000,
2000 (0.9312 at 250), rising 0.012 across them. Front-loading is what
order buys: the best O5 order moves the n just below the 2% widest
(n = 8, 16, 32, 64 at m = 250 to 2000) among the first 5% of
positions, total width 0.138, 0.126, 0.114, 0.104, each below 1/log m
(0.181, 0.161, 0.145, 0.132), the prefix the route spends; (1 - rho_F)
log m reads 0.380, 0.370, 0.367, 0.360 on the grid, but its rate is
unsettled (EDGE below). The O5 optimum sits at the grid's least P at every m and at
its largest n at m = 2000, so rho_F bounds the family's least from
above. R(sort) log m 0.3127, 0.3283, 0.3413, 0.3525 (K 6.63 to 9.12,
max w 0.21 to 0.19; S_m/8 0.0076 to 0.0056, never binding). (Q1)
holds: the best single move lowers R by 0.62, 0.29, 0.13, 0.06%. (Q2)
is wrong: O3 does not stall, 0.963, 0.978, 0.989, 0.995 after 10
rounds; rounds 1 and 10, the ones printed, move a wide interval (rank
about 0.93m) to position 0.13m to 0.27m. (Q3) holds weakly: sawtooth ratios grow with
N at every m (N = 2: 1.0025 to 1.0937; N = 8: 1.22 to 1.46). (Q4) is
wrong as written: the gaining moves carry a wide interval to an early
position, not a narrow one, and no position c is printed. (Q5) holds
at the m read, its limit open (EDGE). Ingredients:
Lambda_H/log m at sorted 1.2143, 1.2039, 1.1965, 1.1889 (|H| 75, 150,
301, 602: with theta = 2/log m above 1/K, H is every interval past
delta), at the O5 order 1.3483 to 1.2821. Run: 160.6 s, peak commit
149.6 MB.
EDGE and E (frozen after the first run, both arms added then; the
second run reprints the first's lines unchanged). EDGE: m = 500 keeps
O5's optimum, rho 0.940393, (1 - rho) log m 0.3704; m = 2000 moves to
(n, t/m, P) = (32, 0.01, 0.02), rho 0.947715, 0.3974: the kill of "the
gain is about c/log m" FIRED (0.027 > 0.02); (Q7) is wrong on its own
terms, m = 2000 moving 0.038 from O5's 0.360; the gain's rate is
unsettled at these m: 1 - rho falls 12% from m = 500 to
2000 where 1/log m falls 18%. E: control C4 1.86e-14. At R's best O5
order E is 2.0010, 1.8708, 1.7819 times E(sort) at m = 120, 250, 500,
and no order of the O5 grid lowers E at all (E2 = 1, its optimum
sorted): front-loading is R's slack. (Q6) holds;
the kill of "front-loading buys E nothing" missed. E(sort) log m
0.4071, 0.4078, 0.4215. Run 191.6 s, peak commit 149.8 MB.
"""
import os
os.environ.setdefault("OPENBLAS_NUM_THREADS", "1")
import math # noqa: E402
import time # noqa: E402
import numpy as np # noqa: E402
from erring_cells import tents # noqa: E402
from erring_order import E, bsum, cuts, lb2p # noqa: E402
CHECKS = []
def check(name, ok, detail=""):
CHECKS.append(bool(ok))
tail = f" ({detail})" if detail else ""
print(f" {'PASS' if ok else 'FAIL'} {name}{tail}")
def bterms(c):
"""B_j per interval j (0-based, I_j = [c[j], c[j+1]))."""
w = np.diff(c)
m = len(w)
idx = np.arange(m)
lag = (idx[None, :] > idx[:, None]) & (c[1:][None, :] <= 2 * c[:-1][:, None])
pair = w[:, None] * np.minimum(w[:, None], w[None, :])
s0 = np.where(lag, pair, 0.0).sum(axis=0)
own = np.where(2 * c[:-1] >= w, w * w / 2, 0.0)
return s0 + own
def diag(c):
"""T_aa per interval, exact: h_a is piecewise linear, so a trapezoid
on its breakpoints integrates it exactly."""
w = np.diff(c)
bd = c[1:]
out = np.empty(len(w))
for a in range(len(w)):
r = w[a] / 2
t0 = 2 * c[a] + r
t1 = t0 + w[a]
lo = np.searchsorted(bd, t0 - r, side="left")
hi = np.searchsorted(bd, t1 + r, side="right")
bs = bd[lo:hi]
if len(bs) == 0:
out[a] = 0.0
continue
pts = np.concatenate([[t0, t1], bs, bs - r, bs + r,
(bs[1:] + bs[:-1]) / 2])
pts = np.unique(np.clip(pts, t0, t1))
h = np.clip(r - np.abs(pts[:, None] - bs[None, :]), 0, None).max(1)
out[a] = float(((h[1:] + h[:-1]) / 2 * np.diff(pts)).sum())
return out
def Rof(w):
c = cuts(w)
return float(diag(c).sum() + bterms(c).sum() / 2)
def profile(m, K):
c = -np.log1p(-np.arange(m + 1) * (-np.expm1(-K)) / m) / K
c[-1] = 1.0
return np.diff(c)
def balanced(m):
best = None
for lr in np.linspace(3 * math.log(m) / 40, 3 * math.log(m), 40):
w = profile(m, lr)
f = max(Rof(w), float((w * w).sum()) / 8)
if best is None or f < best[0]:
best = (f, lr, w)
return best
def moved(order, j, i):
o = list(order)
x = o.pop(j)
o.insert(i, x)
return o
def grid(m, n=16):
return sorted(set(int(round(v)) for v in np.linspace(0, m - 1, n)))
def best_move(ws, order, f0):
pos = grid(len(order))
best = (f0, None)
for j in pos:
for i in pos:
if i == j:
continue
o = moved(order, j, i)
f = Rof(ws[o])
if f < best[0] - 1e-12 * f0:
best = (f, (j, i, o))
return best
def section_controls(rng):
print("CONTROLS")
worst = math.inf
dd = 0.0
for _ in range(200):
m = int(rng.integers(2, 13))
w = rng.dirichlet(np.ones(m) * rng.choice([0.2, 1.0, 5.0]))
r = Rof(w)
worst = min(worst, lb2p(w) - r, E(w) - lb2p(w))
dd = max(dd, abs(diag(cuts(w)).sum() - np.trace(tents(cuts(w)))))
check("C1 R <= LB2' <= E, m = 2..12, and the diagonal against tents",
worst >= -1e-12 and dd < 1e-12,
f"200 partitions, least margin {worst:.2e}, diag diff {dd:.2e}")
bs, br = bsum([0.4, 0.6]), bsum([0.6, 0.4])
check("C2 B at {0.4, 0.6}: 0.18 sorted, 0.08 reversed",
abs(bs - 0.18) < 1e-12 and abs(br - 0.08) < 1e-12,
f"B sorted {bs:.6f}, reversed {br:.6f}; R sorted "
f"{Rof(np.array([0.4, 0.6])):.6f}, reversed "
f"{Rof(np.array([0.6, 0.4])):.6f}")
m = 200
_, _, ws = balanced(m)
base = Rof(ws)
o = moved(list(range(m)), 180, 0)
back = Rof(ws[moved(o, 0, 180)])
check("C3 the move back restores R(sort), m = 200",
abs(back - base) < 1e-15,
f"R(sort) {base:.12f}, moved {Rof(ws[o]):.12f}, back {back:.12f}")
def front(m, n, t, P):
top = m - t
moved_ = list(range(top - n, top))
rest = [k for k in range(m) if not (top - n <= k < top)]
for i, k in enumerate(moved_):
rest.insert(int(round(i * P * m / n)), k)
return rest
def ingredients(w, delta, theta):
c = cuts(w)
w = np.diff(c)
b = bterms(c)
h = (c[:-1] >= delta) & (b < theta * w)
lam = float((w[h] ** 2 / np.maximum(b[h], 1e-300)).sum())
return int(h.sum()), lam
def section_m(m):
t0 = time.time()
fbal, K, ws = balanced(m)
order = list(range(m))
L = math.log(m)
r0 = Rof(ws)
s8 = float((ws * ws).sum()) / 8
nh, lam = ingredients(ws, 1 / L, 2 / L)
print(f"m {m}: K {K:.3f}, max w {ws[-1]:.4f}, R(sort) {r0:.6f}, "
f"S_m/8 {s8:.6f}, R log m {r0 * L:.4f}; sorted |H| {nh}, "
f"Lambda_H/log m {lam / L:.4f}")
for N in (2, 4, 8):
o = [k for r in range(N) for k in range(r, m, N)]
print(f" O4 sawtooth N = {N}: R / R(sort) {Rof(ws[o]) / r0:.4f}")
best = (r0, None, order)
for n in (1, 2, 4, 8, 16, 32, 64):
if n > m // 4:
continue
for tf in (0.0, 0.02, 0.05, 0.1):
for P in (0.05, 0.1, 0.2, 0.3, 0.5):
o = front(m, n, int(round(tf * m)), P)
f = Rof(ws[o])
if f < best[0]:
best = (f, (n, tf, P), o)
fF, par, oF = best
nh, lam = ingredients(ws[oF], 1 / L, 2 / L)
moved_w = 0.0 if par is None else float(
ws[m - int(round(par[1] * m)) - par[0]:
m - int(round(par[1] * m))].sum())
print(f" O5 front-loading: rho_F {fF / r0:.6f} at (n, t/m, P) {par}, "
f"width moved {moved_w:.4f}, R log m {fF * L:.4f}; |H| {nh}, "
f"Lambda_H/log m {lam / L:.4f}")
f, cur = r0, order
for rnd in range(10):
fb, mv = best_move(ws, cur, f)
if mv is None:
print(f" O3 round {rnd + 1}: no move lowers R")
break
j, i, cur = mv
f = fb
if rnd in (0, 9):
print(f" O3 round {rnd + 1}: move {j} -> {i} (width rank "
f"{cur[i]}), R / R(sort) {f / r0:.6f}")
print(f" O3: rho_B {f / r0:.6f} [{time.time() - t0:.1f} s]")
return fF / r0
def o5_best(f, ws, Ps, ns, ts):
m = len(ws)
best = (f(ws), None, list(range(m)))
for n in ns:
if n > m // 4:
continue
for tf in ts:
for P in Ps:
o = front(m, n, int(round(tf * m)), P)
v = f(ws[o])
if v < best[0]:
best = (v, (n, tf, P), o)
return best
def section_edge():
print("EDGE: the O5 family past its grid's edge")
out = []
for m in (500, 2000):
_, _, ws = balanced(m)
r0 = Rof(ws)
f, par, _ = o5_best(Rof, ws, (0.0, 0.01, 0.02, 0.05),
(16, 32, 64, 128, 256), (0.0, 0.01, 0.02))
L = math.log(m)
out.append((1 - f / r0) * L)
print(f" m {m}: rho {f / r0:.6f} at (n, t/m, P) {par}, "
f"(1 - rho) log m {(1 - f / r0) * L:.4f}")
print(f" KILL of 'the gain is about c/log m': "
f"{'FIRED' if out[1] - out[0] > 0.02 else 'missed'}")
def section_e():
print("E: does front-loading lower the error itself")
worst = 0.0
for n in (50, 200):
w = np.concatenate([[0.5], np.ones(n - 1), [0.5]])
worst = max(worst, abs(E(w / w.sum()) - (n * n + 3 * n - 1)
/ (8 * n * n)))
check("C4 E against E(H_n), n = 50 and 200", worst < 1e-12,
f"max |diff| {worst:.2e}")
grid5 = ((0.05, 0.1, 0.2, 0.3, 0.5), (1, 2, 4, 8, 16, 32, 64),
(0.0, 0.02, 0.05, 0.1))
rows = []
for m in (120, 250, 500):
t0 = time.time()
_, _, ws = balanced(m)
e0 = E(ws)
_, parR, oR = o5_best(Rof, ws, *grid5)
e1 = E(ws[oR])
fE, parE, _ = o5_best(E, ws, *grid5)
rows.append((m, e1 / e0, fE / e0))
print(f" m {m}: E(sort) {e0:.6f}, E log m {e0 * math.log(m):.4f}; "
f"E1 at R's O5 order {parR}: {e1 / e0:.6f}; E2 least over "
f"O5 at {parE}: {fE / e0:.6f} [{time.time() - t0:.1f} s]")
kill = all(r[2] <= 0.98 for r in rows if r[0] in (250, 500))
print(f" KILL of 'front-loading buys E nothing': "
f"{'FIRED' if kill else 'missed'}")
def main():
rng = np.random.default_rng(1694)
section_controls(rng)
rows = [(m, section_m(m)) for m in (250, 500, 1000, 2000)]
last = rows[1:]
kill = (all(r[1] <= 0.95 for r in last)
and last[-1][1] - last[0][1] <= 0.01)
print(f"KILL of 'order buys B no constant': "
f"{'FIRED' if kill else 'missed'} (rho_F at m = "
f"{[r[0] for r in last]}: {[round(r[1], 4) for r in last]})")
section_edge()
section_e()
print(f"checks: {sum(CHECKS)}/{len(CHECKS)}")
if __name__ == "__main__":
main()