Spaces:
Sleeping
Sleeping
| import math | |
| import copy | |
| import json | |
| import numpy as np | |
| from fractions import Fraction | |
| from typing import Dict, Tuple, Optional, Any, Callable | |
| from dataclasses import dataclass, field | |
| # ══════════════════════════════════════════════════════════════════════ | |
| # SECTION A: HEISENBERG QUANTUM SIMULATOR (CHP Stabilizer State) | |
| # ══════════════════════════════════════════════════════════════════════ | |
| class StabilizerTableau: | |
| n: int | |
| tableau: np.ndarray = field(default=None) | |
| def __post_init__(self): | |
| if self.tableau is None: | |
| self.tableau = np.zeros((2*self.n, 2*self.n+1), dtype=np.uint8) | |
| for i in range(self.n): | |
| self.tableau[i, i] = 1; self.tableau[self.n+i, self.n+i] = 1 | |
| class HeisenbergSimulator: | |
| def __init__(self, n_qubits: int): | |
| self.n = n_qubits | |
| self.state = StabilizerTableau(n_qubits) | |
| self.gate_count = 0 | |
| def copy(self): | |
| new_sim = HeisenbergSimulator(self.n) | |
| new_sim.state = self.state.copy() | |
| new_sim.gate_count = self.gate_count | |
| return new_sim | |
| def h(self, q: int): | |
| self.gate_count += 1 | |
| t = self.state.tableau | |
| for i in range(2*self.n): | |
| t[i, 2*self.n] ^= t[i, q] & t[i, self.n+q] | |
| t[i, q], t[i, self.n+q] = t[i, self.n+q], t[i, q] | |
| def cnot(self, c: int, target: int): | |
| self.gate_count += 1 | |
| n = self.n; t = self.state.tableau | |
| for i in range(2*n): | |
| t[i, 2*n] ^= t[i,c] & t[i,n+target] & (t[i,target] ^ t[i,n+c] ^ 1) | |
| t[i, target] ^= t[i, c]; t[i, n+c] ^= t[i, n+target] | |
| def get_operator_expectation(self, pauli_string: str) -> float: | |
| n = self.n | |
| x_part, z_part = np.zeros(n, dtype=np.uint8), np.zeros(n, dtype=np.uint8) | |
| phase = 0 | |
| for i, p in enumerate(pauli_string): | |
| if p == 'X': x_part[i] = 1 | |
| elif p == 'Z': z_part[i] = 1 | |
| elif p == 'Y': x_part[i] = 1; z_part[i] = 1; phase += 1 | |
| t = self.state.tableau | |
| for i in range(n, 2*n): | |
| if (np.sum(t[i, :n] * z_part + t[i, n:2*n] * x_part) % 2) == 1: return 0.0 | |
| test_row = np.concatenate([x_part, z_part, [phase % 2]]) | |
| for i in range(n, 2*n): | |
| if np.array_equal(t[i, :2*n], test_row[:2*n]): return 1.0 if t[i, 2*n] == 0 else -1.0 | |
| neg_row = test_row.copy(); neg_row[2*n] = 1 - neg_row[2*n] | |
| if np.array_equal(t[i, :2*n], neg_row[:2*n]): return -1.0 if t[i, 2*n] == 0 else 1.0 | |
| return 1.0 | |
| # ══════════════════════════════════════════════════════════════════════ | |
| # SECTION B: THE TRUTH ORACLE (Ground Truth Verification) | |
| # ══════════════════════════════════════════════════════════════════════ | |
| FACTORIZATION_TEST_CASES = [(15, 3, 5), (35, 5, 7), (85, 5, 17), (143, 11, 13)] | |
| def extract_factors_from_binding(binding: Dict[str, float], N: int) -> Optional[Tuple[int, int]]: | |
| for pname, qname in [('p', 'q'), ('factor_p', 'factor_q'), ('x', 'y')]: | |
| if pname in binding and qname in binding: | |
| p, q = round(binding[pname]), round(binding[qname]) | |
| if p * q == N and p >= 2 and q >= 2: return (min(p, q), max(p, q)) | |
| if 's' in binding and 'd' in binding: | |
| s, d = binding['s'], binding['d'] | |
| p, q = round((s + d) / 2), round((s - d) / 2) | |
| if p * q == N and p >= 2 and q >= 2: return (min(p, q), max(p, q)) | |
| if 't' in binding: | |
| t = round(binding['t']) | |
| if t >= 2 and N % t == 0: return (min(t, N//t), max(t, N//t)) | |
| for v, val in binding.items(): | |
| if not math.isfinite(val): continue | |
| p = round(val) | |
| if p >= 2 and p * p <= N and N % p == 0: return (p, N // p) | |
| return None | |
| def ground_truth_verify(binding: Dict[str, float], N: int) -> Tuple[bool, str]: | |
| result = extract_factors_from_binding(binding, N) | |
| if result is None: return False, "No valid integer factor pair extractable from binding" | |
| p, q = result | |
| if p * q == N: return True, f"{p} × {q} = {N} ✓" | |
| return False, f"Round-trip failed: extracted ({p}, {q}), {p}×{q}={p*q} ≠ {N}" | |
| class AdversaryPreFlight: | |
| def __init__(self, solver_callback: Callable, log_callback: Callable): | |
| self.solver_callback = solver_callback | |
| self.log = log_callback | |
| self.test_cases = FACTORIZATION_TEST_CASES | |
| self.rejection_log = [] | |
| def check_isomorphism(self, axl_def: Any) -> Tuple[bool, str]: | |
| self.log("ADVERSARY PRE-FLIGHT: testing isomorphism on known toy cases...") | |
| passed_count = 0 | |
| for N, true_p, true_q in self.test_cases: | |
| new_axl = copy.deepcopy(axl_def) | |
| if new_axl.observations: | |
| old_N = list(axl_def.observations.values())[0] | |
| for key in list(new_axl.observations.keys()): | |
| if abs(new_axl.observations[key] - old_N) < 1e-6: new_axl.observations[key] = float(N) | |
| for v in new_axl.variables: | |
| if abs(v['lo'] - v['hi']) < 1e-6 and abs(v['lo'] - old_N) < 1e-6: v['lo'] = v['hi'] = float(N) | |
| binding = self.solver_callback(new_axl, N) | |
| if not binding: | |
| self.log(f" N={N}: ✗ FAIL (Exact Integer Solver returned empty binding)") | |
| continue | |
| ok, msg = ground_truth_verify(binding, N) | |
| self.log(f" N={N}: {'✓ PASS' if ok else '✗ FAIL'} | {msg}") | |
| if ok: passed_count += 1 | |
| else: self.rejection_log.append({'N': N, 'reason': msg}) | |
| if passed_count == len(self.test_cases): return True, f"Passed all {passed_count} pre-flight cases" | |
| elif passed_count == 0: return False, "REJECTED: failed all pre-flight cases — isomorphism invalid" | |
| else: return False, f"REJECTED: passed {passed_count}/{len(self.test_cases)} — not reliable enough" | |
| def failure_taxonomy(self, binding: Dict[str, float], N: int, ce: float) -> str: | |
| ok, _ = ground_truth_verify(binding, N) | |
| if ok: return "NONE — solved correctly" | |
| if ce > 100: return "TYPE_A: Structural infeasibility — constraints contradict, redesign isomorphism" | |
| result = extract_factors_from_binding(binding, N) | |
| if result is None: | |
| if any(abs(v - round(v)) > 0.3 for v in binding.values() if math.isfinite(v)): | |
| return "TYPE_D: Integrality failure — optimizer stuck between integers" | |
| return "TYPE_C: Wrong abstraction — inverse map broken, can't extract factors" | |
| p, q = result | |
| if p * q != N: | |
| if abs(p * q - N) / N < 0.1: return "TYPE_D: Precision issue — close but not exact" | |
| return "TYPE_C: Wrong abstraction — extracted integers don't multiply to N" | |
| return "TYPE_B: Optimization landscape — increase rays or change initialization" | |
| # ══════════════════════════════════════════════════════════════════════ | |
| # SECTION C: SHOR DEMO & PHYSICS BENCHMARKS (No missing import bugs!) | |
| # ══════════════════════════════════════════════════════════════════════ | |
| def run_shors_demo(N: int): | |
| print(f"\n=== Shor Demo: N={N} ===") | |
| a = 2 | |
| r = 1 | |
| while pow(a, r, N) != 1: r += 1 | |
| print(f" Period r={r} (a={a} mod N={N})") | |
| sim = HeisenbergSimulator(r) | |
| for i in range(r-1): | |
| sim.h(i) | |
| sim.cnot(i, i+1) | |
| zz = sim.get_operator_expectation('Z' * r) | |
| print(f" ZZ...Z expectation = {zz:.4f}") | |
| print(f" Variables tracked = {r} (vs 2^{r} = {2**r} Schrödinger amplitudes)") | |
| half = pow(a, r//2, N) | |
| f1 = math.gcd(half - 1, N) | |
| f2 = math.gcd(half + 1, N) | |
| print(f" Factors: gcd(a^(r/2)±1, N) = {f1}, {f2}") | |
| print(f" Verify: {f1}×{f2}={f1*f2} {'✓' if f1*f2==N else '✗'}") | |
| def benchmark_heisenberg_vs_schrodinger(n_max: int): | |
| print(f"\n=== Heisenberg vs Schrödinger Scaling ===") | |
| for n in range(2, n_max+1, 2): | |
| heisenberg_vars = 2*n | |
| schrodinger_dim = 2**n | |
| compression = schrodinger_dim // heisenberg_vars | |
| print(f" n={n:3d}: Heisenberg={heisenberg_vars:6d} vars | " | |
| f"Schrödinger=2^{n} | Compression={compression:,}x") | |
| # ══════════════════════════════════════════════════════════════════════ | |
| # SECTION D: SYSTEM PROMPTS & STRUCTURAL TEMPLATES | |
| # ══════════════════════════════════════════════════════════════════════ | |
| NOVEL_ISOMORPHISM_DISCOVERY_PROMPT = """\ | |
| You are an isomorphism designer. Your task is to invent a continuous optimization encoding for integer factorization of N. | |
| HARD REQUIREMENTS: | |
| 1. INVERSE MAP: Write f() and g() explicitly. g() must output exact integers. | |
| 2. MANUAL VERIFICATION: Show your objective evaluates to 0 for N=15, 35, 85 at the exact factor values. | |
| 3. LANDSCAPE: Explain local minima and why descent won't get trapped. | |
| 4. UNIQUENESS: Are the factors the ONLY zeros in [2, sqrt(N)]? | |
| 5. PSL TEMPLATE: Write the full ready-to-run block with tight bounds and GEQ guards. | |
| Do NOT propose an isomorphism if you cannot answer these 5 constraints. | |
| """ | |
| ENCODER_PROMPT = """You are the Encoder. CRITICAL RULE: Every problem MUST include an explicit INVERSE MAP in the description. | |
| OUTPUT SCHEMA (JSON only): | |
| {"name": string, "description": string, "variables": [{"name": string, "lo": number, "hi": number, "is_int": boolean}], "constraints": [{"kind": "EQ"|"GEQ"|"LEQ", "expr": string, "weight": number}], "anchors": [], "observations": {"variable_name": number}} | |
| RULES: Use ** for exponents. Variables > 10000 encode in log2 space.""" | |
| COLLAPSER_PROMPT = """You are the Grounded Hypothesis Collapser. | |
| HARD RULE: For factorization problems, ALWAYS compute the integers (p, q) and verify p*q=N. | |
| OUTPUT SCHEMA (JSON only): | |
| {"scaled_formulas": {"metric": "expr"}, "hypothesis_markdown": "## Markdown output"}""" | |
| STRUCTURAL_HYPOTHESIS_TEMPLATES = { | |
| "Fermat Algebraic (Exact Bijection)": """\ | |
| CONCRETE HYPOTHESIS: Fermat Algebraic Factorization | |
| ==================================================== | |
| N = p*q is equivalent to 4N = s^2 - d^2 | |
| where s = p+q and d = p-q. | |
| Factors: p = (s+d)/2, q = (s-d)/2. | |
| MINIMIZE: d | |
| Variables: | |
| s in [2.0, 172.0] | |
| d in [0.0, 84.0] | |
| sum_sq in [340.0, 340.0] | |
| p in [1.0, 85.0] | |
| q in [1.0, 85.0] | |
| Constraints: | |
| EQ weight=100: sum_sq - (s**2 - d**2) | |
| EQ weight=100: sum_sq - 340.0 | |
| EQ weight=10: p - (s + d) / 2.0 | |
| EQ weight=10: q - (s - d) / 2.0 | |
| GEQ: s - d - 2.0 | |
| GEQ: d - 0.0 | |
| GEQ: s - 2.0 | |
| Anchor: sum_sq - 340.0 = 0 (tolerance 0.001). | |
| Observation: sum_sq = 340.0. | |
| """, | |
| "1D Sine Snap (Fixed Parameterization)": """\ | |
| CONCRETE HYPOTHESIS: 1D Sine-Snap Factorization (Fixed) | |
| ========================================================= | |
| Key fix: only ONE free variable t. y is derived as N/t. | |
| Eliminates y=0 collapse entirely. | |
| MINIMIZE: snap | |
| Variables: | |
| t in [2.0, 9.0] | |
| snap in [0.0, 2.0] | |
| Constraints: | |
| EQ weight=10: snap - (sin(t * 3.14159265)**2 + sin(85.0 / t * 3.14159265)**2) | |
| GEQ: t - 2.0 | |
| LEQ: t - 9.0 | |
| Anchor: snap = 0 (tolerance 0.01). | |
| Observation: snap = 0.0. | |
| """, | |
| "Shor Period-Finding: Heisenberg vs Schrodinger Resource Proof": """\ | |
| CONCRETE HYPOTHESIS: Shor Period-Finding Resource Comparison | |
| ============================================================ | |
| Prove Heisenberg picture uses O(n^2) bits vs Schrodinger 2^n complex numbers. | |
| Variables: | |
| n_qubits in [4.0, 50.0] | |
| log2_heisenberg_bits in [0.0, 14.0] | |
| log2_schrodinger_bits in [4.0, 54.0] | |
| log2_compression in [0.0, 50.0] | |
| clifford_gate_count in [6.0, 1300.0] | |
| n_operators_tracked in [8.0, 100.0] | |
| Constraints: | |
| EQ weight=10: log2_heisenberg_bits - log(2*n_qubits*(2*n_qubits+1))/log(2) | |
| EQ weight=10: log2_schrodinger_bits - n_qubits - log(16)/log(2) | |
| EQ weight=10: log2_compression - log2_schrodinger_bits + log2_heisenberg_bits | |
| EQ weight=10: n_operators_tracked - 2*n_qubits | |
| EQ weight=10: clifford_gate_count - n_qubits*(n_qubits+1)/2 | |
| GEQ: log2_compression - 1.0 | |
| GEQ: n_qubits - 4.0 | |
| Anchor: log2_compression >= 1.0 (mode=geq, tolerance 0.5). | |
| Observation: n_qubits = 20. | |
| """ | |
| } |