""" integral.py -- the carry bound at the d the real support and the pair graph leave open (11, 12, 16, 20), by the carries' integrality. QUESTION. Over every bounded relaxed path from a start (x integral, delta_t = x_t - (x_(t+1) + ... + x_(t+d)) in {-1, 0, 1}, first window sum 0), does a carry x_t or a window sum x_(t+1) + ... + x_(t+d) reach C_d + 1, C_d = max(1, floor((d - 1)/4))? The real support bounds both by L1 norms H_x, H_s, which pass C_d + 1 at d = 11, 12, 16, 20 (3.04, 3.24, 4.11, 5.03); there the integrality of x is put back near t. BY HAND. x_t = g^(t) . delta (real_support.py), g^(t)_k = h_(k - t) + beta^(-t) gam_k, gam = -sum_(s < d) h_(. - s)/(1 + ... + beta^(-(d-1))). For each t a feasibility MILP asks for z_t >= C_d + 1 (or z_t - delta_t, the window sum). Its columns: delta_m integral on [lo_m, F), real in [-1, 1] on [F, F2), F = 10d, F2 past where h's forward tail falls below 1e-3; z_u integral in [-X, X], X = floor(H_x) = C_d + 1 at every open d. Its rows: |row_u . delta - z_u| <= tau_u + 1e-8, tau_u the mass of row u outside the columns (and its coefficients below 1e-9, moved there for conditioning), and z_u - (z_(u+1) + ... + z_(u+d)) = delta_u wherever the window lies inside. h is cut where its terms reach 1e-16 forward and at 200d back, a lost mass below 1e-13, inside the 1e-8. STATIONARY (every t >= B = 5 at once): coordinates relative to t, m in [-B, F2), rows u in [-B + 1, F - d), row u being h_(m - u) plus beta^(-u) w, w = beta^(-t) gam . delta the beta-mode as one real column, cut loose from gam . delta in |w| <= beta^(-B) d L/(1 + ... + beta^(-(d-1))), since ||gam||_1 <= d L/(1 + ... + beta^(-(d-1))) and |delta| <= 1. At t >= B its rows sit at absolute times t + u >= 1 and its deltas at t + m >= 0, so every path's window near t is a point of it. B is a design variable: a B-point restricted to m >= -B', rows u >= -B' + 1 is a B'-point for B' < B (the dropped columns' mass goes to the larger tau, and |w| shrinks with B), so a smaller B is weaker; B = 5 keeps the verdicts of B = 10 and halves the absolute arms. ABSOLUTE (one t < B at a time): rows g^(u), u in [0, F - d), k in [0, F2), the first window sum 0. The mass of g^(u) past F2 is h's forward tail T(F2 - u) plus beta^(-u) gam's, at most beta^(-u) d T(F2 - d + 1), T(n) = sum_(m >= n) |h_m|; both are at most T(F2 - u - d), so tau_u = T(F2 - u - d)(1 + d beta^(-u)) at u >= 0. Only the positive side is asked: x -> -x flips delta, keeps the first window sum 0, and every box is symmetric. Infeasible (status 2) at every arm of d proves the bound at d. PREDICTIONS, CONTROL AND KILL (fixed before the runs). R1 at B = 10, d = 11, 12, the stationary carry and window sum infeasible (TRANSPLANT from the stationary optimum without w at B = 30, 2.0001); R2 every absolute t < 10 infeasible there; R3 likewise at d = 16, 20 (TRANSPLANT from 11, 12: the real support's slack there 0.11 and 0.03). KILL of the route at d: a feasible print at that d. CONTROL, read first: at d = 13, target 3, the stationary carry and window sum, the absolute t = 3 carry and t = 2 window sum all feasible, the triangle wave being a point of each (carry 3 at t = 3 + 12j, window sum 3 at t = 2 + 12j). Fixed after R1..R3's run, for the bare run: T1 at B = 5, presolve off, the stationary arms infeasible at d = 11, 12, 16, 20 (TRANSPLANT from B = 10); T2 the absolute t = 0..4 infeasible there, one solve each; T3 the CONTROL at B = 5, presolve off; T4 at d = 11 the stationary carry at target 3 as an LP (integrality dropped) feasible (property, from B = 10's). KILL of B = 5: a stationary arm feasible or unread at any of the four d (then B = 10 kept); KILL of the fit: the bare run past 900 s. FINDINGS (copied from the runs' prints). - CONTROL holds: all four arms feasible at d = 13, presolve on at B = 10 and presolve off at B = 5 (T3; the bare run's 26.3, 63.3, 30.1, 11.1 s). - R1, R2, R3 hold at B = 10, every arm infeasible on two solves (HiGHS's presolve on and off). T1, T2 hold at B = 5, one solve each. The bound C_d holds at d = 11, 12, 16, 20 (C_d + 1 = 3, 3, 4, 5) over every bounded relaxed path, legal or not, to HiGHS's floating verdicts (below). The kills are silent; the bare run takes 621.4 s. - Integrality refutes what the reals admit: at d = 11 the stationary carry's LP at target 3 is feasible (T4, at B = 5 and B = 10), the MILP is not; at t = 0 the LP alone is infeasible. At B = 10 the same rows at target 2 are feasible (stationary carry and sum, t = 5 carry and sum; at t = 0 not). - A first form, one MILP for every t < 69 with a binary per (t, kind) and their sum >= 1, found no point of d = 13 in 600 s (its big-M relaxation is weak); one t at a time finds it. A search for eventually periodic paths (one MILP per period) found the wave's C_d at d = 9, 13, and found no lasso of 3 at d = 12 to n = 48, its 30 s time-outs unread; the stationary optimum without w (valid at t >= 69) reads 2.000116 at d = 11. Neither is needed by the proof, and both are cut. - Float: at B = 10 HiGHS shifts bounds (values to 1.6e14) in 3 arms at d = 12, each in the presolve-on solve (presolve off: no shift), and in 3 arms at d = 16, the solve and size unread; at B = 5, in the presolve-off solve of d = 16's t = 4 window sum (values to -8.3e10; presolve on: no shift, infeasible), and none at 11, 12 or 20. Every arm solved twice agrees. The verdicts are HiGHS's floating ones (its bound shifts at feasibility 1e-6 times 1 plus a random fraction, as the shift prints read; the rows' tau + 1e-8), no exact certificate. RUN RECORD: the bare run (B = 5, presolve off) 621.4 s, peak commit 349.7 MB, 6/6 checks, each arm within 0.3 s of the runs before it, d = 16's shift as recorded. Before it, at B = 5, one process each under a 512 MB ceiling: d = 13 130.4 s, 214.9 MB; d = 11, 12 246.4 s, 185.8 MB (d = 12's stationary arms 114.5 and 78.2 s); d = 16 110.5 s, 341.9 MB; d = 20 134.0 s, 297.6 MB. At B = 10, each arm solved twice but d = 13's (presolve on), the four runs under a 512 MB ceiling: d = 13 130.1 s, 211.9 MB; TLIM=300 d = 11, 12 343.5 s, 179.0 MB; TLIM=300 d = 16 705.3 s, 300.4 MB; TLIM=600 d = 20 520.4 s, 447.4 MB; the d = 11 LP and target-2 checks 684.3 s, 192.2 MB; the d = 12 presolve-off and B = 15, 20 checks (infeasible, B = 20 without a shift) 272.1 s, 193.8 MB. """ import math import os import sys import time os.environ.setdefault("OPENBLAS_NUM_THREADS", "1") import numpy as np # noqa: E402 from scipy.optimize import Bounds, LinearConstraint, milp # noqa: E402 from scipy.sparse import coo_matrix # noqa: E402 from real_support import by_roots, norms # noqa: E402 TOL = 1e-8 B = 5 CHECKS = [] def check(name, ok, detail=""): CHECKS.append(bool(ok)) print(f" [{'PASS' if ok else 'FAIL'}] {name}" + (f" -- {detail}" if detail else ""), flush=True) def C_of(d): return max(1, (d - 1) // 4) def model(d, t, B=B, F=None, F2=None): """The rows at time t (None: stationary, every t >= B), sparse: a dict with the coordinate list A of (row, col, value), lo, hi, the delta count M (the first F - lo_m integral), the z columns zi, the box X, the row u0 of z_t, and ext, the beta-mode w's column (stationary).""" hf, hb = by_roots(d) L, Ls, nx, ns = norms(d, hf, hb) X = int(math.floor(max(nx.max(), L + 1e-12))) F = F or 10 * d if F2 is None: # forward tail past F2 - F below 1e-3 tail = np.cumsum(np.abs(hf)[::-1])[::-1] F2 = min(len(hf) - 1, F + int(np.argmax(tail < 1e-3))) beta = hb[1] / hb[2] w = t is None if w: # h_(m - u), m in [-B, F2) lo_m = -B def row(u): m = np.arange(lo_m, F2) j = m - u v = np.where(j >= 0, hf[np.clip(j, 0, len(hf) - 1)], hb[np.clip(-j, 0, len(hb) - 1)]) v[j >= len(hf)] = 0 v[-j >= len(hb)] = 0 full = np.abs(hf).sum() + np.abs(hb[1:]).sum() return v, full - np.abs(v).sum() us = range(-B + 1, F - d) u0 = 0 else: # g^(u)_k, k in [0, F2) lo_m = 0 K = np.arange(F2) def shifted(s): j = K - s return np.where(j >= 0, hf[np.clip(j, 0, len(hf) - 1)], hb[np.clip(-j, 0, len(hb) - 1)]) gam = -sum(shifted(s) for s in range(d)) / sum(beta ** -s for s in range(d)) def row(u): v = shifted(u) + beta ** -u * gam # mass of g^(u) beyond F2: h's forward tail past F2 - u, plus # gam's, bounded by d times h's tail past F2 - u - d tail = np.abs(hf[max(0, F2 - u - d):]).sum() * (1 + d * beta ** -u) return v, tail us = range(0, F - d) u0 = t M = F2 - lo_m # delta variables us = list(us) zi = {u: M + i for i, u in enumerate(us)} wc = M + len(zi) if w else None wb = beta ** -B * d * L / sum(beta ** -s for s in range(d)) A, lo, hi = [], [], [] r = 0 for u in us: v, tau = row(u) small = np.abs(v) < 1e-9 # into the slack, for conditioning tau += np.abs(v[small]).sum() A += [(r, int(c), float(v[c])) for c in np.flatnonzero(~small)] A.append((r, zi[u], -1.0)) if w: A.append((r, wc, beta ** -u)) lo.append(-tau - TOL) hi.append(tau + TOL) r += 1 for u in us: # z_u - sum z_(u+i) = delta_u if u + d in zi and lo_m <= u < F: A.append((r, zi[u], 1.0)) A += [(r, zi[u + i], -1.0) for i in range(1, d + 1)] A.append((r, u - lo_m, -1.0)) lo.append(0) hi.append(0) r += 1 if not w: # the first window sum 0 A += [(r, zi[u], 1.0) for u in range(d)] lo.append(0) hi.append(0) r += 1 return dict(A=A, lo=lo, hi=hi, M=M, F=F, lo_m=lo_m, zi=zi, X=X, u0=u0, ext=[(-wb, wb, 0)] if w else []) def solve(mo, c, tlim, presolve, lp): """HiGHS on the model's rows and its ext columns (lp: no integrality).""" M, zi, X = mo["M"], mo["zi"], mo["X"] nz = len(zi) ext = mo["ext"] nv = M + nz + len(ext) rr, cc, vv = zip(*mo["A"]) A = coo_matrix((vv, (rr, cc)), shape=(len(mo["lo"]), nv)).tocsr() integ = np.zeros(nv) if not lp: integ[:mo["F"] - mo["lo_m"]] = 1 integ[M:M + nz] = 1 lb = [-1] * M + [-X] * nz + [e[0] for e in ext] ub = [1] * M + [X] * nz + [e[1] for e in ext] return milp(c=c, integrality=integ, bounds=Bounds(lb, ub), constraints=LinearConstraint(A, mo["lo"], mo["hi"]), options={"time_limit": tlim, "presolve": presolve}) def feasible(d, t, kind, target=None, B=B, tlim=600, presolve=False, lp=False): """A carry (or window sum) of target at t: t None reads every t >= B at once, else the absolute problem at t. Zero objective; HiGHS's status, 2 (infeasible) proving the bound there.""" target = target or C_of(d) + 1 mo = model(d, t, B=B) r = len(mo["lo"]) mo["A"].append((r, mo["zi"][mo["u0"]], 1.0)) if kind == "sum": mo["A"].append((r, mo["u0"] - mo["lo_m"], -1.0)) mo["lo"].append(target) mo["hi"].append(np.inf) nv = mo["M"] + len(mo["zi"]) + len(mo["ext"]) res = solve(mo, np.zeros(nv), tlim, presolve, lp) return res.status def arm(d, t, kind, t0, **kw): t1 = time.time() st = feasible(d, t, kind, **kw) print(f" d = {d}, t {'>= %d' % B if t is None else '= %d' % t}, " f"{kind} >= {kw.get('target') or C_of(d) + 1}" f"{' (LP)' if kw.get('lp') else ''}: status {st} " f"{time.time() - t1:.1f} s [{time.time() - t0:.1f} s]", flush=True) return st def section_control(t0): print("CONTROL d = 13, target 3: the wave a point of every arm", flush=True) sts = [arm(13, t, kind, t0, target=3) for t, kind in ((None, "carry"), (None, "sum"), (3, "carry"), (2, "sum"))] check("every control arm feasible", all(s == 0 for s in sts), f"statuses {sts}") def section_f(t0): print(f"F every arm at d = 11, 12, 16, 20: t >= {B}, then t < {B}", flush=True) for d in (11, 12, 16, 20): sts = [arm(d, t, kind, t0) for t in [None] + list(range(B)) for kind in ("carry", "sum")] check(f"d = {d}: no carry or window sum of {C_of(d) + 1}", all(s == 2 for s in sts), f"statuses {sorted(set(sts))}") if d == 11: st = arm(11, None, "carry", t0, lp=True) check("d = 11: the stationary carry's LP at 3 feasible " "(integrality refutes it, not the box)", st == 0) def main(): t0 = time.time() section_control(t0) section_f(t0) print(f"{sum(CHECKS)}/{len(CHECKS)} checks passed " f"[{time.time() - t0:.1f} s]") sys.exit(0 if all(CHECKS) else 1) if __name__ == "__main__": main()