#!/usr/bin/env python3 """Fixing a real bug in my own filter, and rescuing the conclusion honestly. symn.py reported that both surviving operators have a "rational solution" of sym^4 with exponents (5/3, 4/3, 5/3). Those are not integers, so t^(5/3) is not a rational function and the label was wrong. The ansatz allowed a_s = k(1-e_s) for k = 0..n and never required a_s in Z -- correct for a HYPEREXPONENTIAL solution (which is what the reducibility test wanted, where I first wrote this loop) and wrong for a RATIONAL one. Fourth instance in two days of the same species of bug: answering a nearby question and having it look like the real one. Three things here, in order: 1. The strict test: same search, a_s restricted to integers. This is the real "does sym^n have a rational solution" question. Note it is not vacuous -- for e_s = -2/3, k(1-e_s) = 5k/3 is an integer at k = 3 -- so at n = 4 the integral options are {0, 5} rather than {0, 5/3, 10/3, 5, 20/3}. 2. Independent verification of the hyperexponential hit. I extract w and check sum_k C_k u^(k) = 0 by float Taylor series around a point, which shares no code with the exact linear solver that found it. If the recursion in symn.py were wrong, the exact check would agree with itself and this would not. 3. The conclusion the hit actually supports. If u = prod (t-s)^{a_s} w solves sym^4 L with every a_s in (1/3)Z, then u^3 is a genuine RATIONAL function, and u^3 is a sum of products of twelve solutions of L, so u^3 is a rational solution of sym^12 L. An irreducible order-2 operator with a rational solution of any symmetric power has finite projective monodromy: the zero divisor of that solution is a monodromy-invariant finite set of points of P^1, and a subgroup of PGL_2 preserving one is finite unless it fixes a point (reducible -- excluded by classC2.py) or swaps a pair (dihedral -- excluded because sym^2 has no solution at all, integral or not). So the solutions are algebraic and neither operator is a counterexample. The headline survives; the reason for it changes. """ import json from fractions import Fraction as F from classC import kernel from classC2 import pmul, padd, pscal, pder, parts from symn import binom, sym_power_operator def search(coeffs, exps, lam, n, maxdeg=10, integral_only=False, want_w=False): d, Pn, Qn, (e0, e1, el), L, (t, t1, tl) = parts(exps, lam) dd = pder(d) N = len(coeffs) - 1 opts = [] for e in (e0, e1, el): o = [k * (1 - e) for k in range(n + 1)] if integral_only: o = [a for a in o if a.denominator == 1] opts.append(sorted(set(o))) for a0 in opts[0]: for a1 in opts[1]: for al in opts[2]: Rn = padd(pscal(pmul(t1, tl), a0), pscal(pmul(t, tl), a1), pscal(pmul(t, t1), al)) S = [[F(1)]] for j in range(N): S.append(padd(pmul(pder(S[j]), d), pscal(pmul(S[j], dd), F(-j)), pmul(S[j], Rn))) dpow = [[F(1)]] for _ in range(N + 1): dpow.append(pmul(dpow[-1], d)) for D in range(0, maxdeg + 1): ncols = D + 1 rows = 2 * (max(len(c) for c in coeffs) + len(dpow[N]) + D + 4) M = [[F(0)] * ncols for _ in range(rows)] for k in range(ncols): wk = [F(0)] * ncols wk[k] = F(1) wd = [wk] for _ in range(N): wd.append(pder(wd[-1])) col = [] for kk in range(N + 1): U = [] for i in range(kk + 1): U = padd(U or [F(0)], pscal(pmul(pmul(wd[i], S[kk - i]), dpow[i]), F(binom(kk, i)))) col = padd(col or [F(0)], pmul(pmul(coeffs[kk], U), dpow[N - kk])) for i, v in enumerate(col): if i < rows: M[i][k] = v v = kernel(M, ncols) if v: return ((a0, a1, al), D, list(v)) if want_w else ((a0, a1, al), D) return None def taylor_verify(coeffs, a, lam, w, t0=F(37, 10), order=None): """Independent check by float Taylor series: does u = prod (t-s)^a_s * w satisfy sum_k C_k u^(k) = 0 at t0? Series arithmetic in plain floats, no exact solver.""" N = len(coeffs) - 1 M = N + 2 sings = [0.0, 1.0, float(F(*lam))] x0 = float(t0) def mul(p, q): out = [0.0] * M for i, x in enumerate(p): if x: for j, y in enumerate(q): if y and i + j < M: out[i + j] += x * y return out u = [0.0] * M u[0] = 1.0 for s, aa in zip(sings, a): base = x0 - s av = float(aa) # (base + h)^av = base^av * sum_k C(av,k) (h/base)^k ser = [0.0] * M c = base ** av ser[0] = c for k in range(1, M): c = c * (av - (k - 1)) / (k * base) ser[k] = c u = mul(u, ser) wser = [0.0] * M for i, cw in enumerate(w): # w(x0+h) = sum_i cw * (x0+h)^i for k in range(min(i, M - 1) + 1): wser[k] += float(cw) * binom(i, k) * (x0 ** (i - k)) u = mul(u, wser) fact = [1.0] for k in range(1, M): fact.append(fact[-1] * k) total, scale = 0.0, 0.0 for k, Ck in enumerate(coeffs): cv = sum(float(x) * x0 ** i for i, x in enumerate(Ck)) term = cv * u[k] * fact[k] total += term scale = max(scale, abs(term)) return total, scale if __name__ == "__main__": still = json.load(open('/tmp/stillC.json')) for (exps, lam) in still: exps = tuple(tuple(x) for x in exps) lam = tuple(lam) e = " ".join(f"{a}/{b}" for (a, b) in exps) print(f"\n=== lam={lam[0]}/{lam[1]} [{e}]", flush=True) for n in (2, 4, 6): C, order = sym_power_operator(exps, lam, n) if C is None: print(f" sym^{n}: construction failed: {order}", flush=True) continue strict = search(C, exps, lam, n, maxdeg=8, integral_only=True) print(f" sym^{n} (order {order}) STRICT rational, integer exponents only: " f"{'FOUND ' + str(strict) if strict else 'none up to poly degree 8'}", flush=True) if n == 4: hyp = search(C, exps, lam, n, maxdeg=8, want_w=True) if hyp: (a, D, w) = hyp print(f" sym^4 hyperexponential solution: exponents {tuple(str(x) for x in a)}, " f"w degree <= {D}", flush=True) print(f" w coefficients: {[str(x) for x in w]}", flush=True) tot, sc = taylor_verify(C, a, lam, w) ok = sc > 0 and abs(tot) < 1e-6 * sc print(f" independent float-Taylor check at t=3.7: residual {tot:.3e} " f"vs largest term {sc:.3e} -> {'VERIFIED' if ok else 'DOES NOT VERIFY'}", flush=True) cubed = [3 * x for x in a] print(f" 3*exponents = {tuple(str(x) for x in cubed)} -> u^3 is rational, " f"and u^3 solves sym^12 -> finite projective monodromy, solutions algebraic", flush=True)