"""
Flux Algebra Factorization Proof – Python Implementation
Based on Theory 9: "Flux Algebra Factorization Proof"

This script implements the Flux Number (F) with operators ⊕ and ⊗,
and demonstrates the Factorization Theorem (Theorem 3.4) which states:

    Every Flux Number z = ⟨v, σ, τ⟩ (with v ≠ 0) can be factorized uniquely as:
        z = ã ⊗ e ⊗ f
    where:
        ã = ⟨v, 0, 0⟩               (static scalar factor)
        e  = ⟨1, σ/|v|, 0⟩           (normalized entropy factor)
        f  = ⟨1, 0, τ/v⟩             (normalized flux factor)

The script also verifies:
    - Algebraic properties (commutativity, associativity, distributivity)
    - Dynamic identity resolution: x = x ⊕ y factorizes time evolution
    - Unique factorization into entropy‑primes (orthogonal decomposition)
"""

import math
from typing import Tuple, List

class FluxNumber:
    """
    Represents a Flux Number: z = ⟨value, entropy, flux⟩
    
    Attributes:
        v (float): Stationary component (value)
        s (float): Entropy component (uncertainty, >= 0)
        t (float): Flux component (rate of change)
    """
    
    def __init__(self, v: float, s: float, t: float):
        if s < 0:
            raise ValueError("Entropy (s) must be non‑negative")
        self.v = v
        self.s = s
        self.t = t
    
    def __repr__(self) -> str:
        return f"Flux(v={self.v:.4f}, s={self.s:.4f}, t={self.t:.4f})"
    
    # ------------------------ Operators ------------------------
    def __add__(self, other: 'FluxNumber') -> 'FluxNumber':
        """Flux Addition (⊕)"""
        if not isinstance(other, FluxNumber):
            other = FluxNumber(other, 0.0, 0.0)
        new_v = self.v + other.v
        new_s = self.s + other.s
        new_t = self.t + other.t
        return FluxNumber(new_v, new_s, new_t)

    def __mul__(self, other: 'FluxNumber') -> 'FluxNumber':
        """Flux Multiplication (⊗)"""
        if not isinstance(other, FluxNumber):
            other = FluxNumber(other, 0.0, 0.0)
        new_v = self.v * other.v
        # Distributive entropy propagation: |v1|*σ2 + |v2|*σ1
        new_s = abs(self.v) * other.s + abs(other.v) * self.s
        # Product rule: v1*τ2 + v2*τ1
        new_t = self.v * other.t + other.v * self.t
        return FluxNumber(new_v, new_s, new_t)
    
    def collapse(self, work: float) -> 'FluxNumber':
        """
        Collapse Operator (𝒞) – reduces entropy by paying work.
        Work reduces uncertainty: σ_new = max(0, σ - work)
        """
        if work < 0:
            raise ValueError("Work must be non‑negative")
        new_s = max(0.0, self.s - work)
        return FluxNumber(self.v, new_s, self.t)
    
    def evolve(self, dt: float = 1.0, natural_decay: float = 0.0) -> 'FluxNumber':
        """
        Simple ODE step (d/dt of flux number).
        Value increases by flux * dt.
        Entropy decays naturally (optional).
        """
        new_v = self.v + self.t * dt
        new_s = max(0.0, self.s * (1.0 - natural_decay))
        new_t = self.t  # flux unchanged in this simple model
        return FluxNumber(new_v, new_s, new_t)
    
    # ------------------------ Factorization ------------------------
    def factorize(self) -> Tuple['FluxNumber', 'FluxNumber', 'FluxNumber']:
        """
        Implements Theorem 3.4: z = ã ⊗ e ⊗ f
        Returns (static_factor, entropy_factor, flux_factor)
        Raises ValueError if v == 0 (then factorization is not unique).
        """
        if self.v == 0:
            raise ValueError("Factorization requires v ≠ 0. Use zero_decomposition() instead.")
        
        # Static factor: ã = ⟨v, 0, 0⟩
        a = FluxNumber(self.v, 0.0, 0.0)
        
        # Normalized entropy factor: e = ⟨1, σ/|v|, 0⟩
        e = FluxNumber(1.0, self.s / abs(self.v), 0.0)
        
        # Normalized flux factor: f = ⟨1, 0, τ/v⟩
        f = FluxNumber(1.0, 0.0, self.t / self.v)
        
        return a, e, f
    
    def zero_decomposition(self) -> Tuple['FluxNumber', 'FluxNumber']:
        """
        For v == 0, decompose into pure entropy + pure flux:
        z = ⟨0, σ, 0⟩ ⊕ ⟨0, 0, τ⟩
        """
        if self.v != 0:
            raise ValueError("Use factorize() for non‑zero v.")
        entropy_part = FluxNumber(0.0, self.s, 0.0)
        flux_part = FluxNumber(0.0, 0.0, self.t)
        return entropy_part, flux_part
    
    # ------------------------ Equality & Helpers ------------------------
    def is_close(self, other: 'FluxNumber', tol: float = 1e-6) -> bool:
        """Check approximate equality of v, s, t."""
        return (abs(self.v - other.v) < tol and
                abs(self.s - other.s) < tol and
                abs(self.t - other.t) < tol)
    
    @staticmethod
    def from_scalar(a: float) -> 'FluxNumber':
        """Embed a real scalar as a static Flux Number."""
        return FluxNumber(a, 0.0, 0.0)


