import PhysicsLogic.Papers.KaplanSavageWise1998.ContactChains

set_option backward.isDefEq.respectTransparency false

open PhysicsLogic.Papers.KaplanSavageWise1998.ContactChains

example : forall (g : Nat -> Real) (mu cutoff p : Real) (vertices : List Nat),
    0 < cutoff -> 0 < p -> vertices ≠ [] ->
    p * ‖amplitude g mu p vertices‖ =
      ‖((mu / p : Real) : Complex) + Complex.I‖ ^ (vertices.length - 1) *
        (p / cutoff) ^ vertices.sum *
          (vertices.map (fun n => |cutoff ^ n * p ^ (n + 1) * g n|)).prod :=
  PhysicsLogic.Papers.KaplanSavageWise1998.ContactChains.normalized_chain_magnitude_at_scale

-- The original sharper equal-scale identity must remain unchanged.
example : forall (g : Nat -> Real) (cutoff p : Real) (vertices : List Nat),
    0 < cutoff -> 0 < p -> vertices ≠ [] ->
    p * ‖amplitude g p p vertices‖ =
      ‖(1 : Complex) + Complex.I‖ ^ (vertices.length - 1) *
        (p / cutoff) ^ vertices.sum *
          (vertices.map (fun n => |cutoff ^ n * p ^ (n + 1) * g n|)).prod :=
  PhysicsLogic.Papers.KaplanSavageWise1998.ContactChains.normalized_chain_magnitude

-- Common coefficients and constants precede BOTH mu and p. There are no
-- assigned orders or validity predicates on this literal side of the contract.
example : forall (g : Real -> Nat -> Real) (k : Nat -> Real) (last : Nat)
    (cutoff inverseLength rho bound c : Real),
    0 < cutoff -> |inverseLength| <= rho * cutoff -> rho < 1 -> 0 <= bound -> 1 <= c ->
    (forall j, 0 < j -> j <= last -> |cutoff ^ (2 * j - 1) * k j| <= bound) ->
    (forall mu, 0 < mu -> mu <= cutoff -> |inverseLength| <= rho * mu ->
      (inverseLength - mu) * g mu 0 = 1) ->
    (forall mu, 0 < mu -> mu <= cutoff -> |inverseLength| <= rho * mu ->
      forall n, n < last -> (inverseLength - mu) * g mu (n + 1) =
        (Finset.range (n + 1)).sum (fun j => k (n + 1 - j) * g mu j)) ->
    (exists mu p : Real, 0 < mu ∧ mu <= cutoff ∧ |inverseLength| <= rho * mu ∧
      0 < p ∧ p <= c * mu ∧ mu <= c * p) ∧
    (forall mu p, 0 < mu -> mu <= cutoff -> |inverseLength| <= rho * mu ->
      0 < p -> p <= c * mu -> mu <= c * p ->
      (forall vertices : List Nat, vertices ≠ [] ->
        (forall n, n ∈ vertices -> n <= last) ->
        p * ‖amplitude (g mu) mu p vertices‖ <=
          (c + 1) ^ (vertices.length - 1) * (c / (1 - rho)) ^ vertices.length *
            (c * (1 + bound / (1 - rho)) * p / cutoff) ^ vertices.sum) ∧
      (forall loopCut : Nat, exists vertices : List Nat, vertices ≠ [] ∧
        loopCut < vertices.length - 1 ∧ vertices.sum = 0 ∧
        0 < p * ‖amplitude (g mu) mu p vertices‖ ∧
        (1 / (c * (1 + rho))) ^ (loopCut + 2) <=
          p * ‖amplitude (g mu) mu p vertices‖)) ∧
    (forall vertices : List Nat, vertices ≠ [] ->
      (vertices.sum <= 1 ↔ vertices = List.replicate vertices.length 0 ∨
        exists left right : Nat,
          vertices = List.replicate left 0 ++ [1] ++ List.replicate right 0)) :=
  PhysicsLogic.Papers.KaplanSavageWise1998.ContactChains.comparable_matching_window_controls_actual_chains

