import Mathlib
import Definitions.Def_WeierstrassCurve_Velu
import Definitions.Def_WeierstrassCurve_VeluOrderTwo
import Definitions.Def_WeierstrassCurve_VeluPointMap2
import Definitions.Def_WeierstrassCurve_OddOrderSummingSet

set_option autoImplicit false

noncomputable section

namespace WeierstrassCurve

section XData

variable {R : Type*} [CommRing R] (V : WeierstrassCurve R)

def xVeluT (x : R) : R := 6 * x ^ 2 + V.b₂ * x + V.b₄

def xVeluU (x : R) : R := 4 * x ^ 3 + V.b₂ * x ^ 2 + 2 * V.b₄ * x + V.b₆

def xVeluW (x : R) : R := V.xVeluU x + x * V.xVeluT x

theorem xVeluT_eq_veluT (x y : R) : V.xVeluT x = V.veluT x y := by
  rw [xVeluT, veluT_eq]

theorem xVeluU_eq_veluU {x y : R} (h : V.toAffine.Equation x y) : V.xVeluU x = V.veluU x y := by
  rw [V.veluU_eq_Ψ₂Sq_eval h, xVeluU, Ψ₂Sq]
  simp only [Polynomial.eval_add, Polynomial.eval_mul, Polynomial.eval_pow, Polynomial.eval_C,
    Polynomial.eval_X]

theorem xVeluW_eq_veluW {x y : R} (h : V.toAffine.Equation x y) : V.xVeluW x = V.veluW x y := by
  rw [xVeluW, veluW, V.xVeluU_eq_veluU h, V.xVeluT_eq_veluT x y]

end XData

section Field

variable {L : Type*} [Field L]

def absSum (X : Set L) (g : L → L) : L :=
  haveI := Classical.dec X.Finite
  if X.Finite then ∑ᶠ x ∈ X, g x else 0

theorem absSum_of_finite {X : Set L} (hX : X.Finite) (g : L → L) :
    absSum X g = ∑ x ∈ hX.toFinset, g x := by
  rw [absSum, if_pos hX, finsum_mem_eq_finite_toFinset_sum g hX]

theorem absSum_of_infinite {X : Set L} (hX : X.Infinite) (g : L → L) : absSum X g = 0 := by
  rw [absSum, if_neg hX]

def xVeluCurve (V : WeierstrassCurve L) (X : Set L) : WeierstrassCurve L where
  a₁ := V.a₁
  a₂ := V.a₂
  a₃ := V.a₃
  a₄ := V.a₄ - 5 * absSum X V.xVeluT
  a₆ := V.a₆ - V.b₂ * absSum X V.xVeluT - 7 * absSum X V.xVeluW

def xVeluX (V : WeierstrassCurve L) (X : Set L) (x : L) : L :=
  x + absSum X fun x₀ => V.xVeluT x₀ / (x - x₀) + V.xVeluU x₀ / (x - x₀) ^ 2

@[simp] theorem xVeluCurve_a₁ (V : WeierstrassCurve L) (X : Set L) : (V.xVeluCurve X).a₁ = V.a₁ := rfl
@[simp] theorem xVeluCurve_a₂ (V : WeierstrassCurve L) (X : Set L) : (V.xVeluCurve X).a₂ = V.a₂ := rfl
@[simp] theorem xVeluCurve_a₃ (V : WeierstrassCurve L) (X : Set L) : (V.xVeluCurve X).a₃ = V.a₃ := rfl
theorem xVeluCurve_a₄ (V : WeierstrassCurve L) (X : Set L) :
    (V.xVeluCurve X).a₄ = V.a₄ - 5 * absSum X V.xVeluT := rfl
theorem xVeluCurve_a₆ (V : WeierstrassCurve L) (X : Set L) :
    (V.xVeluCurve X).a₆ = V.a₆ - V.b₂ * absSum X V.xVeluT - 7 * absSum X V.xVeluW := rfl

def twoTorsionY (V : WeierstrassCurve L) (x₀ : L) : L := -(V.a₁ * x₀ + V.a₃) / 2

def xVeluG (V : WeierstrassCurve L) (x₀ : L) : L := V.veluGx x₀ (V.twoTorsionY x₀)