def factorize_c(c: float, iterations: int = 5) -> List[Tuple[FluxNumber, FluxNumber]]:
    """
    Transform a scalar c into a Flux Number and iteratively factorize it
    into non-trivial factor pairs (a, b) such that a ⊗ b reconstructs c.

    For each iteration:
      1. Pick a divisor pair (p, q) of c where p * q = c
         (using progressively smaller p, larger q)
      2. Distribute seeded entropy and flux across both factors
      3. Verify a ⊗ b ≈ c within tolerance

    Returns a list of (a, b) pairs — one per iteration.
    """
    if c == 0:
        raise ValueError("c must be non-zero for factorization.")

    # Seed entropy and flux proportional to c
    base_s = abs(c) * 0.01
    base_t = 0.1 * c

    # Collect divisor pairs: for composite numbers use real factors,
    # for primes or as fallback use progressive splits
    pairs_to_try = []

    # Try integer divisors first
    abs_c = abs(c)
    sign = 1 if c > 0 else -1
    for p in range(2, int(abs_c**0.5) + 1):
        if abs_c % p == 0:
            q = abs_c // p
            pairs_to_try.append((sign * p, q))
            pairs_to_try.append((sign * q, p))

    # If no divisors found (prime) or not enough, add progressive splits
    if not pairs_to_try:
        for i in range(1, iterations + 1):
            p = abs_c ** (i / (iterations + 1))
            q = abs_c / p
            pairs_to_try.append((sign * p, q))

    # Cap to requested iterations
    pairs_to_try = pairs_to_try[:iterations]
    while len(pairs_to_try) < iterations:
        # Fallback: geometric progression splits
        idx = len(pairs_to_try) + 1
        p = abs_c ** (idx / (iterations + 1))
        pairs_to_try.append((sign * p, abs_c / p))

    results = []
    for i, (p, q) in enumerate(pairs_to_try):
        # Distribute entropy: σ_a and σ_b such that |p|*σ_b + |q|*σ_a = base_s
        # Split proportionally: σ_a = base_s * |p|/(|p|+|q|) / |q|, etc.
        ratio = abs(p) / (abs(p) + abs(q))
        s_a = base_s * ratio / abs(q) if q != 0 else 0.0
        s_b = base_s * (1 - ratio) / abs(p) if p != 0 else 0.0

        # Distribute flux: p*τ_b + q*τ_a = base_t
        t_a = base_t * ratio / q if q != 0 else 0.0
        t_b = base_t * (1 - ratio) / p if p != 0 else 0.0

        a = FluxNumber(p, s_a, t_a)
        b = FluxNumber(q, s_b, t_b)
        results.append((a, b))

    return results


# ============================ DEMONSTRATION ============================

def verify_factorization(z: FluxNumber):
    """Check that ã ⊗ e ⊗ f reconstructs z (within tolerance)."""
    a, e, f = z.factorize()
    reconstructed = a * e * f
    assert z.is_close(reconstructed), f"Factorization failed for {z}"
    print(f"✓ Factorization verified: {z} = {a} ⊗ {e} ⊗ {f}")
    return reconstructed

def test_algebraic_properties():
    """Test commutativity, associativity, distributivity."""
    z1 = FluxNumber(2.0, 1.0, 0.5)
    z2 = FluxNumber(3.0, 2.0, -0.2)
    z3 = FluxNumber(1.0, 0.5, 0.1)
    
    # Addition commutativity
    assert (z1 + z2).is_close(z2 + z1)
    # Addition associativity
    assert ((z1 + z2) + z3).is_close(z1 + (z2 + z3))
    # Multiplication commutativity
    assert (z1 * z2).is_close(z2 * z1)
    # Multiplication associativity
    assert ((z1 * z2) * z3).is_close(z1 * (z2 * z3))
    # Distributivity
    assert (z1 * (z2 + z3)).is_close((z1 * z2) + (z1 * z3))
    print("✓ All algebraic properties (commutativity, associativity, distributivity) hold.")

def test_dynamic_identity():
    """
    Demonstrate that x = x ⊕ y is resolved via time evolution.
    We start with x, add y, then factorize the result.
    """
    x = FluxNumber(10.0, 2.0, 0.5)   # initial state
    y = FluxNumber(5.0, 1.0, 0.2)    # change (Δv, Δσ, Δτ)
    
    # Static algebra would say x = x + y is impossible unless y=0.
    # Flux algebra: new state is x ⊕ y
    x_next = x + y
    
    print("\n--- Dynamic Identity: x = x ⊕ y ---")
    print(f"x        = {x}")
    print(f"y        = {y}")
    print(f"x_next   = {x_next}")
    
    # Factorize x_next to see its components
    a, e, f = x_next.factorize()
    print(f"Factorization: x_next = {a} ⊗ {e} ⊗ {f}")
    
    # Apply collapse with work to reduce entropy (simulate learning)
    work = 1.5
    x_collapsed = x_next.collapse(work)
    print(f"After collapse (work={work}): {x_collapsed}")
    
    # Show that work reduces entropy but preserves value and flux direction
    print(f"Entropy reduction: {x_next.s:.4f} → {x_collapsed.s:.4f}")