-- Finite negative scattering length a=-100, nonzero range r0=1, mu=1/20,
-- p=1/10 and Lambda=1 in chosen units. This is an algebraic EFT diagnostic,
-- not a fitted nuclear potential or a physical remainder estimate.
example :
    let g : Real -> Nat -> Real := fun mu n =>
      if n = 0 then 1 / (-1 / 100 - mu) else
      if n = 1 then (1 / 2) / (-1 / 100 - mu) ^ 2 else 0
    g (1 / 20) 0 = -50 / 3 ∧ g (1 / 20) 1 = 1250 / 9 ∧
    (forall vertices : List Nat, vertices ≠ [] ->
      (forall n, n ∈ vertices -> n <= 1) ->
      (1 / 10 : Real) * ‖amplitude (g (1 / 20)) (1 / 20) (1 / 10) vertices‖ <=
        3 ^ (vertices.length - 1) * 4 ^ vertices.length * (2 / 5) ^ vertices.sum) ∧
    (forall loopCut : Nat, exists vertices : List Nat, vertices ≠ [] ∧
      loopCut < vertices.length - 1 ∧ vertices.sum = 0 ∧
      0 < (1 / 10 : Real) * ‖amplitude (g (1 / 20)) (1 / 20) (1 / 10) vertices‖ ∧
      (1 / 3 : Real) ^ (loopCut + 2) <=
        (1 / 10 : Real) * ‖amplitude (g (1 / 20)) (1 / 20) (1 / 10) vertices‖) := by
  dsimp only
  have result := comparable_matching_window_controls_actual_chains
    (fun mu n => if n = 0 then 1 / (-1 / 100 - mu) else
      if n = 1 then (1 / 2) / (-1 / 100 - mu) ^ 2 else 0)
    (fun j => if j = 1 then 1 / 2 else 0) 1 1 (-1 / 100) (1 / 2) (1 / 2) 2
    (by norm_num) (by norm_num) (by norm_num) (by norm_num) (by norm_num)
    (by
      intro j positive retained
      have : j = 1 := by omega
      subst j
      norm_num)
    (by
      intro mu positive _ _
      have nonpole : (-1 / 100 : Real) - mu ≠ 0 := by linarith
      simpa using mul_one_div_cancel nonpole)
    (by
      intro mu positive _ _ n retained
      have nonpole : (-1 / 100 : Real) - mu ≠ 0 := by linarith
      have : n = 0 := by omega
      subst n
      norm_num [Finset.sum_range_succ]
      field_simp)
  have point := result.2.1 (1 / 20) (1 / 10)
    (by norm_num) (by norm_num) (by norm_num) (by norm_num) (by norm_num) (by norm_num)
  refine ⟨by norm_num, by norm_num, ?_⟩
  norm_num at point
  convert point using 1 <;> norm_num

-- With mu held far above p, the claimed uniform leading lower bound fails.
-- The common window and the other comparability inequality still hold.
example :
    ((-1 / 100 : Real) - 1) * (-100 / 101) = 1 ∧
      |(-1 / 100 : Real)| <= (1 / 2) * 1 ∧
      (1 / 100 : Real) <= 2 * 1 ∧ ¬(1 : Real) <= 2 * (1 / 100) ∧
      (1 / 100 : Real) * ‖amplitude (fun _ => -100 / 101) 1 (1 / 100) [0]‖ <
        1 / ((2 : Real) * (1 + 1 / 2)) := by
  norm_num [amplitude, Complex.norm_real, Real.norm_eq_abs]

#print axioms PhysicsLogic.Papers.KaplanSavageWise1998.ContactChains.normalized_chain_magnitude_at_scale
#print axioms PhysicsLogic.Papers.KaplanSavageWise1998.ContactChains.comparable_matching_window_controls_actual_chains