IndisputableMonolith.ClassicalBridge.Fluids.Simulation2D

   1import IndisputableMonolith.ClassicalBridge.Fluids.Galerkin2D
   2import IndisputableMonolith.ClassicalBridge.Fluids.LNALSemantics
   3
   4namespace IndisputableMonolith.ClassicalBridge.Fluids
   5
   6/-!
   7# Simulation2D (Milestone M3)
   8
   9This file states the *one-step simulation* bridge between:
  10- the discrete 2D Galerkin model (`Galerkin2D`), and
  11- a spatial LNAL execution semantics (`LNALSemantics.independent`) applied to an encoded field.
  12
  13At this milestone we do **not** attempt to prove the analytic correctness of the simulation;
  14instead we package the required claim as an explicit `SimulationHypothesis` (no `axiom`, no `sorry`).
  15-/
  16
  17open IndisputableMonolith.ClassicalBridge.Fluids
  18
  19namespace Simulation2D
  20
  21open IndisputableMonolith.LNAL
  22open IndisputableMonolith.ClassicalBridge.Fluids.LNALSemantics
  23open IndisputableMonolith.ClassicalBridge.Fluids.Encoding
  24
  25/-!
  26## A concrete (nontrivial) 1-step simulation instance
  27
  28The current spatial semantics is **independent per voxel**, so the only operations we can
  29fully simulate without additional hypotheses are those that act locally on the encoding.
  30
  31Here we provide a simple, explicit example:
  32- LNAL program: constant `FOLD 1` (increments `nuPhi` by `+1`, with `clampI32`).
  33- Discrete step: update each Galerkin coefficient by replacing it with an **integer** whose
  34  sign matches `coeffSign` and whose magnitude is `clampI32 (coeffMag x + 1)`.
  35
  36This is not yet a Navier–Stokes time step; it’s a proved, quantitative “register-step ↔ encoded-step”
  37example that we can build on once multi-voxel coupling is introduced.
  38-/
  39
  40/-- A nontrivial LNAL program: increment `nuPhi` by `+1` at every instruction pointer. -/
  41@[simp] def foldPlusOneProgram : LProgram :=
  42  fun _ => LInstr.fold 1
  43
  44/-- One-voxel semantics for `foldPlusOneProgram`: increment `nuPhi` by `+1` (clamped); leave `aux5` unchanged. -/
  45lemma voxelStep_foldPlusOneProgram (v : LNALVoxel) :
  46    voxelStep foldPlusOneProgram v = ({ v.1 with nuPhi := clampI32 (v.1.nuPhi + 1) }, v.2) := by
  47  rcases v with ⟨r6, a5⟩
  48  simp [LNALSemantics.voxelStep, foldPlusOneProgram, lStep, -clampI32]
  49
  50/-- The corresponding discrete step on encoded Galerkin coefficients. -/
  51noncomputable def foldPlusOneStep {N : ℕ} (u : GalerkinState N) : GalerkinState N :=
  52  WithLp.toLp 2 (fun i : ((modes N) × Fin 2) =>
  53    let x : ℝ := u i
  54    let m : Int := clampI32 (coeffMag x + 1)
  55    let z : Int := (coeffSign x) * m
  56    (z : ℝ))
  57
  58lemma floor_abs_intCast (z : Int) : Int.floor (|((z : ℝ))|) = |z| := by
  59  have habs : |((z : ℝ))| = ((|z| : Int) : ℝ) := by
  60    simp
  61  rw [habs]
  62  simpa using (Int.floor_intCast (R := ℝ) (z := |z|))
  63
  64lemma cast_lt_zero_iff {z : Int} : ((z : ℝ) < 0) ↔ z < 0 := by
  65  constructor
  66  · intro hz
  67    by_contra h
  68    have : (0 : Int) ≤ z := le_of_not_gt h
  69    have : (0 : ℝ) ≤ (z : ℝ) := by exact_mod_cast this
  70    exact (not_lt_of_ge this) hz
  71  · intro hz
  72    exact_mod_cast hz
  73
  74lemma clampI32_pos_of_ge_one {x : Int} (hx : 1 ≤ x) : 0 < clampI32 x := by
  75  -- From `x ≥ 1`, the negative-saturation branch is impossible; then either we saturate high,
  76  -- or we return `x` itself. In both cases the result is positive.
  77  have hx0 : (0 : Int) ≤ x := le_trans (by decide : (0 : Int) ≤ 1) hx
  78  have hnot : ¬ x ≤ (-i32Bound) := by
  79    have hi : (-i32Bound : Int) < 0 := by
  80      have : (0 : Int) < i32Bound := by decide
  81      linarith
  82    exact not_le_of_gt (lt_of_lt_of_le hi hx0)
  83  -- Convert `hnot` to the numeral form used by definitional unfolding of `i32Bound`.
  84  have hnot' : ¬ x ≤ (-2147483648 : Int) := by simpa using hnot
  85  dsimp [clampI32]
  86  rw [if_neg hnot']
  87  by_cases h₂ : (2147483648 : Int) ≤ x
  88  · rw [if_pos h₂]
  89    have : (1 : Int) < (2147483648 : Int) := by decide
  90    linarith
  91  · rw [if_neg h₂]
  92    exact lt_of_lt_of_le (by decide : (0 : Int) < 1) hx
  93
  94lemma coeffMag_foldPlusOneStep {N : ℕ} (u : GalerkinState N) (i : (modes N) × Fin 2) :
  95    coeffMag ((foldPlusOneStep u) i) = clampI32 (coeffMag (u i) + 1) := by
  96  classical
  97  set x : ℝ := u i
  98  have hx1 : (1 : Int) ≤ coeffMag x + 1 := by
  99    have hx0 : (0 : Int) ≤ coeffMag x := (Int.floor_nonneg (a := |x|)).2 (abs_nonneg x)
 100    linarith
 101  -- Let `m` be the (positive) clamped magnitude and `z` the signed integer we encode as a real.
 102  set m : Int := clampI32 (coeffMag x + 1)
 103  set z : Int := coeffSign x * m
 104  have hmpos : 0 < m := by
 105    simpa [m] using (clampI32_pos_of_ge_one (x := coeffMag x + 1) hx1)
 106  have hm0 : 0 ≤ m := le_of_lt hmpos
 107  -- Expand `coeffMag` and reduce to a floor/abs fact on an integer cast.
 108  -- `foldPlusOneStep u i` is literally `(z : ℝ)`.
 109  have hcore :
 110      coeffMag ((foldPlusOneStep u) i) = Int.floor (|((z : ℝ))|) := by
 111    simp [coeffMag, foldPlusOneStep, x, m, z, -clampI32]
 112  -- Now compute the floor of the abs cast and simplify `|z|`.
 113  calc
 114    coeffMag ((foldPlusOneStep u) i)
 115        = Int.floor (|((z : ℝ))|) := hcore
 116    _ = |z| := floor_abs_intCast z
 117    _ = m := by
 118      by_cases hx : x < 0
 119      · have hz : z = -m := by simp [z, coeffSign, hx]
 120        rw [hz]
 121        -- `| -m | = |m| = m` since `m ≥ 0`.
 122        simp [abs_of_nonneg hm0]
 123      · have hz : z = m := by simp [z, coeffSign, hx]
 124        rw [hz]
 125        simp [abs_of_nonneg hm0]
 126
 127lemma coeffSign_foldPlusOneStep {N : ℕ} (u : GalerkinState N) (i : (modes N) × Fin 2) :
 128    coeffSign ((foldPlusOneStep u) i) = coeffSign (u i) := by
 129  classical
 130  set x : ℝ := u i
 131  have hx1 : (1 : Int) ≤ coeffMag x + 1 := by
 132    have hx0 : (0 : Int) ≤ coeffMag x := (Int.floor_nonneg (a := |x|)).2 (abs_nonneg x)
 133    linarith
 134  have hmpos : 0 < clampI32 (coeffMag x + 1) :=
 135    clampI32_pos_of_ge_one (x := coeffMag x + 1) hx1
 136  by_cases hx : x < 0
 137  ·
 138    -- `x < 0` ⇒ `coeffSign x = -1`, and the step coefficient is negative (since the magnitude is positive).
 139    have hzInt : (coeffSign x) * clampI32 (coeffMag x + 1) < 0 := by
 140      have hneg : -(clampI32 (coeffMag x + 1)) < 0 := neg_neg_of_pos hmpos
 141      calc
 142        (coeffSign x) * clampI32 (coeffMag x + 1)
 143            = (-1) * clampI32 (coeffMag x + 1) := by simp [coeffSign, hx]
 144        _ = -(clampI32 (coeffMag x + 1)) := by ring
 145        _ < 0 := hneg
 146    have hzReal :
 147        (((coeffSign x) * clampI32 (coeffMag x + 1) : Int) : ℝ) < 0 :=
 148      (cast_lt_zero_iff (z := (coeffSign x) * clampI32 (coeffMag x + 1))).2 hzInt
 149    have hy : (foldPlusOneStep u) i < 0 := by
 150      simpa [foldPlusOneStep, x] using hzReal
 151    simp [coeffSign, hx, hy, x]
 152  ·
 153    -- `¬ x < 0` ⇒ `coeffSign x = 1`, and the step coefficient is nonnegative.
 154    have hzInt : (0 : Int) ≤ (coeffSign x) * clampI32 (coeffMag x + 1) := by
 155      have hm0 : (0 : Int) ≤ clampI32 (coeffMag x + 1) := le_of_lt hmpos
 156      simpa [coeffSign, hx] using hm0
 157    have hzReal :
 158        (0 : ℝ) ≤ (((coeffSign x) * clampI32 (coeffMag x + 1) : Int) : ℝ) := by
 159      exact_mod_cast hzInt
 160    have hy : ¬ (foldPlusOneStep u) i < 0 := by
 161      have hy0 : (0 : ℝ) ≤ (foldPlusOneStep u) i := by
 162        simpa [foldPlusOneStep, x] using hzReal
 163      exact not_lt_of_ge hy0
 164    simp [coeffSign, hx, hy, x]
 165
 166lemma decide_lt_zero_foldPlusOneStep {N : ℕ} (u : GalerkinState N) (i : (modes N) × Fin 2) :
 167    decide ((foldPlusOneStep u) i < 0) = decide (u i < 0) := by
 168  classical
 169  set x : ℝ := u i
 170  have hx1 : (1 : Int) ≤ coeffMag x + 1 := by
 171    have hx0 : (0 : Int) ≤ coeffMag x := (Int.floor_nonneg (a := |x|)).2 (abs_nonneg x)
 172    linarith
 173  have hmpos : 0 < clampI32 (coeffMag x + 1) :=
 174    clampI32_pos_of_ge_one (x := coeffMag x + 1) hx1
 175  by_cases hx : x < 0
 176  ·
 177    have hzInt : (coeffSign x) * clampI32 (coeffMag x + 1) < 0 := by
 178      have hneg : -(clampI32 (coeffMag x + 1)) < 0 := neg_neg_of_pos hmpos
 179      calc
 180        (coeffSign x) * clampI32 (coeffMag x + 1)
 181            = (-1) * clampI32 (coeffMag x + 1) := by simp [coeffSign, hx]
 182        _ = -(clampI32 (coeffMag x + 1)) := by ring
 183        _ < 0 := hneg
 184    have hzReal :
 185        (((coeffSign x) * clampI32 (coeffMag x + 1) : Int) : ℝ) < 0 :=
 186      (cast_lt_zero_iff (z := (coeffSign x) * clampI32 (coeffMag x + 1))).2 hzInt
 187    have hy : (foldPlusOneStep u) i < 0 := by
 188      simpa [foldPlusOneStep, x] using hzReal
 189    simp [x, hx, hy]
 190  ·
 191    have hzInt : (0 : Int) ≤ (coeffSign x) * clampI32 (coeffMag x + 1) := by
 192      have hm0 : (0 : Int) ≤ clampI32 (coeffMag x + 1) := le_of_lt hmpos
 193      simpa [coeffSign, hx] using hm0
 194    have hzReal :
 195        (0 : ℝ) ≤ (((coeffSign x) * clampI32 (coeffMag x + 1) : Int) : ℝ) := by
 196      exact_mod_cast hzInt
 197    have hy : ¬ (foldPlusOneStep u) i < 0 := by
 198      have hy0 : (0 : ℝ) ≤ (foldPlusOneStep u) i := by
 199        simpa [foldPlusOneStep, x] using hzReal
 200      exact not_lt_of_ge hy0
 201    simp [x, hx, hy]
 202
 203lemma voxelStep_foldPlusOne_encodeIndex {N : ℕ} (u : GalerkinState N) (i : (modes N) × Fin 2) :
 204    voxelStep foldPlusOneProgram (encodeIndex u i) = encodeIndex (foldPlusOneStep u) i := by
 205  classical
 206  -- use the closed-form voxel semantics for `FOLD 1`, then rewrite encoding fields
 207  simp [voxelStep_foldPlusOneProgram, encodeIndex, coeffMag_foldPlusOneStep, coeffSign_foldPlusOneStep,
 208    decide_lt_zero_foldPlusOneStep, -clampI32]
 209
 210/-!
 211## Decoding the encoding (for weaker, decode-based simulation)
 212
 213Exact commutation at the *encoding* level is too strict for decay steps that can hit `0` and flip
 214sign bits (`sigma/phiE`) when re-encoded. A practical workaround is to compare **decoded**
 215coefficients after stepping.
 216-/
 217
 218/-- Decode a nonnegative magnitude from a register `nuPhi`.
 219
 220We interpret negative values as `0` (via `Int.toNat`) so that a decrement step can “saturate at 0”
 221even if the underlying register becomes negative. -/
 222def decodeMag (z : Int) : Int :=
 223  Int.ofNat (Int.toNat z)
 224
 225/-- Decode a single real coefficient from a voxel: sign from `sigma`, magnitude from `nuPhi`. -/
 226noncomputable def decodeCoeff (v : LNALVoxel) : ℝ :=
 227  ((v.1.sigma * decodeMag v.1.nuPhi : Int) : ℝ)
 228
 229/-- Decode an `LNALField` (of the expected Galerkin encoding length) back to a Galerkin state.
 230
 231This is a left-inverse of `encodeGalerkin2D` only up to the coarse quantization used by `coeffMag`.
 232-/
 233noncomputable def decodeGalerkin2D {N : ℕ} (field : LNALField)
 234    (hsize : field.size = Fintype.card ((modes N) × Fin 2)) : GalerkinState N :=
 235  WithLp.toLp 2 (fun i : ((modes N) × Fin 2) =>
 236    let j : Fin (Fintype.card ((modes N) × Fin 2)) := (Fintype.equivFin ((modes N) × Fin 2)) i
 237    decodeCoeff (field[(Fin.cast hsize.symm j)]))
 238
 239/-- Hypothesis: one LNAL spatial step simulates one discrete Galerkin step (exactly). -/
 240structure SimulationHypothesis (N : ℕ) where
 241  /-- The LNAL program used for the simulation. -/
 242  P : LProgram
 243  /-- The discrete (time-stepping) map on Galerkin states. -/
 244  step : GalerkinState N → GalerkinState N
 245  /-- One-step commutation: execute then encode = encode then step. -/
 246  comm :
 247      ∀ u : GalerkinState N,
 248        (independent.step P (encodeGalerkin2D u)) = encodeGalerkin2D (step u)
 249
 250/-- Trivial simulation: use the `LISTEN noop` LNAL program and take the discrete step as `id`. -/
 251@[simp] def SimulationHypothesis.trivial (N : ℕ) : SimulationHypothesis N :=
 252  { P := listenNoopProgram
 253    step := id
 254    comm := by
 255      intro u
 256      simp }
 257
 258/-- A concrete, nontrivial simulation instance: `FOLD 1` corresponds to `foldPlusOneStep`. -/
 259noncomputable def SimulationHypothesis.foldPlusOne (N : ℕ) : SimulationHypothesis N :=
 260  { P := foldPlusOneProgram
 261    step := foldPlusOneStep
 262    comm := by
 263      intro u
 264      classical
 265      -- Prove array equality by pointwise equality.
 266      refine Array.ext (by simp [LNALSemantics.independent, encodeGalerkin2D]) ?_
 267      intro j hj₁ hj₂
 268      have hj : j < Fintype.card ((modes N) × Fin 2) := by
 269        simpa [encodeGalerkin2D] using hj₂
 270      let jFin : Fin (Fintype.card ((modes N) × Fin 2)) := ⟨j, hj⟩
 271      -- Reduce to the single-voxel commutation lemma.
 272      simpa [LNALSemantics.independent, encodeGalerkin2D, foldPlusOneProgram] using
 273        (show voxelStep foldPlusOneProgram (encodeIndex u ((Fintype.equivFin ((modes N) × Fin 2)).symm jFin))
 274            = encodeIndex (foldPlusOneStep u) ((Fintype.equivFin ((modes N) × Fin 2)).symm jFin) from by
 275            simpa using
 276              (voxelStep_foldPlusOne_encodeIndex (u := u)
 277                (i := (Fintype.equivFin ((modes N) × Fin 2)).symm jFin))) }
 278
 279/-- The one-step simulation lemma (directly from `SimulationHypothesis.comm`). -/
 280theorem simulation_one_step {N : ℕ} (H : SimulationHypothesis N) (u : GalerkinState N) :
 281    independent.step H.P (encodeGalerkin2D u) = encodeGalerkin2D (H.step u) :=
 282  H.comm u
 283
 284/-!
 285## Decode-based simulation (weaker notion) and a diagonal decay step
 286
 287For some local operations (e.g. sign-preserving decay of a magnitude register), exact commutation
 288`encode(step u) = step(encode u)` can fail due to sign-bit conventions at `0`.
 289
 290We therefore also provide a **decode-based** simulation notion:
 291we compare decoded coefficients after executing the LNAL step.
 292-/
 293
 294/-- Hypothesis: one LNAL spatial step simulates one discrete Galerkin step **after decoding**. -/
 295structure DecodedSimulationHypothesis (N : ℕ) where
 296  /-- The LNAL program used for the simulation. -/
 297  P : LProgram
 298  /-- The discrete (time-stepping) map on Galerkin states (on decoded/quantized coefficients). -/
 299  step : GalerkinState N → GalerkinState N
 300  /-- One-step commutation after decoding. -/
 301  comm :
 302      ∀ u : GalerkinState N,
 303        decodeGalerkin2D (N := N)
 304          (field := independent.step P (encodeGalerkin2D u))
 305          (hsize := by simp [LNALSemantics.independent, encodeGalerkin2D])
 306          = step u
 307
 308/-- The one-step decoded simulation lemma (directly from `DecodedSimulationHypothesis.comm`). -/
 309theorem decoded_simulation_one_step {N : ℕ} (H : DecodedSimulationHypothesis N) (u : GalerkinState N) :
 310    decodeGalerkin2D (N := N)
 311      (field := independent.step H.P (encodeGalerkin2D u))
 312      (hsize := by simp [LNALSemantics.independent, encodeGalerkin2D])
 313      = H.step u :=
 314  H.comm u
 315
 316/-- A constant LNAL program: decrement `nuPhi` by `1` (via `FOLD (-1)`). -/
 317@[simp] def foldMinusOneProgram : LProgram :=
 318  fun _ => LInstr.fold (-1)
 319
 320/-- One-voxel semantics for `foldMinusOneProgram`: decrement `nuPhi` by `1` (clamped); leave `aux5` unchanged. -/
 321lemma voxelStep_foldMinusOneProgram (v : LNALVoxel) :
 322    voxelStep foldMinusOneProgram v = ({ v.1 with nuPhi := clampI32 (v.1.nuPhi + (-1)) }, v.2) := by
 323  rcases v with ⟨r6, a5⟩
 324  simp [LNALSemantics.voxelStep, foldMinusOneProgram, lStep, -clampI32]
 325
 326/-- The corresponding discrete step on **decoded** Galerkin coefficients.
 327
 328We mimic the VM update on `nuPhi` (via `clampI32`) and then decode a nonnegative magnitude using `decodeMag`. -/
 329noncomputable def foldMinusOneDecodedStep {N : ℕ} (u : GalerkinState N) : GalerkinState N :=
 330  WithLp.toLp 2 (fun i : ((modes N) × Fin 2) =>
 331    let x : ℝ := u i
 332    let m : Int := decodeMag (clampI32 (coeffMag x + (-1)))
 333    let z : Int := (coeffSign x) * m
 334    (z : ℝ))
 335
 336/-- Decoding after a single `FOLD (-1)` step matches the corresponding decoded discrete step. -/
 337lemma decodeCoeff_voxelStep_foldMinusOne_encodeIndex {N : ℕ} (u : GalerkinState N) (i : (modes N) × Fin 2) :
 338    decodeCoeff (voxelStep foldMinusOneProgram (encodeIndex u i)) = (foldMinusOneDecodedStep u) i := by
 339  classical
 340  -- both sides reduce to the same signed-integer expression
 341  simp [decodeCoeff, foldMinusOneDecodedStep, decodeMag, voxelStep_foldMinusOneProgram, encodeIndex, -clampI32]
 342
 343/-- A concrete, proved decay-step simulation instance (after decoding). -/
 344noncomputable def DecodedSimulationHypothesis.foldMinusOne (N : ℕ) : DecodedSimulationHypothesis N :=
 345  { P := foldMinusOneProgram
 346    step := foldMinusOneDecodedStep
 347    comm := by
 348      intro u
 349      classical
 350      -- Prove equality of Galerkin states by pointwise equality of coefficients.
 351      ext i
 352      -- Reduce to the single-voxel lemma `decodeCoeff_voxelStep_foldMinusOne_encodeIndex`.
 353      simp [decodeGalerkin2D, LNALSemantics.independent, encodeGalerkin2D,
 354        decodeCoeff_voxelStep_foldMinusOne_encodeIndex] }
 355
 356end Simulation2D
 357
 358end IndisputableMonolith.ClassicalBridge.Fluids
 359