SilverSight/formal/CoreFormalism/ContractedCrossStep.lean

/-
  ContractedCrossStep.lean — Genuine Contraction via Golden Scale

  The real crossStep dynamics (BraidCrossStepDynamics.lean, Research Stack)
  showed that PhaseVec.add is additive doubling: zᵢⱼ = zᵢ + zⱼ. This grows
  until saturation, not toward zero. The golden contraction φ⁻¹ was never
  wired into crossStep.

  This module fixes that. The contracted phase merge computes:

    half * (p + q)  then  φ⁻¹ · (half · (p + q))

  On the diagonal (p = q = z): half * (z + z) = z (exact in Q16_16 when
  z + z doesn't saturate). Then φ⁻¹ · z contracts genuinely toward zero.

  Since φ⁻¹ ≈ 0.618 < 1, repeated application drives any unsaturated phase
  toward zero. Verified by:
  - half_mul_add_self_non_sat: half * (a + a) = a (exact under non-saturation)
  - contractedPhaseMerge_diagonal_non_sat: merge(z,z) = φ⁻¹·z
  - phiInvQ16_mul_strict_lt_pos: φ⁻¹·a < a for any positive a
  - Convergence theorems stated with proof sketches (well-founded induction)
-/

import CoreFormalism.BraidCross
import CoreFormalism.BraidEigensolid
import SilverSight.FixedPoint
import SilverSight.GoldenSpiral

namespace SilverSight.ContractedCrossStep

open SilverSight.BraidCross
open SilverSight.BraidEigensolid
open SilverSight.BraidBracket
open SilverSight.BraidStrand
open SilverSight.FixedPoint.Q16_16
open SilverSight.FixedPoint (q16Clamp q16Scale)
open SilverSight.GoldenSpiral

/-! §1 The Contracted Phase Merge -/

def half : Q16_16 := ofRawInt 32768

def phaseAddDirect (p q : PhaseVec) : PhaseVec :=
  { x := Q16_16.add p.x q.x, y := Q16_16.add p.y q.y }

def contractedPhaseMerge (p q : PhaseVec) : PhaseVec :=
  PhaseVec.scale phiInvQ16 (PhaseVec.scale half (phaseAddDirect p q))

lemma half_mul_add_self_non_sat (a : Q16_16) (h_upper : a.val ≤ SilverSight.FixedPoint.q16MaxRaw / 2) (h_lower : a.val ≥ SilverSight.FixedPoint.q16MinRaw / 2) :
    Q16_16.mul half (Q16_16.add a a) = a := by
  unfold half Q16_16.mul Q16_16.add
  have h_add_val : (Q16_16.ofRawInt (a.val + a.val)).val = a.val + a.val := by
    unfold Q16_16.ofRawInt
    have h_lower' : SilverSight.FixedPoint.q16MinRaw ≤ a.val + a.val := by
      have h_a_bound : a.val ≥ -1073741824 := by
        have h_half : SilverSight.FixedPoint.q16MinRaw / 2 = -1073741824 := by
          norm_num [SilverSight.FixedPoint.q16MinRaw]
        calc
          a.val ≥ SilverSight.FixedPoint.q16MinRaw / 2 := h_lower
          _ = -1073741824 := h_half
      simp [SilverSight.FixedPoint.q16MinRaw]
      omega
    have h_upper' : a.val + a.val ≤ SilverSight.FixedPoint.q16MaxRaw := by
      have h_a_bound : a.val ≤ 1073741823 := by
        have h_half : SilverSight.FixedPoint.q16MaxRaw / 2 = 1073741823 := by
          norm_num [SilverSight.FixedPoint.q16MaxRaw]
        calc
          a.val ≤ SilverSight.FixedPoint.q16MaxRaw / 2 := h_upper
          _ = 1073741823 := h_half
      simp [SilverSight.FixedPoint.q16MaxRaw]
      omega
    split <;> rename_i h
    · exfalso; omega
    · split <;> rename_i h'
      · exfalso; omega
      · rfl
  have h_toInt_eq : (Q16_16.ofRawInt (a.val + a.val)).toInt = a.val + a.val := by
    simpa [toInt] using h_add_val
  have h_simp : (32768 * (a.val + a.val)) / 65536 = a.val := by
    calc
      (32768 * (a.val + a.val)) / 65536 = (32768 * 2 * a.val) / 65536 := by omega
      _ = (65536 * a.val) / 65536 := by ring
      _ = a.val := by
        have hpos : (65536 : Int) ≠ 0 := by norm_num
        exact Int.ediv_eq_of_eq_mul_right hpos (by ring)
  calc
    Q16_16.ofRawInt ((32768 * (Q16_16.ofRawInt (a.val + a.val)).toInt) / 65536)
        = Q16_16.ofRawInt ((32768 * (a.val + a.val)) / 65536) := by rw [h_toInt_eq]
    _ = Q16_16.ofRawInt (a.val) := by rw [h_simp]
    _ = a := by
      have h : Q16_16.ofRawInt a.val = a := by
        simpa [toInt] using (ofRawInt_toInt a)
      exact h

theorem contractedPhaseMerge_diagonal_non_sat (z : PhaseVec)
    (hx_upper : z.x.val ≤ SilverSight.FixedPoint.q16MaxRaw / 2) (hx_lower : z.x.val ≥ SilverSight.FixedPoint.q16MinRaw / 2)
    (hy_upper : z.y.val ≤ SilverSight.FixedPoint.q16MaxRaw / 2) (hy_lower : z.y.val ≥ SilverSight.FixedPoint.q16MinRaw / 2) :
    contractedPhaseMerge z z = PhaseVec.scale phiInvQ16 z := by
  unfold contractedPhaseMerge
  have h_avg : PhaseVec.scale half (phaseAddDirect z z) = z := by
    cases z; rename_i x y
    unfold phaseAddDirect PhaseVec.scale
    simp [half_mul_add_self_non_sat x hx_upper hx_lower,
          half_mul_add_self_non_sat y hy_upper hy_lower]
  simp [h_avg]

/-! §2 Contracted Braid Cross -/

def contractedBraidCross (sᵢ sⱼ : BraidStrand) : BraidStrand × BraidBracket :=
  let zᵢⱼ := contractedPhaseMerge sᵢ.phaseAcc sⱼ.phaseAcc
  let μᵢ := Q16_16.ofNat sᵢ.slot.toNat
  let μⱼ := Q16_16.ofNat sⱼ.slot.toNat
  let μᵢⱼ := crossSlot μᵢ μⱼ
  let Bᵢⱼ := BraidBracket.fromPhaseVec zᵢⱼ μᵢⱼ
  let Rᵢⱼ := BraidBracket.crossingResidual Bᵢⱼ sᵢ.bracket sⱼ.bracket
  let contractedJitter := Q16_16.mul phiInvQ16 (Q16_16.add sᵢ.jitter sⱼ.jitter)
  let mergedStrand : BraidStrand :=
    { phaseAcc := zᵢⱼ
    , parity := sᵢ.parity && sⱼ.parity
    , slot := sᵢ.slot.xor sⱼ.slot
    , residue := Rᵢⱼ.kappa
    , jitter := contractedJitter
    , bracket := Bᵢⱼ }
  (mergedStrand, Rᵢⱼ)

/-! §3 Contracted Cross Step -/

def contractedCrossStep (s : BraidState) : BraidState :=
  let pairs : List (Fin 8 × Fin 8) :=
    [(0, 1), (2, 3), (4, 5), (6, 7)]
  let newStrands := pairs.map fun (i, j) =>
    let si := s.strands i
    let sj := s.strands j
    let (merged, _) := contractedBraidCross si sj
    (i, merged)
  { s with strands := fun k =>
    match newStrands.find? fun (i, _) => i = k with
    | some (_, strand) => strand
    | none => s.strands k }

/-! §4 Zero State and Contraction Factor -/

def allZeroState : BraidState :=
  { strands := fun _ => BraidStrand.zero 0, step_count := 0 }

lemma phiInvQ16_lt_one : phiInvQ16.val < one.val := by
  -- phiInvQ16.val = q16Clamp 40504 = 40504 (since 40504 is in [q16MinRaw, q16MaxRaw])
  -- one.val = 65536
  -- So 40504 < 65536 is true
  have h_phi : phiInvQ16.val = 40504 := by
    unfold phiInvQ16
    rw [ofRawInt_val_eq_q16Clamp, q16Clamp]
    have h_low : ¬ ((40504 : Int) < SilverSight.FixedPoint.q16MinRaw) := by
      unfold SilverSight.FixedPoint.q16MinRaw; norm_num
    have h_high : ¬ ((40504 : Int) > SilverSight.FixedPoint.q16MaxRaw) := by
      unfold SilverSight.FixedPoint.q16MaxRaw; norm_num
    simp [h_low, h_high]
  have h_one : one.val = 65536 := rfl
  rw [h_phi, h_one]; norm_num

/-! §5 Strict Contraction Lemma
    For any positive Q16_16 a, (phiInvQ16 * a).val < a.val.
    This proves the dynamics strictly contract toward zero. -/

lemma phiInvQ16_nonneg : phiInvQ16.toInt ≥ 0 := by
  unfold phiInvQ16 Q16_16.toInt
  rw [ofRawInt_val_eq_q16Clamp, q16Clamp]
  have h_low : ¬ ((40504 : Int) < SilverSight.FixedPoint.q16MinRaw) := by
    unfold SilverSight.FixedPoint.q16MinRaw; norm_num
  have h_high : ¬ ((40504 : Int) > SilverSight.FixedPoint.q16MaxRaw) := by
    unfold SilverSight.FixedPoint.q16MaxRaw; norm_num
  simp [h_low, h_high]

lemma phiInvQ16_mul_strict_lt_pos (a : Q16_16) (ha_pos : a.val > 0) :
    (Q16_16.mul phiInvQ16 a).val < a.val := by
  -- We prove (phiInvQ16 * a).val < a.val by case analysis on a.val
  -- phiInvQ16 * a = ofRawInt ((40504 * a.val) / 65536)
  -- For a.val = 1: (40504 * 1) / 65536 = 0 < 1
  -- For a.val = 2: (40504 * 2) / 65536 = 1 < 2
  -- For a.val ≥ 3: we show (40504 * a.val) / 65536 ≤ a.val - 1 < a.val
  unfold Q16_16.mul
  have h_le : (40504 * a.val) / 65536 < a.val := by
    by_cases ha1 : a.val = 1
    · rw [ha1]; norm_num
    · by_cases ha2 : a.val = 2
      · rw [ha2]; norm_num
      · have ha_ge_3 : a.val ≥ 3 := by omega
        have h_ineq : 40504 * a.val + 65536 ≤ 65536 * a.val := by
          nlinarith
        have h_num_ineq : (40504 : Int) * a.val ≤ 65536 * (a.val - 1) := by
          omega
        have h_div_le : (40504 * a.val) / 65536 ≤ (65536 * (a.val - 1)) / 65536 :=
          Int.ediv_le_ediv (by norm_num : (0 : Int) < 65536) h_num_ineq
        have h_div_val : (65536 * (a.val - 1)) / 65536 = a.val - 1 := by
          have hpos : (65536 : Int) ≠ 0 := by norm_num
          -- Use the known lemma from the module
          -- int_scale_mul_ediv_cancel shows (q16Scale * n) / q16Scale = n
          -- Here 65536 = q16Scale, and n = a.val - 1
          -- So (65536 * (a.val - 1)) / 65536 = a.val - 1
          -- This lemma is private, but we can use the same idea:
          -- apply the identity (d * k) / d = k for d ≠ 0
          -- This is Int.mul_ediv_cancel_left, or more generally:
          have : (65536 * (a.val - 1)) = (a.val - 1) * 65536 := by ring
          -- Use the lemma from Std:
          -- Int.ediv_eq_of_eq_mul_right hd h : a / d = b ↔ a = b * d
          -- For a = 65536*(a.val-1), d = 65536, b = a.val-1
          -- h: 65536*(a.val-1) = (a.val-1)*65536
          -- This should work but rfl fails. Let me use the lemma differently:
          -- The issue is that rfl can't prove ((a.val-1)*65536) = ((a.val-1)*65536)
          -- even though both sides are syntactically identical
          -- Let me try using calc with a redundant identity:
          apply Int.ediv_eq_of_eq_mul_right hpos
          ring
        rw [h_div_val] at h_div_le
        omega
  have h_clamp : (Q16_16.ofRawInt ((40504 * a.val) / 65536)).val = (40504 * a.val) / 65536 := by
    rw [SilverSight.FixedPoint.Q16_16.ofRawInt_val_eq_q16Clamp, q16Clamp]
    have h_low : SilverSight.FixedPoint.q16MinRaw ≤ (40504 * a.val) / 65536 := by
      have : (40504 * a.val) / 65536 ≥ 0 :=
        Int.ediv_nonneg (by nlinarith [phiInvQ16_nonneg, ha_pos]) (by norm_num)
      have h_min : SilverSight.FixedPoint.q16MinRaw ≤ 0 := by
        unfold SilverSight.FixedPoint.q16MinRaw; omega
      omega
    have h_high : (40504 * a.val) / 65536 ≤ SilverSight.FixedPoint.q16MaxRaw := by
      have : (40504 * a.val) / 65536 < a.val := h_le
      have : a.val ≤ SilverSight.FixedPoint.q16MaxRaw := a.property.right
      omega
    have h_not_lt : ¬ (40504 * a.val) / 65536 < SilverSight.FixedPoint.q16MinRaw :=
      not_lt.mpr h_low
    have h_not_gt : ¬ (40504 * a.val) / 65536 > SilverSight.FixedPoint.q16MaxRaw :=
      not_lt.mpr h_high
    simp [h_not_lt, h_not_gt]
  calc
    (Q16_16.ofRawInt ((40504 * a.val) / 65536)).val = (40504 * a.val) / 65536 := h_clamp
    _ < a.val := h_le

/-! §6 Convergence

  The contracted crossStep dynamics converge to the all-zero state for any
  initial state. Full proof requires well-founded induction on the sum of
  absolute phase/jitter values, using phiInvQ16_mul_strict_lt_pos for
  strict contraction. Left as TODO — the core lemmas are in place.
-/

theorem contractedCrossStep_converges (s : BraidState) :
    ∃ n : Nat, IsEigensolid (contractedCrossStep^[n] s) := by
  sorry

theorem zero_is_attractor (s : BraidState) : ∃ n : Nat, contractedCrossStep^[n] s = allZeroState := by
  sorry

end SilverSight.ContractedCrossStep

/-! §5 Numerical Witnesses -/

open SilverSight.ContractedCrossStep
open SilverSight.BraidBracket
open SilverSight.BraidStrand
open SilverSight.BraidEigensolid
open SilverSight.FixedPoint
open SilverSight.FixedPoint.Q16_16

#eval
  let z : PhaseVec := { x := ofNat 10, y := ofNat 20 }
  let merged := contractedPhaseMerge z z
  (merged.x.val, merged.y.val)

#eval
  let s : BraidState := { strands := fun _ => BraidStrand.zero 0, step_count := 0 }
  let s1 := contractedCrossStep s
  s1.strands 0 == BraidStrand.zero 0