From b80dd57a9b03e5494a1a879546d7af4a2a2596f7 Mon Sep 17 00:00:00 2001 From: Shay Date: Tue, 12 May 2026 19:12:48 -0700 Subject: [PATCH] =?UTF-8?q?init:=20algebra=20layer=20=E2=80=94=20cl41,=20v?= =?UTF-8?q?ersor,=20cga,=20holonomy?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- algebra/__init__.py | 4 + algebra/cga.py | 70 ++++++++++++++++++ algebra/cl41.py | 174 ++++++++++++++++++++++++++++++++++++++++++++ algebra/holonomy.py | 71 ++++++++++++++++++ algebra/versor.py | 50 +++++++++++++ 5 files changed, 369 insertions(+) create mode 100644 algebra/__init__.py create mode 100644 algebra/cga.py create mode 100644 algebra/cl41.py create mode 100644 algebra/holonomy.py create mode 100644 algebra/versor.py diff --git a/algebra/__init__.py b/algebra/__init__.py new file mode 100644 index 00000000..94df920d --- /dev/null +++ b/algebra/__init__.py @@ -0,0 +1,4 @@ +from .cl41 import geometric_product, reverse, grade_project, scalar_part, norm_squared, basis_vector +from .versor import versor_apply, normalize_to_versor, versor_condition +from .cga import cga_inner, outer_product, is_null, null_project, embed_point +from .holonomy import holonomy_encode, holonomy_similarity diff --git a/algebra/cga.py b/algebra/cga.py new file mode 100644 index 00000000..e142572b --- /dev/null +++ b/algebra/cga.py @@ -0,0 +1,70 @@ +""" +Conformal Geometric Algebra geometry on Cl(4,1). + +Key identity: for null vectors X, Y on the horosphere, + cga_inner(X, Y) = -d(X, Y)^2 / 2 +where d is Euclidean distance. + +This is the ONLY distance metric in CORE-AI. +No cosine similarity. No L2 norm. No approximate indexing. +""" + +import numpy as np +from .cl41 import geometric_product, reverse, scalar_part, basis_vector, N_COMPONENTS + + +def cga_inner(X: np.ndarray, Y: np.ndarray) -> float: + """ + Symmetric inner product: 0.5 * scalar_part(X*Y + Y*X). + For null vectors representing conformal points: equals -d^2 / 2. + """ + XY = geometric_product(X, Y) + YX = geometric_product(Y, X) + return 0.5 * scalar_part(XY + YX) + + +def outer_product(X: np.ndarray, Y: np.ndarray) -> np.ndarray: + """ + Outer (wedge) product: X ^ Y. + For a prompt versor X_p and response versor X_r, + X_p ^ X_r is a grade-2 object encoding their geometric relationship. + A real (non-imaginary) result means the dialogue is coherent. + """ + XY = geometric_product(X, Y) + YX = geometric_product(Y, X) + return 0.5 * (XY - YX) + + +def is_null(X: np.ndarray, tol: float = 1e-6) -> bool: + """Check if X lies on the null cone: X*X = 0.""" + return abs(cga_inner(X, X)) < tol + + +def null_project(X: np.ndarray) -> np.ndarray: + """ + Re-project X onto the null cone. + Call on vault entries periodically to correct floating-point null-cone drift. + This is numerical maintenance, not a heat shield. + Method: extract Euclidean part, re-embed via standard CGA point map. + """ + euclidean = X[1:4].copy().astype(np.float32) + x_sq = float(np.dot(euclidean, euclidean)) + result = np.zeros(N_COMPONENTS, dtype=np.float32) + result[1:4] = euclidean + result[4] = 0.5 * x_sq # e+ coefficient + result[5] = 1.0 # e- coefficient + return result + + +def embed_point(x: np.ndarray) -> np.ndarray: + """ + Embed a Euclidean point x in R^3 into the CGA null cone. + Standard map: X = x + (1/2)|x|^2 * e+ + e- + """ + x = np.asarray(x, dtype=np.float32) + assert len(x) == 3, "embed_point expects a 3D vector" + result = np.zeros(N_COMPONENTS, dtype=np.float32) + result[1:4] = x + result[4] = 0.5 * float(np.dot(x, x)) + result[5] = 1.0 + return result diff --git a/algebra/cl41.py b/algebra/cl41.py new file mode 100644 index 00000000..9d6d6f48 --- /dev/null +++ b/algebra/cl41.py @@ -0,0 +1,174 @@ +""" +Cl(4,1) multivector arithmetic. + +Signature: (+,+,+,+,-). Basis e1..e5. +Multivectors are float32 arrays of shape (32,) ordered by grade: + grade-0: index 0 (1 component) + grade-1: indices 1-5 (5 components) + grade-2: indices 6-15 (10 components) + grade-3: indices 16-25 (10 components) + grade-4: indices 26-30 (5 components) + grade-5: index 31 (1 component) +""" + +from itertools import combinations +from math import comb +import numpy as np + +N_DIMS = 5 +N_COMPONENTS = 32 +SIGNATURE = np.array([1, 1, 1, 1, -1], dtype=np.float64) + +# --- Grade offset table --- + +def _grade_offsets(): + offsets = [] + start = 0 + for k in range(N_DIMS + 1): + count = comb(N_DIMS, k) + offsets.append((start, count)) + start += count + return offsets + +_GRADE_OFFSETS = _grade_offsets() + +def grade_start(k: int) -> int: + return _GRADE_OFFSETS[k][0] + +def grade_count(k: int) -> int: + return _GRADE_OFFSETS[k][1] + +# --- Blade index maps --- + +def _all_blades(): + """Return ordered list of blade tuples (one per component, ordered by grade).""" + blades = [] + for k in range(N_DIMS + 1): + for combo in combinations(range(N_DIMS), k): + blades.append(combo) + return blades + +_BLADES = _all_blades() # index -> tuple of basis vector indices +_BLADE_TO_IDX = {b: i for i, b in enumerate(_BLADES)} + +def _blade_product(blade_a, blade_b): + """Compute geometric product of two basis blades. Returns (sign, result_blade_tuple).""" + # Concatenate and bubble-sort, tracking sign flips and metric contractions + seq = list(blade_a) + list(blade_b) + sign = 1 + # Bubble sort to canonical order, tracking swaps + n = len(seq) + for i in range(n): + for j in range(n - i - 1): + if seq[j] > seq[j + 1]: + seq[j], seq[j + 1] = seq[j + 1], seq[j] + sign *= -1 + elif seq[j] == seq[j + 1]: + # Metric contraction: e_i^2 = signature[i] + metric = int(SIGNATURE[seq[j]]) + sign *= metric + seq.pop(j) + seq.pop(j) # second element now at same index after first pop + n -= 2 + break + else: + continue + break + # After contraction there may still be duplicates — recurse + result = tuple(sorted(set(seq))) # this is wrong for multi-contraction; use proper loop + return sign, tuple(seq) + +def _build_multiplication_table(): + """Precompute full 32x32 geometric product table.""" + table_idx = np.zeros((N_COMPONENTS, N_COMPONENTS), dtype=np.int32) + table_sign = np.zeros((N_COMPONENTS, N_COMPONENTS), dtype=np.float32) + + for i, blade_a in enumerate(_BLADES): + for j, blade_b in enumerate(_BLADES): + sign, result_blade = _compute_blade_product(blade_a, blade_b) + result_idx = _BLADE_TO_IDX.get(result_blade, 0) + table_idx[i, j] = result_idx + table_sign[i, j] = sign + + return table_idx, table_sign + +def _compute_blade_product(blade_a, blade_b): + """Compute geometric product of two basis blades via bubble sort + metric.""" + seq = list(blade_a) + list(blade_b) + sign = 1 + i = 0 + while i < len(seq) - 1: + j = i + while j < len(seq) - 1: + if seq[j] == seq[j + 1]: + # Contract: e_k^2 = signature[k] + sign *= int(SIGNATURE[seq[j]]) + seq.pop(j) + seq.pop(j) + if j > 0: + i = max(0, j - 1) + break + elif seq[j] > seq[j + 1]: + seq[j], seq[j + 1] = seq[j + 1], seq[j] + sign *= -1 + j += 1 + else: + j += 1 + else: + i += 1 + result_blade = tuple(seq) + if result_blade not in _BLADE_TO_IDX: + return 0, () + return sign, result_blade + +_TABLE_IDX, _TABLE_SIGN = _build_multiplication_table() + +# --- Core operations --- + +def geometric_product(A: np.ndarray, B: np.ndarray) -> np.ndarray: + """Full geometric product in Cl(4,1).""" + A = np.asarray(A, dtype=np.float32) + B = np.asarray(B, dtype=np.float32) + result = np.zeros(N_COMPONENTS, dtype=np.float32) + for i in range(N_COMPONENTS): + if A[i] == 0.0: + continue + for j in range(N_COMPONENTS): + if B[j] == 0.0: + continue + result[_TABLE_IDX[i, j]] += _TABLE_SIGN[i, j] * A[i] * B[j] + return result + +def reverse(A: np.ndarray) -> np.ndarray: + """ + Reverse (main anti-automorphism). + Grade-k blades pick up sign (-1)^(k*(k-1)/2). + Grade 0,1: +1. Grade 2,3: -1. Grade 4,5: +1. + """ + A = np.asarray(A, dtype=np.float32).copy() + # Grade 2: indices 6-15 + A[6:16] *= -1.0 + # Grade 3: indices 16-25 + A[16:26] *= -1.0 + return A + +def grade_project(A: np.ndarray, k: int) -> np.ndarray: + """Extract grade-k part of A.""" + result = np.zeros(N_COMPONENTS, dtype=np.float32) + start, count = _GRADE_OFFSETS[k] + result[start:start + count] = A[start:start + count] + return result + +def scalar_part(A: np.ndarray) -> float: + """Return grade-0 component.""" + return float(A[0]) + +def norm_squared(A: np.ndarray) -> float: + """||A||^2 = scalar_part(A * reverse(A)).""" + return scalar_part(geometric_product(A, reverse(A))) + +def basis_vector(i: int) -> np.ndarray: + """Return the i-th basis vector (0-indexed) as a 32-component multivector.""" + v = np.zeros(N_COMPONENTS, dtype=np.float32) + v[1 + i] = 1.0 + return v diff --git a/algebra/holonomy.py b/algebra/holonomy.py new file mode 100644 index 00000000..c4edb2a7 --- /dev/null +++ b/algebra/holonomy.py @@ -0,0 +1,71 @@ +""" +Holonomy prompt encoding. + +A prompt w1, w2, ..., wn is encoded as the geometric holonomy of its +forward+reverse versor walk. The walk closes, producing a versor that +is bounded by construction and invariant to global phase. + +The holonomy IS a versor — it drops directly into versor_apply with +no bridging code. The fuel and the engine are the same substance. +""" + +import numpy as np +from .cl41 import geometric_product, reverse as cl_reverse +from .versor import normalize_to_versor +from .cga import cga_inner + + +def holonomy_encode( + word_versors: list, + alpha: float = 0.5, + weights: list = None, +) -> np.ndarray: + """ + Compute the holonomy of the word versor sequence. + + Forward walk: F = w1 * w2 * ... * wn (weighted by word frequency inverse) + Reverse walk: R = (1-alpha) * reverse(wn) * ... * reverse(w1) + Holonomy: H = geometric_product(F, R) + + H is a versor. For alpha=0.5, the holonomy captures the geometric + curvature of the prompt path. Prompts with different semantic content + produce geometrically distinct holonomies even at the same length. + + weights: optional list of float scalars (e.g. inverse token frequency). + Rare content words rotate more than common function words. + If None, uniform weights are used. + """ + if not word_versors: + raise ValueError("Cannot encode empty prompt.") + + n = len(word_versors) + if weights is None: + weights = [1.0] * n + assert len(weights) == n + + # Forward accumulation + F = word_versors[0].copy() * weights[0] + F = normalize_to_versor(F) + for k in range(1, n): + w = word_versors[k] * weights[k] + w = normalize_to_versor(w) + F = geometric_product(F, w) + + # Reverse accumulation with alpha damping + R = cl_reverse(word_versors[-1]) * (1.0 - alpha) + R = normalize_to_versor(R) + for k in range(n - 2, -1, -1): + r = cl_reverse(word_versors[k]) + r = normalize_to_versor(r) + R = geometric_product(r, R) + + H = geometric_product(F, R) + return normalize_to_versor(H) + + +def holonomy_similarity(H1: np.ndarray, H2: np.ndarray) -> float: + """ + Compare two holonomies via CGA inner product. + Used for prompt-level semantic similarity without embedding lookup. + """ + return cga_inner(H1, H2) diff --git a/algebra/versor.py b/algebra/versor.py new file mode 100644 index 00000000..f644c152 --- /dev/null +++ b/algebra/versor.py @@ -0,0 +1,50 @@ +""" +The three versor primitives. + +These are the ONLY normalization/transition/check functions in the system. +Do not add correction, monitoring, or grade-guard functions here. +If you think you need something else, you have an unclosed operation upstream. +""" + +import numpy as np +from .cl41 import geometric_product, reverse, scalar_part, norm_squared + + +def versor_apply(V: np.ndarray, F: np.ndarray) -> np.ndarray: + """ + Sandwich product: V * F * reverse(V). + + The ONLY allowed field transition in the system. + Algebraically closed on the versor manifold: + if V and F are versors, V*F*reverse(V) is a versor. + No pre/post normalization. No grade projection. No guards. + """ + return geometric_product(V, geometric_product(F, reverse(V))) + + +def normalize_to_versor(F: np.ndarray) -> np.ndarray: + """ + Project F onto the versor manifold: F / sqrt(|F * reverse(F)|). + + Call this ONCE per input at the injection gate (ingest/gate.py). + Never call mid-propagation, mid-generation, or in the vault. + If you feel the urge to call this elsewhere, fix the upstream operation. + """ + n2 = norm_squared(F) + if abs(n2) < 1e-12: + raise ValueError("Cannot normalize a null multivector to a versor.") + return F / np.sqrt(abs(n2)) + + +def versor_condition(F: np.ndarray) -> float: + """ + Returns ||F * reverse(F) - 1||_F. + + Zero means F is on the versor manifold. + Use in tests and at the injection gate only. + Never call in the generation hot path. + """ + product = geometric_product(F, reverse(F)) + product = product.copy() + product[0] -= 1.0 + return float(np.linalg.norm(product))