def twoVeluCurve (V : WeierstrassCurve L) (X : Set L) : WeierstrassCurve L where
  a₁ := V.a₁
  a₂ := V.a₂
  a₃ := V.a₃
  a₄ := V.a₄ - 5 * absSum X V.xVeluG
  a₆ := V.a₆ - V.b₂ * absSum X V.xVeluG - 7 * absSum X fun x => x * V.xVeluG x

def twoVeluX (V : WeierstrassCurve L) (X : Set L) (x : L) : L :=
  x + absSum X fun x₀ => V.xVeluG x₀ / (x - x₀)

@[simp] theorem twoVeluCurve_a₁ (V : WeierstrassCurve L) (X : Set L) : (V.twoVeluCurve X).a₁ = V.a₁ := rfl
@[simp] theorem twoVeluCurve_a₂ (V : WeierstrassCurve L) (X : Set L) : (V.twoVeluCurve X).a₂ = V.a₂ := rfl
@[simp] theorem twoVeluCurve_a₃ (V : WeierstrassCurve L) (X : Set L) : (V.twoVeluCurve X).a₃ = V.a₃ := rfl
theorem twoVeluCurve_a₄ (V : WeierstrassCurve L) (X : Set L) :
    (V.twoVeluCurve X).a₄ = V.a₄ - 5 * absSum X V.xVeluG := rfl
theorem twoVeluCurve_a₆ (V : WeierstrassCurve L) (X : Set L) :
    (V.twoVeluCurve X).a₆ = V.a₆ - V.b₂ * absSum X V.xVeluG - 7 * absSum X fun x => x * V.xVeluG x :=
  rfl

variable [DecidableEq L]