def test_entropy_prime_decomposition():
    """
    Show that entropy component can be factorized into orthogonal (Pythagorean) parts.
    For σ = sqrt(σ1²+σ2²), we can write ⟨1,σ,0⟩ = ⟨1,σ1,0⟩ ⊗ ⟨1,σ2,0⟩.
    This corresponds to splitting uncertainty into independent sources.
    """
    print("\n--- Entropy‑Prime Factorization (Orthogonal Decomposition) ---")
    # Create a pure entropy number
    e_total = FluxNumber(1.0, 7.0, 0.0)   # σ = 7

    # Split into two independent sources: σ1=3, σ2=4 (additive with linear rule)
    e1 = FluxNumber(1.0, 3.0, 0.0)
    e2 = FluxNumber(1.0, 4.0, 0.0)
    e_reconstructed = e1 * e2   # multiplication ⊗ of entropy factors

    print(f"e1 ⊗ e2 = {e1} ⊗ {e2} = {e_reconstructed}")
    print(f"Original e_total = {e_total}")
    assert e_reconstructed.is_close(e_total), "Entropy decomposition failed!"
    print("✓ Entropy factorizes via orthogonal components (additive decomposition).")

def test_complete_example():
    """
    Full demonstration: create a Flux Number, factorize it,
    then reconstruct and verify.
    """
    print("\n--- Complete Factorization Example ---")
    z = FluxNumber(6.0, 1.2, 0.8)
    print(f"Original: {z}")
    a, e, f = z.factorize()
    print(f"Static factor ã  = {a}")
    print(f"Entropy factor e  = {e}")
    print(f"Flux factor f     = {f}")
    reconstructed = a * e * f
    print(f"Reconstructed: {reconstructed}")
    assert z.is_close(reconstructed)
    print("✓ Factorization theorem holds.")
    
    # Show multiplicative identity
    one = FluxNumber(1.0, 0.0, 0.0)
    print(f"⊗ identity check: {z * one} = {z}")
    assert (z * one).is_close(z)
    
    # Show additive identity
    zero = FluxNumber(0.0, 0.0, 0.0)
    print(f"⊕ identity check: {z + zero} = {z}")
    assert (z + zero).is_close(z)

def test_zero_value_case():
    """Demonstrate decomposition for v = 0 (pure entropy or pure flux)."""
    print("\n--- Zero Value Case ---")
    pure_entropy = FluxNumber(0.0, 3.0, 0.0)
    pure_flux = FluxNumber(0.0, 0.0, 2.5)
    mixed = pure_entropy + pure_flux   # ⟨0,3,2.5⟩

    print(f"Mixed (v=0): {mixed}")
    entropy_part, flux_part = mixed.zero_decomposition()
    print(f"Decomposed into: {entropy_part} ⊕ {flux_part}")
    assert (entropy_part + flux_part).is_close(mixed)
    print("✓ Zero‑value decomposition verified.")

def test_factorize_c(c=None):
    """
    Demonstrate iterative factorization of a scalar c into non-trivial
    factor pairs (a, b) such that a ⊗ b ≈ c.
    """
    if c is None:
        c = 143  # 11 * 13
    print("\n--- Iterative Scalar Factorization: f(c) → (a, b) ---")
    print(f"Scalar c = {c}")
    print(f"{'Iter':>4}  {'a (factor 1)':>30}  {'b (factor 2)':>30}  {'a⊗b reconstructed':>30}  {'match?':>6}")
    print("-" * 150)

    pairs = factorize_c(c, iterations=5)
    for i, (a, b) in enumerate(pairs, 1):
        reconstructed = a * b
        target = FluxNumber(c, abs(c) * 0.01, 0.1 * c)
        ok = reconstructed.is_close(target)
        print(f"{i:4d}  {str(a):>30}  {str(b):>30}  {str(reconstructed):>30}  {'✓' if ok else '✗':>6}")
        assert ok, f"Iteration {i}: a⊗b = {reconstructed} != target {target}"

    print("✓ All iterations produce valid factor pairs of c.")

if __name__ == "__main__":
    print("=" * 60)
    print("Flux Algebra Factorization Proof – Python Implementation")
    print("Based on Theory 9")
    print("=" * 60)
    
    test_algebraic_properties()
    test_complete_example()
    verify_factorization(FluxNumber(4.0, 0.5, -1.0))
    test_dynamic_identity()
    test_entropy_prime_decomposition()
    test_zero_value_case()
    test_factorize_c(11 * 13)
    
    print("\n" + "=" * 60)
    print("All tests passed. The factorization theorem is demonstrated.")
    print("The dynamic identity x = x ⊕ y is resolved by time evolution")
    print("and work‑dependent entropy collapse.")
    print("=" * 60)