|
| 1 | +""" |
| 2 | +M15 core -- class/representative machinery. |
| 3 | +
|
| 4 | +Builds on milestone13's identity_complement (complement spectra, orbits, |
| 5 | +deformation). Adds: |
| 6 | + - ADE + affine-A + unicyclic graph builders |
| 7 | + - the complement-eigenvector CONNECTION (Procrustes transport over shared |
| 8 | + support) and its cycle holonomy -- the genuinely non-exact object. |
| 9 | + (Scalar/vector spectral differences are potentials -> exact -> zero |
| 10 | + holonomy identically; M15 founding journal, exp_01 registration.) |
| 11 | +
|
| 12 | +Gauge note: each vertex's frame is the eigenvector matrix of its complement |
| 13 | +subgraph, computed once (deterministic eigh). Holonomy is defined up to |
| 14 | +conjugation by the start vertex's frame; its eigenvalue ANGLES and the |
| 15 | +Frobenius deficit ||H - I|| are conjugation-invariant up to the registered |
| 16 | +tolerance, and labeling-invariance is tested explicitly, not assumed. |
| 17 | +""" |
| 18 | + |
| 19 | +import sys |
| 20 | +import numpy as np |
| 21 | +from pathlib import Path |
| 22 | + |
| 23 | +_M13_CORE = Path(__file__).resolve().parent.parent.parent / "milestone13" / "core" |
| 24 | +sys.path.insert(0, str(_M13_CORE)) |
| 25 | +from identity_complement import ( # noqa: E402 |
| 26 | + PHI, INV_PHI, LN_PHI, |
| 27 | + complement_spectrum, vertex_orbits, |
| 28 | + complement_deformation_rate, max_deformation_rate, |
| 29 | + find_shortest_path, _convert_numpy, |
| 30 | +) |
| 31 | + |
| 32 | +RESULTS_DIR = Path(__file__).resolve().parent.parent / "results" |
| 33 | + |
| 34 | + |
| 35 | +def save_m15_results(experiment_name, data): |
| 36 | + import json |
| 37 | + from datetime import datetime |
| 38 | + RESULTS_DIR.mkdir(exist_ok=True) |
| 39 | + ts = datetime.now().strftime("%Y%m%d_%H%M%S") |
| 40 | + out = RESULTS_DIR / f"{experiment_name}_{ts}.json" |
| 41 | + with open(out, 'w') as f: |
| 42 | + json.dump(data, f, indent=2, default=str) |
| 43 | + print(f"\n Results saved: {out}") |
| 44 | + return out |
| 45 | + |
| 46 | + |
| 47 | +# ============================================================ |
| 48 | +# Graph builders |
| 49 | +# ============================================================ |
| 50 | + |
| 51 | +def build_path(n): |
| 52 | + """A_n Dynkin diagram: path on n vertices.""" |
| 53 | + a = np.zeros((n, n)) |
| 54 | + for i in range(n - 1): |
| 55 | + a[i, i + 1] = a[i + 1, i] = 1.0 |
| 56 | + return a |
| 57 | + |
| 58 | + |
| 59 | +def build_d(n): |
| 60 | + """D_n: path on n-1 vertices with an extra leaf on vertex 1.""" |
| 61 | + a = np.zeros((n, n)) |
| 62 | + for i in range(n - 2): |
| 63 | + a[i, i + 1] = a[i + 1, i] = 1.0 |
| 64 | + a[1, n - 1] = a[n - 1, 1] = 1.0 |
| 65 | + return a |
| 66 | + |
| 67 | + |
| 68 | +def build_cycle(m): |
| 69 | + """Affine A_{m-1} extended Dynkin diagram: cycle on m vertices.""" |
| 70 | + a = np.zeros((m, m)) |
| 71 | + for i in range(m): |
| 72 | + a[i, (i + 1) % m] = a[(i + 1) % m, i] = 1.0 |
| 73 | + return a |
| 74 | + |
| 75 | + |
| 76 | +def build_tadpole(cycle_len, tail_len): |
| 77 | + """Cycle with a path tail attached (unicyclic, non-transitive).""" |
| 78 | + m = cycle_len + tail_len |
| 79 | + a = np.zeros((m, m)) |
| 80 | + for i in range(cycle_len): |
| 81 | + a[i, (i + 1) % cycle_len] = a[(i + 1) % cycle_len, i] = 1.0 |
| 82 | + prev = 0 |
| 83 | + for t in range(tail_len): |
| 84 | + j = cycle_len + t |
| 85 | + a[prev, j] = a[j, prev] = 1.0 |
| 86 | + prev = j |
| 87 | + return a |
| 88 | + |
| 89 | + |
| 90 | +def random_unicyclic(n, rng): |
| 91 | + """Random connected unicyclic graph on n vertices (tree + one extra edge).""" |
| 92 | + a = np.zeros((n, n)) |
| 93 | + nodes = list(rng.permutation(n)) |
| 94 | + for i in range(1, n): # random tree (random attachment) |
| 95 | + j = nodes[rng.randint(0, i)] |
| 96 | + a[nodes[i], j] = a[j, nodes[i]] = 1.0 |
| 97 | + while True: # add one non-edge -> single cycle |
| 98 | + u, v = rng.randint(0, n), rng.randint(0, n) |
| 99 | + if u != v and a[u, v] == 0: |
| 100 | + a[u, v] = a[v, u] = 1.0 |
| 101 | + return a |
| 102 | + |
| 103 | + |
| 104 | +def cycle_basis_single(adjacency): |
| 105 | + """For a unicyclic graph, return the unique cycle as a vertex list.""" |
| 106 | + n = adjacency.shape[0] |
| 107 | + deg = adjacency.sum(axis=1).astype(int) |
| 108 | + a = adjacency.copy() |
| 109 | + # iteratively strip leaves |
| 110 | + changed = True |
| 111 | + alive = set(range(n)) |
| 112 | + while changed: |
| 113 | + changed = False |
| 114 | + for v in list(alive): |
| 115 | + if a[v].sum() == 1: |
| 116 | + u = int(np.argmax(a[v])) |
| 117 | + a[v, u] = a[u, v] = 0 |
| 118 | + alive.discard(v) |
| 119 | + changed = True |
| 120 | + cyc_nodes = sorted(alive) |
| 121 | + # order the cycle by walking |
| 122 | + start = cyc_nodes[0] |
| 123 | + cycle = [start] |
| 124 | + prev, cur = None, start |
| 125 | + while True: |
| 126 | + nbrs = [j for j in np.nonzero(a[cur])[0] if j != prev] |
| 127 | + nxt = int(nbrs[0]) |
| 128 | + if nxt == start: |
| 129 | + break |
| 130 | + cycle.append(nxt) |
| 131 | + prev, cur = cur, nxt |
| 132 | + return cycle |
| 133 | + |
| 134 | + |
| 135 | +# ============================================================ |
| 136 | +# The complement-eigenvector connection |
| 137 | +# ============================================================ |
| 138 | + |
| 139 | +def complement_frame(adjacency, vertex): |
| 140 | + """Eigen-decomposition of the complement subgraph G \\ vertex. |
| 141 | +
|
| 142 | + Returns (eigvals ascending, eigvecs columns, kept_vertices list).""" |
| 143 | + n = adjacency.shape[0] |
| 144 | + keep = [i for i in range(n) if i != vertex] |
| 145 | + sub = adjacency[np.ix_(keep, keep)] |
| 146 | + vals, vecs = np.linalg.eigh(sub) |
| 147 | + return vals, vecs, keep |
| 148 | + |
| 149 | + |
| 150 | +def edge_transport(adjacency, u, v, k, frames=None): |
| 151 | + """Orthogonal transport (Procrustes) from u's complement frame to v's, |
| 152 | + over the shared support V \\ {u, v}, using the top-k eigenvectors |
| 153 | + (largest eigenvalues). Returns (T [k x k orthogonal], min_eigengap).""" |
| 154 | + if frames is None: |
| 155 | + frames = {} |
| 156 | + for w in (u, v): |
| 157 | + if w not in frames: |
| 158 | + frames[w] = complement_frame(adjacency, w) |
| 159 | + vals_u, vecs_u, keep_u = frames[u] |
| 160 | + vals_v, vecs_v, keep_v = frames[v] |
| 161 | + common = [w for w in keep_u if w != v] # = V \ {u, v} |
| 162 | + rows_u = [keep_u.index(w) for w in common] |
| 163 | + rows_v = [keep_v.index(w) for w in common] |
| 164 | + Vu = vecs_u[rows_u, :][:, -k:] # top-k by eigenvalue |
| 165 | + Vv = vecs_v[rows_v, :][:, -k:] |
| 166 | + gap_u = float(vals_u[-k] - vals_u[-k - 1]) if len(vals_u) > k else np.inf |
| 167 | + gap_v = float(vals_v[-k] - vals_v[-k - 1]) if len(vals_v) > k else np.inf |
| 168 | + M = Vv.T @ Vu |
| 169 | + U, _, Wt = np.linalg.svd(M) |
| 170 | + T = U @ Wt # orthogonal k x k |
| 171 | + return T, min(gap_u, gap_v) |
| 172 | + |
| 173 | + |
| 174 | +def cycle_holonomy(adjacency, cycle, k): |
| 175 | + """Holonomy of the connection around an ordered vertex cycle. |
| 176 | +
|
| 177 | + Returns dict: deficit ||H - I||_F, sorted |rotation angles| (conjugation |
| 178 | + invariants), min eigengap encountered (degeneracy guard).""" |
| 179 | + frames = {} |
| 180 | + H = np.eye(k) |
| 181 | + min_gap = np.inf |
| 182 | + m = len(cycle) |
| 183 | + for i in range(m): |
| 184 | + u, v = cycle[i], cycle[(i + 1) % m] |
| 185 | + T, gap = edge_transport(adjacency, u, v, k, frames) |
| 186 | + min_gap = min(min_gap, gap) |
| 187 | + H = T @ H |
| 188 | + eig = np.linalg.eigvals(H) |
| 189 | + angles = np.sort(np.abs(np.angle(eig))) |
| 190 | + deficit = float(np.linalg.norm(H - np.eye(k))) |
| 191 | + return {'deficit': deficit, |
| 192 | + 'angles': [float(a) for a in angles], |
| 193 | + 'det': float(np.linalg.det(H)), |
| 194 | + 'min_eigengap': float(min_gap)} |
| 195 | + |
| 196 | + |
| 197 | +def relabeled(adjacency, perm): |
| 198 | + """Apply vertex permutation: perm[i] = new label of old vertex i.""" |
| 199 | + n = adjacency.shape[0] |
| 200 | + P = np.zeros((n, n)) |
| 201 | + for i in range(n): |
| 202 | + P[perm[i], i] = 1.0 |
| 203 | + return P @ adjacency @ P.T |
| 204 | + |
| 205 | + |
| 206 | +# ============================================================ |
| 207 | +# Scalar potential (the exact part -- for exp_01 T1) |
| 208 | +# ============================================================ |
| 209 | + |
| 210 | +def spectral_potential(adjacency, vertex): |
| 211 | + """g(v) = sum of complement spectrum -- a vertex potential. Signed edge |
| 212 | + differences of g are exact by construction (telescoping).""" |
| 213 | + return float(np.sum(complement_spectrum(adjacency, vertex))) |
0 commit comments