def kernelXSet (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (ℓ : ℕ) : Set L :=
  (fun P => P.coordsOrZero.1) '' {P | P ∈ H ∧ P ≠ 0 ∧ ℓ • P = 0}

def coKernelXSet (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (ℓ : ℕ) : Set L :=
  (fun P => P.coordsOrZero.1) '' {P | P ∈ H ∧ ℓ • P ≠ 0}

def stepCurve (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (ℓ : ℕ) :
    WeierstrassCurve L :=
  if ℓ = 2 then V.twoVeluCurve (V.kernelXSet H 2) else V.xVeluCurve (V.kernelXSet H ℓ)

def stepX (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (ℓ : ℕ) (x : L) : L :=
  if ℓ = 2 then V.twoVeluX (V.kernelXSet H 2) x else V.xVeluX (V.kernelXSet H ℓ) x

theorem stepCurve_two (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) :
    V.stepCurve H 2 = V.twoVeluCurve (V.kernelXSet H 2) :=
  if_pos rfl

theorem stepCurve_of_ne_two (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) {ℓ : ℕ}
    (hℓ : ℓ ≠ 2) : V.stepCurve H ℓ = V.xVeluCurve (V.kernelXSet H ℓ) :=
  if_neg hℓ

theorem stepX_two (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (x : L) :
    V.stepX H 2 x = V.twoVeluX (V.kernelXSet H 2) x :=
  if_pos rfl

theorem stepX_of_ne_two (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) {ℓ : ℕ}
    (hℓ : ℓ ≠ 2) (x : L) : V.stepX H ℓ x = V.xVeluX (V.kernelXSet H ℓ) x :=
  if_neg hℓ

def subgroupOfX (V' : WeierstrassCurve L) (Y : Set L) : AddSubgroup V'.toAffine.Point :=
  AddSubgroup.closure {P | P ≠ 0 ∧ P.coordsOrZero.1 ∈ Y}

def stepSubgroup (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (ℓ : ℕ) :
    AddSubgroup (V.stepCurve H ℓ).toAffine.Point :=
  subgroupOfX (V.stepCurve H ℓ) (V.stepX H ℓ '' V.coKernelXSet H ℓ)

variable (L) in

abbrev CQJState : Type _ := (Σ V : WeierstrassCurve L, AddSubgroup V.toAffine.Point) × ℕ

def cqjStep (s : CQJState L) : CQJState L :=
  ⟨⟨s.1.1.stepCurve s.1.2 s.2.minFac, s.1.1.stepSubgroup s.1.2 s.2.minFac⟩, s.2 / s.2.minFac⟩

def cqjIterate (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (N : ℕ) : CQJState L :=
  cqjStep^[N.primeFactorsList.length] ⟨⟨V, H⟩, N⟩

def cyclicQuotientCurve (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (N : ℕ) :
    WeierstrassCurve L :=
  (V.cqjIterate H N).1.1

def cyclicQuotientJ (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (N : ℕ) : L :=
  (V.cyclicQuotientCurve H N).c₄ ^ 3 / (V.cyclicQuotientCurve H N).Δ

theorem cqjIterate_def (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (N : ℕ) :
    V.cqjIterate H N = cqjStep^[N.primeFactorsList.length] ⟨⟨V, H⟩, N⟩ := rfl

theorem cyclicQuotientCurve_def (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point)
    (N : ℕ) : V.cyclicQuotientCurve H N = (V.cqjIterate H N).1.1 := rfl

theorem cyclicQuotientJ_def (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (N : ℕ) :
    V.cyclicQuotientJ H N =
      (V.cyclicQuotientCurve H N).c₄ ^ 3 / (V.cyclicQuotientCurve H N).Δ := rfl

theorem cqjStep_apply (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (n : ℕ) :
    cqjStep ⟨⟨V, H⟩, n⟩ =
      ⟨⟨V.stepCurve H n.minFac, V.stepSubgroup H n.minFac⟩, n / n.minFac⟩ := rfl

theorem cqjIterate_one (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) :
    V.cqjIterate H 1 = ⟨⟨V, H⟩, 1⟩ := by
  simp [cqjIterate, Nat.primeFactorsList_one]

theorem cyclicQuotientCurve_one (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) :
    V.cyclicQuotientCurve H 1 = V := by
  rw [cyclicQuotientCurve_def, cqjIterate_one]

theorem cyclicQuotientJ_one (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) :
    V.cyclicQuotientJ H 1 = V.c₄ ^ 3 / V.Δ := by
  rw [cyclicQuotientJ_def, cyclicQuotientCurve_one]

theorem length_primeFactorsList_eq_succ {N : ℕ} (hN : 2 ≤ N) :
    N.primeFactorsList.length = (N / N.minFac).primeFactorsList.length + 1 := by
  obtain ⟨k, rfl⟩ : ∃ k, N = k + 2 := ⟨N - 2, by omega⟩
  rw [Nat.primeFactorsList, List.length_cons]

theorem cqjIterate_eq_of_two_le (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point)
    {N : ℕ} (hN : 2 ≤ N) :
    V.cqjIterate H N =
      (V.stepCurve H N.minFac).cqjIterate (V.stepSubgroup H N.minFac) (N / N.minFac) := by
  rw [cqjIterate, length_primeFactorsList_eq_succ hN, Function.iterate_succ_apply]
  rfl

theorem cyclicQuotientCurve_eq_of_two_le (V : WeierstrassCurve L)
    (H : AddSubgroup V.toAffine.Point) {N : ℕ} (hN : 2 ≤ N) :
    V.cyclicQuotientCurve H N =
      (V.stepCurve H N.minFac).cyclicQuotientCurve (V.stepSubgroup H N.minFac) (N / N.minFac) := by
  rw [cyclicQuotientCurve_def, cqjIterate_eq_of_two_le V H hN, cyclicQuotientCurve_def]

theorem cyclicQuotientJ_eq_of_two_le (V : WeierstrassCurve L)
    (H : AddSubgroup V.toAffine.Point) {N : ℕ} (hN : 2 ≤ N) :
    V.cyclicQuotientJ H N =
      (V.stepCurve H N.minFac).cyclicQuotientJ (V.stepSubgroup H N.minFac) (N / N.minFac) := by
  rw [cyclicQuotientJ_def, cyclicQuotientCurve_eq_of_two_le V H hN, cyclicQuotientJ_def]

theorem cyclicQuotientJ_eq_j (V : WeierstrassCurve L) (H : AddSubgroup V.toAffine.Point) (N : ℕ)
    [(V.cyclicQuotientCurve H N).IsElliptic] :
    V.cyclicQuotientJ H N = (V.cyclicQuotientCurve H N).j := by
  rw [cyclicQuotientJ_def, WeierstrassCurve.j, div_eq_mul_inv, mul_comm, Units.val_inv_eq_inv_val,
    WeierstrassCurve.coe_Δ']

end Field

end WeierstrassCurve

end